坚持实干、坚持一线、坚持积累、坚持思考，坚持创新。
接着博文《[动力学] 模态分析基础 (Basics of modal analysis)》继续介绍桁架单元用于模态分析的方法。该部分内容也是 书本 《有限单元法：编程与软件应用》模态分析章节的部分内容节选。
- 10.2 2D桁架结构模态分析
模态分析需要的参数主要有结构的刚度矩阵和质量矩阵，其中的刚度矩阵与静力分析中的刚度矩阵相同，这里仅给出进行模态分析所需的质量矩阵等其他条件。
- 10.2.1 桁架单元质量矩阵
集中质量矩阵（局部坐标）
设杆材料密度为 ，单元长度为 ，截面积为 ，每个节点分担单元1/2质量，则单元质量矩阵
(10.2‑1)
整体坐标下的质量矩阵转换
整体坐标下单元一致质量矩阵为
(10.2‑2)
【推导】：质点动能与坐标系无关，则有
(10.2‑3)
即
(10.2‑4)
从而得到式(10.2‑2)。
- 10.2.2 算例：2D桁架模态分析
算例采用与2D桁架静力分析中相同的结构，其中材料密度为7.85e-9t/mm3，采用集中质量矩阵，将单元质量集中于两端节点。由于算例模态分析的Matlab代码与前面章节中静力分析的Matlab代码大体相同，因此这里仅给进行模态分析所需的新增的 MATLAB 代码。
% Truss 2D modal analysis
% Author ： JiDong Cui（崔济东），Xuelong Shen（沈雪龙）
% Website : www.jdcui.com
% 20151108
……
此前代码与静力分析代码相同，不再赘述。
% Properties definition
……
p=7.85E-9 % Density
……
定义单元属性时，增加单元材料密度ρ=7.85e-9，用于计算单元质量。
% 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;
-S,C, 0,0;
0,0, C,S;
0,0,-S,C
];
eleK = EA/L*[ C*C C*S -C*C -C*S;
C*S S*S -C*S -S*S;
-C*C -C*S C*C C*S;
-C*S -S*S C*S S*S];
eleM0=p*A*L/2*diag([1,1,1,1]);
eleM=T'*eleM0*T;
% Corresponding freedom of the element nodes
eleDof = [noindex(1)*2-1 noindex(1)*2 noindex(2)*2-1 noindex(2)*2];
% Integration of element stiffness matrix and mass matrix
stiffness( eleDof , eleDof) = stiffness( eleDof , eleDof) + eleK;
mass( eleDof , eleDof) = mass( eleDof , eleDof) + eleM;
end
……
本段代码的作用是集成结构整体刚度矩阵和整体质量矩阵。本结构共有numDOF个自由度，因此首先采用zeros()方法定义了一个numDOF×numDOF的零矩阵mass。单元的质量矩阵eleM0采用集中质量形式求得，并采用与刚度矩阵相同的“对号入座法”将单元质量矩阵组装到结构整体质量矩阵mass中。
% Get the eigenvalues and eigenvectors
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))';
f=sort(fo)
……
本段代码通过求解标准特征值问题，求得特征值和特征向量。Matlab内置函数eig([A])可直接求解特征值和特征向量，本例中 [A]=[K]-1[M]，求得的特征值放在向量{Do}中，相应的特征向量放在矩阵[Vo]中。其后的代码用于处理求解得到的数据，其中代码“fo=”用于根据特征值求得频率，代码“V=”用于将特征向量排序，代码“f=”用于将结构自振频率进行排序。图 10‑1所示是该桁架结构的前三阶振型。
（a）第一振型
（b）第二振型
（c）第三振型
图 10‑1 桁架结构振型
- 相关内容 ( 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)
欢迎关注 “结构之旅” 微信公众号