Error due to “zero” eigenvalues
I’m following this documentation: Model an Excavator Dipper Arm as a Flexible Body – MATLAB & Simulink (mathworks.com).
I copied the entire code and it works with the Dipper.STL file, but for every other STL file i import, I get the ‘zero eigenvalue" error mentioned in the title. The simplest STL file I used was a cube modeled in spaceclaim (10 x 10 x 60 mm). I’ve varied the mesh fineness to be on a spectrum similar to the dipper arm.
I’ve looked at all settings and not figured out what might be wrong. I also added a fixed support for the first frame; didn’t make a difference. It’s also worthy noting that the stiffness matrix seemed to have non-zero values just the first few rows, but the rest seemd to be absolutely zero. Here’s the entire thing.
clear
clc
cd(‘G:My DriveStudiesM.Sc. Mechanical Engineering24 SMMatlabBeam_Atli’)
stlFile = ‘beam_Atli2.stl’;
figure
trisurf(stlread(stlFile))
axis equal
E = 10e9; % Young’s modulus in Pa
nu = 0.26; % Poisson’s ratio (nondimensional)
rho = 7800; % Mass density in kg/m^3
origins = [60 0.5 -5 % Frame 1: Cylinder connection point
30 0.5 0 % Frame 2: Bucket connection point
0 0.5 -5]; % Frame 3: Fulcrum point
numFrames = size(origins,1);
feModel = createpde(‘structural’,’modal-solid’);
importGeometry(feModel,stlFile);
structuralProperties(feModel, …
‘YoungsModulus’,E, …
‘PoissonsRatio’,nu, …
‘MassDensity’,rho);
generateMesh(feModel, …
‘GeometricOrder’,’quadratic’, …
‘Hmax’,2.5, …
‘Hmin’,0.25);
figure
pdegplot(feModel,’FaceLabels’,’on’,’FaceAlpha’,0.5)
faceIDs = [1,2,3]; % List in the same order as the interface frame origins
figure
pdemesh(feModel,’FaceAlpha’,0.5)
hold on
colors = [‘rgb’ repmat(‘k’,1,numFrames-3)];
assert(numel(faceIDs) == numFrames);
for k = 1:numFrames
nodeIdxs = findNodes(feModel.Mesh,’region’,’Face’,faceIDs(k));
scatter3( …
feModel.Mesh.Nodes(1,nodeIdxs), …
feModel.Mesh.Nodes(2,nodeIdxs), …
feModel.Mesh.Nodes(3,nodeIdxs), …
‘ok’,’MarkerFaceColor’,colors(k))
scatter3( …
origins(k,1), …
origins(k,2), …
origins(k,3), …
80,colors(k),’filled’,’s’)
end
hold off
structuralBC(feModel, ‘Face’, faceIDs(1), ‘Constraint’, ‘fixed’);
for k = 1:numFrames
structuralBC(feModel, …
‘Face’,faceIDs(k), …
‘Constraint’,’multipoint’, …
‘Reference’,origins(k,:));
end
rom = reduce(feModel,’FrequencyRange’,[0 1e4]);
arm.P = rom.ReferenceLocations’; % Interface frame locations (n x 3 matrix)
arm.K = rom.K; % Reduced stiffness matrix
arm.M = rom.M; % Reduced mass matrix
dampingRatio = 0.05;
arm.C = computeModalDampingMatrix(dampingRatio,rom.K,rom.M);
frmPerm = zeros(numFrames,1); % Frame permutation vector
dofPerm = 1:size(arm.K,1); % DOF permutation vector
assert(size(arm.P,1) == numFrames);
for i = 1:numFrames
for j = 1:numFrames
if isequal(arm.P(j,:),origins(i,:))
frmPerm(i) = j;
dofPerm(6*(i-1)+(1:6)) = 6*(j-1)+(1:6);
continue;
end
end
end
assert(numel(frmPerm) == numFrames);
assert(numel(dofPerm) == size(arm.K,1));
arm.P = arm.P(frmPerm,:);
arm.K = arm.K(dofPerm,:);
arm.K = arm.K(:,dofPerm);
arm.M = arm.M(dofPerm,:);
arm.M = arm.M(:,dofPerm);
arm.C = arm.C(dofPerm,:);
arm.C = arm.C(:,dofPerm);
function C = computeModalDampingMatrix(dampingRatio,K,M)
% To avoid numerical issues (such as complex eigenvalues with very small
% imaginary parts), make the matrices exactly symmetric.
K = (K+K’)/2; % Stiffness matrix
M = (M+M’)/2; % Mass matrix
% Compute the eigen-decomposition associated with the mass and stiffness
% matrices, sorting the eigenvalues in ascending order and permuting
% the corresponding eigenvectors.
[V,D] = eig(K,M);
[d,sortIdxs] = sort(diag(D));
V = V(:,sortIdxs);
% Due to small numerical errors, the six eigenvalues associated with the
% rigid-body modes may not be exactly zero. To avoid numerical issues,
% check that the first six eigenvalues are close enough to zero. Then
% replace them with exact 0 values.
assert(all(abs(d(1:6))/abs(d(7)) < 1e-9),’Error due to "zero" eigenvalues.’);
d(1:6) = 0;
% Vectors of generalized masses and natural frequencies
MV = M*V;
generalizedMasses = diag(V’*MV);
naturalFrequencies = sqrt(d);
% Compute the modal damping matrix associated with K and M
C = MV * diag(2*dampingRatio*naturalFrequencies./generalizedMasses) * MV’;
% Print the eigenmodes
disp(‘Eigenmodes (eigenvectors) V:’);
disp(V);
disp(‘Eigenvalues D:’);
disp(d);
endI’m following this documentation: Model an Excavator Dipper Arm as a Flexible Body – MATLAB & Simulink (mathworks.com).
I copied the entire code and it works with the Dipper.STL file, but for every other STL file i import, I get the ‘zero eigenvalue" error mentioned in the title. The simplest STL file I used was a cube modeled in spaceclaim (10 x 10 x 60 mm). I’ve varied the mesh fineness to be on a spectrum similar to the dipper arm.
I’ve looked at all settings and not figured out what might be wrong. I also added a fixed support for the first frame; didn’t make a difference. It’s also worthy noting that the stiffness matrix seemed to have non-zero values just the first few rows, but the rest seemd to be absolutely zero. Here’s the entire thing.
clear
clc
cd(‘G:My DriveStudiesM.Sc. Mechanical Engineering24 SMMatlabBeam_Atli’)
stlFile = ‘beam_Atli2.stl’;
figure
trisurf(stlread(stlFile))
axis equal
E = 10e9; % Young’s modulus in Pa
nu = 0.26; % Poisson’s ratio (nondimensional)
rho = 7800; % Mass density in kg/m^3
origins = [60 0.5 -5 % Frame 1: Cylinder connection point
30 0.5 0 % Frame 2: Bucket connection point
0 0.5 -5]; % Frame 3: Fulcrum point
numFrames = size(origins,1);
feModel = createpde(‘structural’,’modal-solid’);
importGeometry(feModel,stlFile);
structuralProperties(feModel, …
‘YoungsModulus’,E, …
‘PoissonsRatio’,nu, …
‘MassDensity’,rho);
generateMesh(feModel, …
‘GeometricOrder’,’quadratic’, …
‘Hmax’,2.5, …
‘Hmin’,0.25);
figure
pdegplot(feModel,’FaceLabels’,’on’,’FaceAlpha’,0.5)
faceIDs = [1,2,3]; % List in the same order as the interface frame origins
figure
pdemesh(feModel,’FaceAlpha’,0.5)
hold on
colors = [‘rgb’ repmat(‘k’,1,numFrames-3)];
assert(numel(faceIDs) == numFrames);
for k = 1:numFrames
nodeIdxs = findNodes(feModel.Mesh,’region’,’Face’,faceIDs(k));
scatter3( …
feModel.Mesh.Nodes(1,nodeIdxs), …
feModel.Mesh.Nodes(2,nodeIdxs), …
feModel.Mesh.Nodes(3,nodeIdxs), …
‘ok’,’MarkerFaceColor’,colors(k))
scatter3( …
origins(k,1), …
origins(k,2), …
origins(k,3), …
80,colors(k),’filled’,’s’)
end
hold off
structuralBC(feModel, ‘Face’, faceIDs(1), ‘Constraint’, ‘fixed’);
for k = 1:numFrames
structuralBC(feModel, …
‘Face’,faceIDs(k), …
‘Constraint’,’multipoint’, …
‘Reference’,origins(k,:));
end
rom = reduce(feModel,’FrequencyRange’,[0 1e4]);
arm.P = rom.ReferenceLocations’; % Interface frame locations (n x 3 matrix)
arm.K = rom.K; % Reduced stiffness matrix
arm.M = rom.M; % Reduced mass matrix
dampingRatio = 0.05;
arm.C = computeModalDampingMatrix(dampingRatio,rom.K,rom.M);
frmPerm = zeros(numFrames,1); % Frame permutation vector
dofPerm = 1:size(arm.K,1); % DOF permutation vector
assert(size(arm.P,1) == numFrames);
for i = 1:numFrames
for j = 1:numFrames
if isequal(arm.P(j,:),origins(i,:))
frmPerm(i) = j;
dofPerm(6*(i-1)+(1:6)) = 6*(j-1)+(1:6);
continue;
end
end
end
assert(numel(frmPerm) == numFrames);
assert(numel(dofPerm) == size(arm.K,1));
arm.P = arm.P(frmPerm,:);
arm.K = arm.K(dofPerm,:);
arm.K = arm.K(:,dofPerm);
arm.M = arm.M(dofPerm,:);
arm.M = arm.M(:,dofPerm);
arm.C = arm.C(dofPerm,:);
arm.C = arm.C(:,dofPerm);
function C = computeModalDampingMatrix(dampingRatio,K,M)
% To avoid numerical issues (such as complex eigenvalues with very small
% imaginary parts), make the matrices exactly symmetric.
K = (K+K’)/2; % Stiffness matrix
M = (M+M’)/2; % Mass matrix
% Compute the eigen-decomposition associated with the mass and stiffness
% matrices, sorting the eigenvalues in ascending order and permuting
% the corresponding eigenvectors.
[V,D] = eig(K,M);
[d,sortIdxs] = sort(diag(D));
V = V(:,sortIdxs);
% Due to small numerical errors, the six eigenvalues associated with the
% rigid-body modes may not be exactly zero. To avoid numerical issues,
% check that the first six eigenvalues are close enough to zero. Then
% replace them with exact 0 values.
assert(all(abs(d(1:6))/abs(d(7)) < 1e-9),’Error due to "zero" eigenvalues.’);
d(1:6) = 0;
% Vectors of generalized masses and natural frequencies
MV = M*V;
generalizedMasses = diag(V’*MV);
naturalFrequencies = sqrt(d);
% Compute the modal damping matrix associated with K and M
C = MV * diag(2*dampingRatio*naturalFrequencies./generalizedMasses) * MV’;
% Print the eigenmodes
disp(‘Eigenmodes (eigenvectors) V:’);
disp(V);
disp(‘Eigenvalues D:’);
disp(d);
end I’m following this documentation: Model an Excavator Dipper Arm as a Flexible Body – MATLAB & Simulink (mathworks.com).
I copied the entire code and it works with the Dipper.STL file, but for every other STL file i import, I get the ‘zero eigenvalue" error mentioned in the title. The simplest STL file I used was a cube modeled in spaceclaim (10 x 10 x 60 mm). I’ve varied the mesh fineness to be on a spectrum similar to the dipper arm.
I’ve looked at all settings and not figured out what might be wrong. I also added a fixed support for the first frame; didn’t make a difference. It’s also worthy noting that the stiffness matrix seemed to have non-zero values just the first few rows, but the rest seemd to be absolutely zero. Here’s the entire thing.
clear
clc
cd(‘G:My DriveStudiesM.Sc. Mechanical Engineering24 SMMatlabBeam_Atli’)
stlFile = ‘beam_Atli2.stl’;
figure
trisurf(stlread(stlFile))
axis equal
E = 10e9; % Young’s modulus in Pa
nu = 0.26; % Poisson’s ratio (nondimensional)
rho = 7800; % Mass density in kg/m^3
origins = [60 0.5 -5 % Frame 1: Cylinder connection point
30 0.5 0 % Frame 2: Bucket connection point
0 0.5 -5]; % Frame 3: Fulcrum point
numFrames = size(origins,1);
feModel = createpde(‘structural’,’modal-solid’);
importGeometry(feModel,stlFile);
structuralProperties(feModel, …
‘YoungsModulus’,E, …
‘PoissonsRatio’,nu, …
‘MassDensity’,rho);
generateMesh(feModel, …
‘GeometricOrder’,’quadratic’, …
‘Hmax’,2.5, …
‘Hmin’,0.25);
figure
pdegplot(feModel,’FaceLabels’,’on’,’FaceAlpha’,0.5)
faceIDs = [1,2,3]; % List in the same order as the interface frame origins
figure
pdemesh(feModel,’FaceAlpha’,0.5)
hold on
colors = [‘rgb’ repmat(‘k’,1,numFrames-3)];
assert(numel(faceIDs) == numFrames);
for k = 1:numFrames
nodeIdxs = findNodes(feModel.Mesh,’region’,’Face’,faceIDs(k));
scatter3( …
feModel.Mesh.Nodes(1,nodeIdxs), …
feModel.Mesh.Nodes(2,nodeIdxs), …
feModel.Mesh.Nodes(3,nodeIdxs), …
‘ok’,’MarkerFaceColor’,colors(k))
scatter3( …
origins(k,1), …
origins(k,2), …
origins(k,3), …
80,colors(k),’filled’,’s’)
end
hold off
structuralBC(feModel, ‘Face’, faceIDs(1), ‘Constraint’, ‘fixed’);
for k = 1:numFrames
structuralBC(feModel, …
‘Face’,faceIDs(k), …
‘Constraint’,’multipoint’, …
‘Reference’,origins(k,:));
end
rom = reduce(feModel,’FrequencyRange’,[0 1e4]);
arm.P = rom.ReferenceLocations’; % Interface frame locations (n x 3 matrix)
arm.K = rom.K; % Reduced stiffness matrix
arm.M = rom.M; % Reduced mass matrix
dampingRatio = 0.05;
arm.C = computeModalDampingMatrix(dampingRatio,rom.K,rom.M);
frmPerm = zeros(numFrames,1); % Frame permutation vector
dofPerm = 1:size(arm.K,1); % DOF permutation vector
assert(size(arm.P,1) == numFrames);
for i = 1:numFrames
for j = 1:numFrames
if isequal(arm.P(j,:),origins(i,:))
frmPerm(i) = j;
dofPerm(6*(i-1)+(1:6)) = 6*(j-1)+(1:6);
continue;
end
end
end
assert(numel(frmPerm) == numFrames);
assert(numel(dofPerm) == size(arm.K,1));
arm.P = arm.P(frmPerm,:);
arm.K = arm.K(dofPerm,:);
arm.K = arm.K(:,dofPerm);
arm.M = arm.M(dofPerm,:);
arm.M = arm.M(:,dofPerm);
arm.C = arm.C(dofPerm,:);
arm.C = arm.C(:,dofPerm);
function C = computeModalDampingMatrix(dampingRatio,K,M)
% To avoid numerical issues (such as complex eigenvalues with very small
% imaginary parts), make the matrices exactly symmetric.
K = (K+K’)/2; % Stiffness matrix
M = (M+M’)/2; % Mass matrix
% Compute the eigen-decomposition associated with the mass and stiffness
% matrices, sorting the eigenvalues in ascending order and permuting
% the corresponding eigenvectors.
[V,D] = eig(K,M);
[d,sortIdxs] = sort(diag(D));
V = V(:,sortIdxs);
% Due to small numerical errors, the six eigenvalues associated with the
% rigid-body modes may not be exactly zero. To avoid numerical issues,
% check that the first six eigenvalues are close enough to zero. Then
% replace them with exact 0 values.
assert(all(abs(d(1:6))/abs(d(7)) < 1e-9),’Error due to "zero" eigenvalues.’);
d(1:6) = 0;
% Vectors of generalized masses and natural frequencies
MV = M*V;
generalizedMasses = diag(V’*MV);
naturalFrequencies = sqrt(d);
% Compute the modal damping matrix associated with K and M
C = MV * diag(2*dampingRatio*naturalFrequencies./generalizedMasses) * MV’;
% Print the eigenmodes
disp(‘Eigenmodes (eigenvectors) V:’);
disp(V);
disp(‘Eigenvalues D:’);
disp(d);
end matlab, simulink, stl, error, eigenvalue, flexible MATLAB Answers — New Questions
