Why does the PDE toolbox construct a Stiffness Matrix with negative Eigenvalues?

Question

Lukas C 2023-5-11

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1961469-why-does-the-pde-toolbox-construct-a-stiffness-matrix-with-negative-eigenvalues

回答： Karan Singh 2023-8-30

Dear Community,

I created a simple geometry using the PDE toolbox, set the boundary conditions and meshed it. When trying to assemble the finite element matrices, using assembleFEMatrices, I get a stiffness matrix K, of which the first few eigenvalues - depending on the size of the matrix - are small negative numbers. This, of course, means that K is not positive definite and thus cannot be factorized using Cholesky, which is what I ultimately need. The assembled mass matrix has no such problem.

My question is: Is there a problem with the pde toolbox or is there something I am missing? From my understanding a stiffness matrix should always at least be semipositive definite. Is there a way to force the matrix to be changed in such a way that those small negative eigenvalues will be small positive eigenvalues?

Thank you in advance.

beam = createpde("structural", "static-solid");

%% Create Geometry

gm = multicuboid(5,1,1);

beam.Geometry = gm;

%Material Properties

structuralProperties(beam, "Cell", 1, "YoungsModulus",210e9, "PoissonsRatio",0.3, "MassDensity",7800);

%Boundary Conditions

structuralBC(beam, "Face", 5 ,"Constraint","fixed");

force = 1000; %apply force

structuralBoundaryLoad(beam, "Face",3, "SurfaceTraction", [0;0;-force]);

%% Generate Mesh

generateMesh(beam, "GeometricOrder","linear","Hmax",1);

%compute the pde solutions

results = solve(beam);

%Assemble FEM-Matrices

FEM = assembleFEMatrices(beam,"MK");

pdeplot3D(beam)

eigenK = eig(FEM.K);

eigenM = eig(FEM.M);

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Karan Singh 2023-8-30

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1961469-why-does-the-pde-toolbox-construct-a-stiffness-matrix-with-negative-eigenvalues#answer_1297141

Hey Lukas, you have correctly pointed out that a stiffness matrix should always be at least semipositive definite. Regarding your question about whether there is a problem with the “PDE Toolbox” or if something is missing, it is difficult to determine without further information. However, it is more likely that the issue lies in the specific problem being modelled or the boundary conditions applied, rather than a problem with the “PDE Toolbox” itself.

You have also asked if there is a way to force the “stiffness matrix” to be changed so that the small negative eigenvalues become small positive eigenvalues. It is generally not recommended to arbitrarily modify the “stiffness matrix” to change the eigenvalues. Instead, it is better to investigate and address the underlying issues causing the negative eigenvalues.

Here are a few suggestions to address this problem:

Check mesh quality: Ensure that the mesh generated by the “generateMesh” function is of good quality. Poor mesh quality, such as distorted elements or highly skewed elements, can lead to inaccurate or unstable results. You can try refining the mesh using a smaller “Hmax” value or using a different meshing algorithm.
Verify boundary conditions: Double-check the boundary conditions you have applied to the problem. Incorrect or inconsistent boundary conditions can lead to non-physical or unstable solutions. Make sure that the fixed constraint and applied force are correctly defined on the appropriate boundaries.
Consider different element types: The default element type used by the “PDE Toolbox” is linear tetrahedral elements. In some cases, using higher-order elements or a different element type (e.g., quadratic elements) can lead to more accurate and stable results. You can try changing the element type using the “generateMesh” function's options.

Hope it helps!