坚持实干、坚持一线、坚持积累、坚持思考，坚持创新。
接着博文《[动力学] 模态分析基础 (Basics of modal analysis)》继续介绍框架单元用于模态分析的方法。该部分内容也是 书本 《有限单元法：编程与软件应用》模态分析章节的部分内容节选。
10.3.1 梁单元质量矩阵
集中质量矩阵（局部坐标）
设杆材料密度为 ，单元长度为 ，截面积为 ，每个节点分担单元1/2的平动质量，无转动惯量，则单元质量矩阵
(10.3‑1)
10.3.2 算例：2D框架模态分析
算例采用与2D框架结构静力分析中相同的结构，采用欧拉梁单元，材料密度为2.5493e-9t/mm3，采用集中质量矩阵，将单元质量集中于两端节点，且不考虑集中质量后的节点转动质量。
由于算例模态分析的Matlab代码与前面章节中静力分析的Matlab代码大体相同，因此这里仅给进行模态分析所需的新增代码。
% Truss 2D modal analysis
% Author ： JiDong Cui（崔济东），Xuelong Shen（沈雪龙）
% Website : www.jdcui.com
% 20151114
……
此前代码与静力分析代码相同，不再赘述。
% Properties definition
……
p=2.5493 %Density
……
定义单元属性时，增加单元密度ρ=2.5493，即重度25kN/m3，用于计算单元质量。
% Global stiffness matrix and mass matrix
stiffness = zeros(numDOF);
mass=zeros(numDOF);
for i=1:numEle;
% Index of the element nodes
noindex = EleNode(i,:);
% Element length
deltax = xx( noindex(2)) - xx( noindex(1));
deltay = yy( noindex(2)) - yy( noindex(1));
L = sqrt( deltax*deltax + deltay*deltay );
% Element stiffness matrix
C = deltax / L;
S = deltay / L;
T = [ C S 0 0 0 0;
-S C 0 0 0 0;
0 0 1 0 0 0;
0 0 0 C S 0;
0 0 0 -S C 0;
0 0 0 0 0 1
];
% Section properties
E = EList(EleSect(i));
A = AList(EleSect(i));
I = IList(EleSect(i));
EAL = E*A/L;
EIL1 = E*I/L;
EIL2 = 6.0*E*I/(L*L);
EIL3 = 12.0*E*I/(L*L*L);
% Element stiffness matrix in local coordinate system
eleKo = [ EAL 0 0 -EAL 0 0;
0 EIL3 EIL2 0 -EIL3 EIL2;
0 EIL2 4*EIL1 0 -EIL2 2*EIL1;
-EAL 0 0 EAL 0 0;
0 -EIL3 -EIL2 0 EIL3 -EIL2;
0 EIL2 2*EIL1 0 -EIL2 4*EIL1
];
% Element mass matrix in local coordinate system
eleMo=p*A*L/2*diag([1,1,0,1,1,0]);
% Corresponding freedom of the element nodes
eleDof = [noindex(1)*3-2 noindex(1)*3-1 noindex(1)*3 noindex(2)*3-2 noindex(2)*3-1 noindex(2)*3];
% Integration of element stiffness matrix and mass matrix
stiffness( eleDof , eleDof) = stiffness( eleDof,eleDof) + T'*eleKo*T;
mass( eleDof , eleDof) = mass( eleDof , eleDof) + T'*eleMo*T;
end
本段代码的作用是集成结构整体刚度矩阵和整体质量矩阵。本结构共有numDOF个自由度，因此首先采用zeros()方法定义了一个numDOF×numDOF的零矩阵mass。单元的质量矩阵eleM0采用集中质量形式求得，并采用与刚度矩阵相同的“对号入座法”将单元质量矩阵组装到结构整体质量矩阵mass中。
% Get the eigenvalues and eigenvectors
activeDof = setdiff( [1:numDOF]',restrainedDof);
if det(mass(activeDof,activeDof))>0
A=mass(activeDof,activeDof)\stiffness(activeDof,activeDof);
[Vo,Do] = eig(A);
fo=diag(sqrt(Do)/(2*pi));
Vf=sortrows([Vo',fo],size(Vo,1)+1);
V=Vf(:,1:size(Vo,1))'; % Sort eigenvectors
elseif det(stiffness(activeDof,activeDof))>0
A=stiffness(activeDof,activeDof)\mass(activeDof,activeDof);
[Vo,Do] = eig(A);
fo=power(diag(sqrt(Do)),-1)/(2*pi);
Vf=sortrows([Vo',fo],size(Vo,1)+1);
V=Vf(:,1:size(Vo,1))'; % Sort eigenvectors
end
f=sort(fo)
本段代码通过求解标准特征值问题，求得特征值和特征向量。Matlab内置函数eig([A])可直接求解特征值和特征向量。由于本例不考虑节点的转动质量，因此 det(mass(activeDof,activeDof)) 等于零，所以本例做了一个判断“if…else if…”，最终采用的是 [A]=[K]-1[M] 。求得的特征值放在向量{Do}中，相应的特征向量放在矩阵[Vo]中。其后的代码用于处理求解得到的数据，其中代码“fo=”用于根据特征值求得频率，代码“V=”用于将特征向量排序，代码“f=”用于将结构自振频率进行排序。图 10‑7所示是该桁架结构的前三阶振型。
（a）第1振型
（b）第2振型
（c）第3振型
- 相关内容 ( Related Topics )
[00]. [Book][书] 有限单元法：编程与软件应用
[01]. 2D Truss by Matlab [Matlab 2D桁架有限元] (code by myself)
[02]. 3D Truss by MATLAB[MATLAB 3D 桁架静力分析] (code by myself)
[03]. [FEM][有限元][编程][Matlab][Code by myself] FEM Analysis: 2D Truss Element [有限元分析: 2D桁架单元]
[04]. [FEM][有限元][编程][Matlab][Code by myself] FEM Analysis: 2D Euler Beam Element [有限元分析: 2D欧拉梁单元]
[05]. [FEM][有限元][编程][Matlab][Code by myself] 2D剪切梁单元
[06]. [FEM][有限元][编程][Matlab][Code by myself] 2D Timoshenko梁单元
[07]. [FEM][有限元][编程][Matlab][Code by myself] 三角形常应变单元(CST)
[08]. [FEM][有限元][编程][Matlab][Code by myself] 平面4变形单元(Q4)
[09]. [FEM][有限元][编程][Matlab][Code by myself] 4节点四面体单元(TET4)
[10]. [FEM][有限元][编程][Matlab][Code by myself] 8节点六面体单元(C3D8)
[11]. [FEM][有限元][编程][Matlab][Code by myself] 平面8节点二次“完全积分”单元(CPS8)
[12]. [FEM][有限元][编程][Matlab][Code by myself] 平面6节点二次“完全积分”单元(CPS6)
[13]. [FEM][MATLAB][有限元] FEM Modal Analysis Programming with MATLAB (Truss Element) (桁架单元模态分析编程)
[14]. [FEM][MATLAB][有限元] FEM Modal Analysis Programming with MATLAB (Frame Elements) (框架单元模态分析编程)
[15]. [力学][有限元]Basics of Linear Buckling Analysis [线性曲屈分析基础]
[16]. [FEM][MATLAB][有限元][编程] FEM Buckling Analysis Programming with MATLAB (Truss Elements) (桁架单元曲屈分析编程)
[17]. [FEM][MATLAB][有限元][编程] FEM Buckling Analysis Programming with MATLAB (Frame Elements) (框架单元曲屈分析编程)
- 注释 ( Comments )
( 如果您发现有错误，欢迎批评指正。邮箱：jidong_cui@163.com . 如果您喜欢这篇博文，请在上面给我 点个赞 吧！ 🙂 🙂
( If you found any mistakes in the post, please let me know. Email : jidong_cui@163.com. If you like this posts, please give me a “thumbs up“ rating on the above button! )
- 微信公众号 ( Wechat Subscription)
欢迎关注 “结构之旅” 微信公众号
