基于MATLAB的振动模态分析

摘要
振动系统是研究机械振动的运动学和动力学，研究单自由系统的振动有着实际意义，因为工程上有许多问题通过简化，用单自由度系统的振动理论就能得到满意的结果。

模态是振动系统的一种固有振动特性，模态一般包含频率、振型、阻尼。

振动系统问题是个比较虚拟的问题，比较抽象的理论分析，对于问题的分析可以实体化建立数学模型，通过MATLAB可以转化成为图像。

单自由度频率、阻尼、振型的分析，我们可以建立数学模型，最后通过利用MATLAB编程实现数据图形；多自由度主要研究矩阵的迭代求解，我们在分析抽象的理论的同时根据MATLAB编程实现数据的迭代最后可以得到所要的数据，使我们的计算更加简便。

利用MATLAB编程并验证程序的正确性。

通过程序的运行，能快速获得多自由度振动系统的固有频率以及主振型，为设计人员提供了防止系统共振的理论依据，也为初步分析各构件的振动情况以及解耦分析系统响应奠定了基础。

关键词：振动系统；单自由度；MATLAB；多自由度
Abstract
Vibration system is to study the kinematics and dynamics of mechanical vibration, the vibration of a single free system has practical significance, because there are many engineering problems by simplifying, using the vibration theory of a single degree of freedom system can be satisfied with the results.
Vibration system problems is a relatively virtual problems, more abstract and theoretical analysis, problem analysis for a mathematical model can be materialized by MATLAB can be converted into images. Single degree of freedom frequency, damping, mode shape analysis, we can create mathematical models, the final program data through the use of MATLAB graphics; many degrees of freedom main matrix iterative solution, our analysis based on abstract theory, while MATLAB programming The last iteration of data can be the desired data, so our calculations easier
Using MATLAB programming and verify the correctness of the program.Through the process of operation, can quickly obtain multiple degrees of freedom vibration system and the main vibration mode natural frequency for the design to prevent resonance provide the theoretical basis for the preliminary analysis of the vibration of each component, and laid the decoupling of system response basis.
Key words:vibrating system; Single Degree of Freedom ;MATLAB; multiple degree of
freedom
辽宁工程技术大学毕业设计（论文）
1 绪论
1.1问题的提出
机械振动是一门既古老又年轻的科学，随着人类科学技术的不断进步振动理论得到不断的发展和完善。

机械振动在许多情况下是有害的，人们想方设法避免它：另一方面，人们利用机械振动原理制造了各种机械或仪表来为人类服务。

振动机械是20世纪后半期得到迅速发展的一类机械，它是利用振动原理来完成各种工艺过程的机械设备。

其中，Mathorks 公司推出的MATLAB以其强大的功能和易用性受到越来越多科技工作者的欢迎。

它把计算、可视化、程序设计融合到了一个交互的工作环境中，可以实现工程计算、算法研究、建模和仿真、数据分析及可视化、科学和工程绘图、应用程序开发（包括图形用户界面程序设计）等功能。

它在美国等发达国家的大学里已经成为一种必须掌握的基本编程语言，而在国外的研究设计单位和工业部门，更是早己成为研究和解决工程计算问题的一种标准软件。

在国内也有越来越多的科学技术工作者参加到学习和倡导这种语言的行列中来。

应用MATLAB软件对选矿用振动筛的振动特性进行研究，可以充分发挥计算机技术的优势，为选矿用振动筛振动特性研究探索新的途径。

在工程振动中，确定系统固有频率与主振型是非常重要的。

固有频率是决定系统振动特性的重要物理量，它既是防止系统共振的依据，又是多自由度系统解耦分析（模态分析）的前提，因此研究某系统振动时，首先要求出系统的固有频率。

主振型则为初步分析各构件的振动情况以及解耦分析奠定了基础。

对于多自由度振动系统，计算系统固有频率与主振型主要有2种方法[1]：(1)利用特征矩阵方程式与特征方程式求解；(2)矩阵迭代法求解。

2种方法各有各的特色。

对于低自由度的振动系统，方法一容易、快捷。

但是在实际工程中，大多数振动
系统都是自由度较多，用特征矩阵方程式与特征方程式求解系统固有频率与主振型这种传缆的计算方法虽然从原则上可行，但当自由度增加时，惯性、刚度阵的阶数增高，计算量也急剧加大，这显然很不方便。

