# Why does matrix inversion cause extremely large values to appear in the inverted matrix?

10 views (last 30 days)
raymond bryant on 18 Dec 2021
Commented: raymond bryant on 18 Dec 2021
I'm working on making an fea program where I have a stiffness matrix thats 158 x 158.
When I do \ to find the displacements, I get a sparse vector with several values around 1E16.
I inspected my sitffness matrix and found that most values were in the range of 0.001 to .500
So I inverted the matrix to see what kind of values would appear. I found that my last 4 rows and columns of the inverted stiffness matrix
are all around 1E16.
What could be causing this behavior?

Chunru on 18 Dec 2021
Edited: Chunru on 18 Dec 2021
Consider ax = b where all variables are scaler, what will happen if a is very samll? Then x wil be large
In linear system Ax = b, if A is not well conditioned, in a sense that its inverse can be very large, x = A\b may have some very large elements. You can use "cond(A)" to check if A is well conditioned or not in matlab.
Chunru on 18 Dec 2021
A non-singulare square matrix is full rank. The rank should be size of the matrix.

John D'Errico on 18 Dec 2021
As others have said, a singular stiffness matrix typically means one of several things in a Finite Element problem.
It might be that you made a mistake in the creation of your global stiffness matrix. It happens. Since we cannot possibly know what you did wrong there, then only you can go through your code to fix that.
But perhaps the most common error made in FEA is to insufficiently constrain the problem. That is, if your entire system can freely rotate or translate in some direction as a RIGID body? Then the stiffness matrix will be singular. Remember that FEA works by minimizing the potential energy of the object, so if an enrgy free rotation of translatino exists, then there will be infinitely many solutions. You will then have a singularity. The clue is often a simple one, where the rank is typically one less than the size of the matrix.
Other possibilities are less simple, sometimes reducing to a divide by zero. For example, suppose I designed a very simple truss that connnects nodes a-b, b-c, c-d, but where nodes b and c are coincident? Then the truss member between b amd c will have zero length. Since the stiffness of that truss is found by dividing by the length of that truss member, you have a divide by zero. Other cases might be a 2-d FEA problem composed of triangular elements. But if some of those triangles are constrained to have all nodes in a straight line, then they have zero ara. And again, the stiffness matrix will divide by that area.
There are other ways I could surely come up with to create a numerically singular stiffness matrix, even though it is still mathematically well-posed. But they might involve systems where some elements were made of titanium, and others of polystyrene, so extremely varying stiffnesses.
Regardless, you CANNOT invert a singular stiffness matrix.
raymond bryant on 18 Dec 2021
Okay, thank you very much. I think part of my problem is I've only work on truss FEA systems and this is a beam. So I need to add in the rotation dof's for each node. I might also use a different method for adding the boundary conditions to my K matrix.

R2021b

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by