但采用矩阵迭代法，即使是自由度很大的振动系统，计算量也只不过是多进行矩阵迭代而已，而且假设的初始矩阵愈接近实际状况，迭代的次数愈少，相应的计算量也愈少。

王超：基于MA TLAB的振动系统编程分析
1.2 国内外研究现状
1.2.1机械振动理论的发展状况及应用现状
振动理论是力学的一个重要组成部分[2]，人类对振动现象的认识有悠久的历史。

振动力学的物理基础在17世纪已经奠定，到了18世纪，振动力学已从物理学中独立出来。

最主要的成就为线性振动理论的形成，它是与数学中的常微分方程和偏微分方程同步发展的。

目前，振动及系统按运动微分方程的形式分为以下两种。

线性振动：描述其运动的方程为线性微分方程，相应的系统称为线性系统。

线性振动的一个重要特性是线性叠加原理成立。

非线性振动[3]:描述其运动的方程为非绒性微分方程，相应的系统称为非线性系统。

对于非线性振动叠加原理不再成立。

在实际的振动机械或振动系统中，严格的讲，都是非线性的。

但是，建立振动系统的非线性力学模型难度大，求解困难，有些问题甚至无解可求。

在实际的工程应用中，很多情况下在误差允许的范围之内用线性的方法解决复杂的近线性问题。

线性振动有确定的力学模型一一线性微分方程，可以求得准确的解，能够描述出振动系统的主要特征。

由于用线性振动的方法能够解决众多的工程实际问题，线性振动的理论一直倍受关注，并且在理论和实验方面已经得到很大的发展和成熟。

特别是多自由度系统的振动的理论，可以说既是振动力学的核心又是应用得最广泛的振动理论。

线性振动在当今不仅是作为基础科学的力学的一个重要组成部分，而且正走上向工程科学发展的道路，它在航空、机械、船舶、车辆、建筑、水利等工业技术部门中占有愈来愈重要的地位。

线性振动的应用可分为两个方面：一个方面是减少由于振动而造成的危害，目的在于减振甚至于避免有害的振动；另一个方面利用振动，如工业上常采用的振动筛选、振动沉桩、振动输送以及按振动理论设计的测量传感器、地震仪等等就是这方面的典型例子。

选矿用振动筛是振动筛选设各中的—种，线性振动理论在选矿用振动筛的设计制造及生产运行中有着广泛的应用，有关这方面的内容将在下一节中详细介绍。

线性振动的理论在发展过程中产生了一个重要分支，那就是模态分析理论。

在对选矿用振动筛进行分析时，需要通过实验来验证理论的正确性，振动实验则需要用到模态分析技术。

模态分析技术从20世纪60年代后期发展至今已趋成熟[4]。

它和有限元分析技术一起，已成为结构动力学中两大支柱。

模态分析是结构动力学中的～种“逆问题”分析方法，它与传统的“正问题”方法（主要是指有限元方法）不同，是建立在实验（或实测）的基础上，采用实验与理论相结合的方法来处理工程中的振动问题。

目前这一技术已发展成为解决工程中振动问题的重要手段，在机械、航空、航天、土木、建筑、造船、化工等
辽宁工程技术大学毕业设计（论文）
工程领域被广泛应用[5]。

近十年来，模态分析理论吸取了振动理论、信号处理、信号分析、数据处理、数理统计及自动控制理论中的有关“营养”，结合自身内容的发展，形成了一套独特的理论为模态分析及参数识别技术的发展奠定了理论基础。

模态分析的基础理论概念主要包括；机械阻抗、导纳、传递函数（或频响函数）、实模态、复模态等。

模态测试技术主要采用同时测量输入及输出的方法，对一个振动系统来说，可以表示成图l-1所示的框图
图1.1 模态分析框图
Fig. 1.1 Modal Analysis Diagram
通过测量激励和响应，进行模念分析可以确定系统。

自从FFT问世以来，目前广泛采用宽频带激振技术。

其中主要有脉冲、阶跃激励，快速正弦扫描等瞬态激励和纯随机、伪随机、周期随机、瞬态随机等激励方法。

此外，由于F弦慢扫描技术测试精度高，它仍不失为重要激励手段。

模态参数辨识的频域方法有：分量分析法、导纳圆辨识方法、正交多项式曲线拟合、非线性优化辨识方法等。

模态参数辨识的时域方法与模态参数辨识的频域方法不同，它无需将所测得的响应与激励的时间历程信号变换到频域中去，而是直接在时域中进行参数辨识。

它与频域法相比，两者所采取的分析路线不同，如图1.2所示。

图1.2 模态参数辨识分析路线框图
Fig.1.2 Modal parameter identification of line diagram
时域法比频域法发展较晚，但近几年来有长足的进展。

自70年代以来主要有：
Ibrahim时域法（简称LTD法）、最小二乘复指数法(LSCE法)、多参考点复指数法(PRCE 法)、特征系统实现算法(ERA)。

模态分析技术在动态载荷识别、模型修正与结构动力修改中有广泛的应用，结构动态特征灵敏度分析是非常重要的方法之一。

模态综合技术主要有组合系统法和模态综合法。

随着电子技术与计算机技术的迅速发展，模态分析已成为解决复杂结构振动问题的主要工具，并与计算机辅助设计(CAD)．计算机辅助实验(CAT)相结合，
王超：基于MA TLAB的振动系统编程分析
进入产品设计阶段，作为计算机辅助工程中的重要环节，有着广泛的应用[6]。

1.2.2 MATLAB软件的发展状况及应用现状
MATLAB软件概述：MATLAB的名称源自Matrix Laboratory，是一门计算语言口[7]。

在工程计算领域，计算机技术的应用正逐步将科技人员从繁重的计算工作中解放出来。

在科学计算和工程应用的过程中，一些技术人员尝试用Basic，Fortran以及C语言编制程序来减轻计算的工作量，但编制程序不仅要掌握所用语-的语法，还要对有关算法进行深入分析。

为了满足用户对工程数学计算的要求，MATLAB的功能、特点、应用范围：MATLAB 越来越广地被人们应用是源于它在求解方程、数值计算、程序编写上的优点，而它的这些优点是由它的功能和特点决定的。

MATLAB的主要功能:(1)数值计算功能，一条MATLAB 语句相当于几十条C语言或Fortran语言的语句。

(2)符号计算功能，利用MATLAB的符号计算功能可以清晰地获得解的表达式，对于避免出错和提高程序的可读性均有很大的帮助。

(3)数据分析和可视化功能，在科学计算和研究工作中，技术人员经常会遇到大量的原始数据，而对数据的分析往往难于入手。

MATLAB能将这些数据以图形的方式显示出来，不仅使数据间的关系清晰明了，而且对于揭示其内在本质往往有着非常重要的作用。

MATLAB提供了良好的用户界面，许多函数本身会自动绘制出图形，而且会自动选取坐标刻度，绘制出直角坐标、极坐标、对数坐标下的二维和三维图形，以及条形图、直方图、等高线图、饼形图、离散数据图和瀑布图等专用图形。

(4)文字处理功能。

MATLAB的主要特点：(1)功能强大，MATLAB不但在数值计算和符号计算方面具有强大的功能，而且在计算结果的分析和数据可视化方面也有着其它类似软件难以匹敌的优势[9]。

Notebook，Simulink功能以及各种专业工具箱将MATLAB的应用扩展到非常广的领域。

(2)界面友好、编程效率高，MATLAB的指令表达方式与标准教科书的数学表达式非常相近，用户不需要有较高的计算机编程基础，只要按照计算要求输入表达式，
MATLAB将为用户计算出结果。

同时使用MATLAB语言设计的程序，其编译和执行速度都超过了传统c和Fortran语言设计的程序，在工程计算方面的编程效率也高于其它编程语言。

(3)扩展性强，MATLAB的最重要特点之一就是它的可扩展性。

这个特点使得用户能够自由地开发自己的应用程序。

这些年来，许多使用MATLAB的科学家、工程师和技术人员已经开发出相当多的不同领域的应用程序。

MATLAB的应用范围：MATLAB由主包和各种工具箱组成。

主包是MATLAB的核心，工具箱是扩展的有专门功能的函数。

例如，控制系统工具箱应用于连续和离散系统设计、
辽宁工程技术大学毕业设计（论文）
频域和时域响应等控制领域；信号处理工具箱应用于自适应去噪和压缩、谱分析和估计等信号处理领域；通信工具箱应用于信号编码、调制解调等通信领域。

应用MATLAB
的各种工具箱可以在很大程度上减小用户编程时的复杂度，因此MATLAB在很广的领域内得到了应用，其典型应用有；自动控制、图像信号处理、生物医学工程、语音处理、雷达工程、信号分析、振动理论、时序分析与建模、化学、统计学、经济学等
1.3 MATLAB语言的优点
MATLAB作为一个以矩阵和数组为核心计算的软件，对矩阵迭代法中的矩阵迭代计算尤其适合[10]。

就所查的资料看，以前的学者和研究人员迭代求解系统固有频率与主振型时，大部分都是用Visiual Basic或Fortran语言来编写程序[11]。

限于Visiual Basic或Fortran本身语句以及语法的局限性，用这种高级语言编写的程序涉及到选择合适的算法和编写冗长的语言代码以及键入和调试等一系列问题。

即使有现成的标准予程序可供调用，要在一些较复杂的、科研问题中编写一个完整的程序仍然是一个复杂的、技巧性很强的工作。

因此，用高级语言编写的程序一般代码段较长，需要调用的子程序较多，整个程序的通读性较差。

相反，MATLAB则有简洁、可读性强等优点。

1.4本文研究的内容
振动机械在国民经济中占有重要的位置，振动筛是振动机械中的重要一员。

一直以来有许多人对振动筛进行设计和研究，但是，振动筛的动态设计和计算机辅助设计近年来刚刚起步。

振动特性是振动筛非常重要的有别于非振动机械的一个本质特点，却往往被设计者和制造者简单化。

客观的说，一般的振动都是非线性的，但在许多情况下可以近似看作线性来处理。

线性振动理论不论从基础理论还是实验技术方面近年来都有很大的发展，特别是应用现代化振动测试仪器测量振动信号以及应用计算机软件来分析处理振动信号，为从事振动研究的科技人员带来了极大的方便。

把振动的理论应用到工程实际中去，切实解决工程中遇到的实际的振动问题是研究振动理论的根本目的。

需要对该力学模型进行深入的分析（借助MATLAB软件进行仿真分析）。

本文主要利用MATLAB对振动系统进行模拟分析对于虚拟抽象的理论图像化，处理单自由度振动的3个阻尼和强迫单自由度阻尼振动，多自由度系统振动矩阵迭代求解。

王超：基于MA TLAB 的振动系统编程分析
2单自由度系统的振动
2.1 单自由度振动系统数学模型的建立
[12] 建立和分析有粘性阻尼时的自由度振动微分方程。

以静平衡位置为原点建立如图坐
标，由牛顿定律得运动方程为[13]:
0=++kx x c x
m (2-1) 令
m
k m c n n ==2,2ω 其中n 称为衰减系数，单位为s 1；n ω是相应的无阻尼时的固有频率，式(2-1)可以写为：
022=++x x n x n ω (2-2)
如果进一步令
n n
ως= (2-3)
其中无量纲的ς称为相对阻尼系数，则式(2-2)可写为：
022=++x x x n n ωςω
(2-4) 为了求解，令
x x c
辽宁工程技术大学毕业设计（论文）
st e x = (2-6)
代入(2-4)后得到特征方程：
0222=++n n s s ωςω (2-7)
他的两个特征根为： 122,1-±-=ςωςωn n s (2-8)
根据相对阻尼系数ς的不同大小，可以将阻尼分为三种状态：1>ς时为过阻尼，1=ς时为临界阻尼，10<<ς时为欠阻尼。

1）过阻尼状态
1>ς，1s 与2s 是两个不等的负实根，令
12*-=ςωωn （2-9）
初始条件
00)0(,)0(x x
x x == （2-10） 系统初始条件响应为
)()(**00*0t sh x x t ch x e t x n t n ωωςωωςω++
=- （2-11）
临界阻尼状态 n s ως-==,1是二重根，方程（2-4）的通解为系统对式（2-10）的初始条件的响应为
])([)(000t x x
x e t x n t n ωω++=- （2-12） 欠阻尼状态
1<ς，其中
21ςωω-=n d （2-13）
初始条件响应
)sin cos ()(000t x x
t x e t x d d n d t n ωωςωωςω++=- （2-14）
2.2 参数设定与求解
阻尼比ς分别取；应用Matlab 对式(2-11)和(2-12),（2-14）求解。

程序如下：
王超：基于MA TLAB的振动系统编程分析
clear,format compact;
a=0.5;t=0:0.1:18;;w0=1;
k=1;x0=1;
wd=w0.*sqrt(1-a*a);x1=wd
y=exp(-a*w0.*t).*(x0.*cos(wd.*t)+((x1+a*wd*x0)./wd)*sin(wd.*t))
figure(1),plot(t,y,'r');hold on
a=1.0;t=0:0.1:18;
w0=1;wd=1;x1=wd;
y=exp(-wd.*t).*(x0+(x1+wd*x0).*t);
figure(1),plot(t,y,'d');hold on
a=2.0;t=0:0.1:18;w0=1;wd=w0*sqrt(a*a-1);
y=exp(-a*w0.*t).*(x0.*cosh(wd.*t)+(x1+a*w0*x0)/w0.*sinh(t));
figure(1),plot(t,y,'v');hold on
结论：
图2-2为Matlab计算后给出的响应曲线，从中可以得到一些重要的结论[14]：
在1
<ς的情况下，阶跃信号输入时，输出信号为衰减振荡，其振荡角频率（阻尼0<
ω，幅值按指数衰减越大，阻尼越大，衰减越快。

振荡角频率）为
d
ς时，振荡系统等同于两个一阶系统串联。

此时虽然不产生振荡，但也需要经过较1
>
长时间才能达到稳态。

在一定的ς之下，欠阻尼系统能够更快地达到稳态值；而过阻尼系统反应迟饨，动作缓慢，所以系统通常设计成欠阻尼系统，ς取值为2
t(s)
x (t )
图2-2
算例绘制无阻尼单自由度系统的固有频率和周期随静变形的变化曲线。

固有频率n ω和周期n τ
st
n g
δω=
，g
st n δπ
τ2=
取2/81.9s m g =。

可以利用下列MATLAB 程序画出st δ在0~0.5范围内n ω和n τ的变换曲线：
%Ex2_17.m g=9.81;
for i=1:101 t(i)=0.01+(0.5-0.01)*(i-1)/100;w(i)=(g/t(i))^0.5; tao(i)=2*pi*(t(i)/g)^0.5; end
plot(t,w);gtext('w_n'); hold on;plot(t,tao);gtext('T_n');
xlabel('delta_s_t'); title('Example2.1');
0.05
0.1
0.15
0.2
0.250.3
0.35
0.4
0.45
0.5
05
10
15
20
25
30
35
delta st
Example2.1
2.3单自由度系统的强迫振动
[15]
简谐激励是激励形式中最简单的一种，虽然它在实际中存在的场合比较少但掌握系统对于简谐激励的响应的规律，是理解系统对周期激励或更一般形式激励的响应基础。

图所示的弹簧质量系统中，质量块上作用有简谐激振力
t P t P ωsin )(0= (2-15)
其中0P 为激振力幅，ω为激振频率。

以静平衡位置为坐标原点建立图示的坐标系。

从图的受力分析，得到运动微分方程为：
t p kx x c x
m ωsin 0=++ (2-16) 由常微分方程理论知道，方程(3.2)的通解x 由相应的齐次方程的通解h x 和非齐次方程的任意特解p x 两部分组成，即
)()()(t x t x t x p h += (2-17)
当欠阻尼时，式中)(t x h 为有阻尼自由振动，它的特点是振动频率为阻尼固有频率，振幅按指数规律衰减，称为瞬态振动或瞬态响应；)(t x p 是一种持续的等幅振动，它是由于简谐激励振力的持续作用而产生的，称之为稳态强迫振动或稳态振动，在间隔充分长时间考虑的振动就是这种稳态振动，而在刚受到外界激励时，系统的响应则是上述两种振动之和。

可见，系统受简谐激励后的响应可以分为两个阶段，一开始的过程称为过渡阶段，经充分长时间后，瞬态响应消失这时进入过渡阶段，经充分长时间后，瞬态响应消失，这是进入稳态阶段。

将方程（2-15）的两端同除以质量m ，并且令
2
2n m
c ςω= （2-18） 其中ς为相对阻尼系数，n ω为相应的无阻尼系统的固有频率，则方程（2-15）成为
x
m
t m
p x x x
n n ωωςωsin 20
2=++ （2-19） 上述方程特解可以通过)sin(ϕω-=t B x 或者t B t A x ωωsin cos +=来求得，这里介绍用复数方法求式（2-19）的特解。

先将式（2-19）写为下列的复数形式
t
i n n e m
p x x x
ωωςω022=++ （2-20） 其中x 是复数设复数形式的特解为
t i Be x ω= （2-21）
其中B 称为复振幅，其意义是包含有相位的振幅。

将式（2-21）代入（2-20），解得
ω
ςωωωn n i m P B 21
2
20+-=
（2-22） 记λ为频率比，它定义为
n
ωω
λ=
（2-23） 则式（2-22）可以写成
ϕςλλςλλi i Be e k
p i k P B --=+-=
+-=
2
22020)2()1(1211
（2-24）
式中
2
2
20
)
2()1(1ςλλ+-=
k
p B （2-25）
2
1
12λςλ
ϕ-=-tg （2-26） 将式（2-24）代入（2-21），得到复数形式的特解为
)(ϕω-=t i Be x （2-27）
比较方程（2-17）与（2-18），可知（2-19）中的位移x 是（2-20）中复数x 的虚部，因此（2-25）的虚部就是方程（2-12）的特解，即有
)sin(ϕω-=t B x （2-28）
其中B 为振幅，ϕ为相位差。

由式（2-26）、2-23）及（2-24）得出稳态强迫振动有如下的基本特点：
1)线性系统对简谐激励的稳态响应是频率等同于激励频率而相位滞后与激振力的简谐振动；
2)稳态响应的振幅及相位差只取决于系统本身的物理性质和激振的频率及力幅，而与系统本身进入运动的方式无关。

无阻尼系统对简谐激励的稳态响应可以从式（2-26）得出。

当n ωω<时，得到1<λ，0=ϕ，这时
t k P x ωλsin 11
2
0-=
（2-29）
当n ωω>时，得到1>λ，πϕ=，这时
)sin(1
1
2
0πωλ--=
t k P x （2-23） 式（2-21）也可以写成（2-22）的形式，这时相位差反映在振幅2
011
λ-k P 的符号中。

上述结果也可以由直接设t B x ωsin =并代入下列方程而得到：
t P kx x
m ωsin 0=+ （2-24） 为了具体讨论影响稳定响应的振幅和相位差的各种因素，记
k
P B 0
0=
（2-25） 0B 实际是质量块在激振力幅静作用下的最大位移。

再引入无量纲的振幅放大因子
β，它定义为
2220)
2()1(1
ςλλβ+-=
=
B B （2-26） 由式（2-26）和（2-19）可以分别画出以相对阻尼系数ς为参数的曲线——λβ-曲线与λϕ-曲线，前者称为幅频响应曲线，后者称为相频响应曲线如图所示
程序如下
for kesai=[0.05,0.10,0.15,0.25,0.375,0.50,1.0] lamda=0:0.01:5.0;
beta=1./(sqrt((1-lamda.^2).^2+(2*kesai*lamda).^2)); plot(lamda,beta)
hold on end
axis([0 5 0 3]);
频率比
振幅放大因子
偏心质量引起的强迫振动振幅放大因子
2222
)
2()1(ςλλλ+-=me MB (2-27)
程序如下：
for kesai=[0.05,0.10,0.15,0.25,0.375,0.50] lamda=0:0.01:5.0;
beta=lamda./(sqrt((1-lamda.^2).^2+(2*kesai*lamda).^2)); plot(lamda,beta) hold on end
axis([0 5 0 3]);
频率比
M B /m e
支撑运动引起的强迫振动振幅放大因子
2
222
)2()1()2(1ςλλςλβ+-+=
=a B (2-28) 程序如下：
for kesai=[0.05,0.10,0.15,0.25,0.375,0.50,1.0] lamda=0:0.01:5.0;
beta=sqrt((1+(2*kesai*lamda).^2)./((1-lamda.^2).^2+(2*kesai*lamda).^2)); plot(lamda,beta) hold on end
axis([0 5 0 3]);
00.51 1.52
2.53
3.54
4.55
0.5
1
1.5
2
2.5
3
频率比
振幅放大因子
算例利用MATLAB ，绘制弹簧-质量系统在简谐力作用下的响应曲线。

已知数据
如下：s m x
m x tN t F m N k kg m /1.0,1.0,30cos 100)(,/2000,500===== 。

系统全解形式如下：
t f t f x t x
t x n n n n n
ωω
ωωωωωωcos cos )(sin )(2
20
22000-+--
+=
式中，
s rad s rad m
k
m F f n /30,/20,205
100
00======
ωω 利用MATLAB 绘制解曲线上式的程序如下： %Ex3_11.m F0=100; wn=20; m=5; w=30;
x0=0.1; x0_dot=0.1; f_0=F0/m; for i=1:101
t(i)=2*(i-1)/100;
x(i)=x0_dot*sin(wn*t(i))/wn+(x0-f_0/(wn^2-w^2))*cos(wn*t(i)) +f_0/(wn^2-w^2)*cos(w*t(i)); end plot(t,x); xlabel('t'); ylabel('x(t)'); title('Ex3.11')
t
x (t )
Ex3.11
2.4本章小结
基于MATLAB 对单自由度自由振动绘制振动图像，进行粘性阻尼，强迫振动振幅放大因子绘图进行数据分析，使振动数据更加明显。

3 基于MATLAB 的多自由度系统编程分析
3.1 多自由度系统
[16]
多自由度振动系统的数学模型：
[]{}[]{}[]{}{}M x C x K x f ++= （3-1）
其中[]M 、[]C 、[]K 、{}f 和{}x 分别为质量矩阵、阻尼矩阵、刚度矩阵、力向量和响应向量。

把这个时域矩阵方程变换到拉氏域（变数为p ），并假定初始位移和初始速度为零，则得：
[][][]{}{}2()()()p M p C K X p F p ++= （3-2）
或 []{}{}()()()Z p X p F p = （3-3）
式中 []()Z p ：动刚度矩阵。

由式3-2）或（3-3）可以得出传递函数矩阵[]()H p ：
{}[]{}()()()X p H p F p = （3-4）
借助矩阵相关理论计算出来：
[][]
[]1
(())()()()
adj Z p H p Z p Z p -==
（3-5）
式中 [](())adj Z p ：为伴随矩阵； ()Z p ：为[]()Z p 的行列式。

式（3-5）的分母，叫做系统的特征方程。

类似单自由度系统，特征方程的根，即系统极点，决定系统的共振频率。

应为[]C 一般粘性阻尼矩阵，不一定满足矩阵可对角化条件，为了把系统方程（3-2）转化为一般特征值问题公式，需引入恒等式：
[][]{}{}()0p M p M X -= （3-6）
将此式与（3-2）合并：
[][]{}{}()p A B Y F '+= （3-7）
其中
[][]
[][][]0M A M C ⎡⎤=⎢
⎥
⎣⎦, [][][][][]00M B K ⎡⎤-=⎢⎥⎣⎦
,
{}{}{}p X Y X ⎧⎫⎪⎪=⎨⎬⎪⎪⎩⎭, {}{}{}0F F ⎧⎫⎪⎪
'=⎨⎬⎪⎪⎩⎭
令{}F '={}0, 则（3-7）的特征值满足下列方程：
[][]0p A B += （3-8）
对于N 自由度系统，此方程有2N 个复共轭对出现的特征根：
i i i
i
i i j j λσωλσω*
=-+⎧⎨=--⎩ 其中i σ阻尼因子；i ω为阻尼固有频率。

3.2 第一阶固有频率及主振型[17,18]
在求解系统动力响应时，系统较低的前几阶固有频率及相应的主振型占有重要的地位，为计算它们而采用下面的矩阵迭代法是比较简单的。

将i φφλλ==和i 带入公式中，得
i i i φλφ=A （3-9）
若将上式左端看作新列阵，上式表示：对于精确的主振型。

新列阵(i A φ) 与原来的列阵i φ的各个对应元素之间都相差同一常倍数，这个常倍数即特征值1λ。

记1X 为初始迭代列阵，由展开定理，1X 可以表示为
111φa X =n n a a φφ+++ 22 （3-10）
对上式左乘矩阵A ，由式（3-9）得知第一次迭代后所得的列阵为
n n n a a a AX X φλφλφλ+++== 22211112=⎪⎪⎭⎫ ⎝⎛+++n n n a a a φλλφλλφλ1212
2111 （3-11） 如果特征值1λ不是特征方程的重根，那么上式中的
1
n 1312λλ
λλλλ、、、 都小于1，因此比起其他主振型1φ在2X 内占的比重相对地比在1X 中占的比重大，换句话说，用矩阵A 迭代计算一次后，扩大了迭代列阵中第一阶主振型的优势。

经第二次迭代后，得。