毕业设计弹簧质量阻尼系统的建模与控制系统设计线性系统理论结课报告.doc

华北电力大学（保定）线性系统理论结课报告
分 数: ___________
任课教师签字：___________

华北电力大学研究生结课作业

学 年 学 期：第一学年第一学期
课 程 名 称：线性系统理论
学 生 姓 名：
学 号：
提 交 时 间：2014.11.27

目录
目录 2
1 研究背景及意义 3
2 弹簧-质量-阻尼模型 3
2.1 系统的建立 4
2.1.1 系统传递函数的计算 5
2.2 系统的能控能观性分析 7
2.2.1 系统能控性分析 8
2.2.2 系统能观性分析 9
2.3 系统的稳定性分析 10
2.3.1 反馈控制理论中的稳定性分析方法 10
2.3.2 利用Matlab分析系统稳定性 10
2.3.3 Simulink仿真结果 12
2.4 系统的极点配置 15
2.4.1 状态反馈法 15
2.4.2 输出反馈法 16
2.4.2 系统极点配置 16
2.5系统的状态观测器 18
2.6 利用离散的方法研究系统的特性 20
2.6.1 离散化定义和方法 20
2.6.2 零阶保持器 22
2.6.3 一阶保持器 24
2.6.4 双线性变换法 26
3.总结 28
4.参考文献 28

弹簧-质量-阻尼系统的建模与控制系统设计
1 研究背景及意义
弹簧、阻尼器、质量块是组成机械系统的理想元件。由它们组成的弹簧-质量-阻尼系统是最常见的机械振动系统，在生活中具有相当广泛的用途，缓冲器就是其中的一种。缓冲装置是吸收和耗散过程产生能量的主要部件，其吸收耗散能量的能力大小直接关系到系统的安全与稳定。缓冲器在生活中处处可见，例如我们的汽车减震装置和用来消耗碰撞能量的缓冲器，其缓冲系统的性能直接影响着汽车的稳定与驾驶员安全；另外，天宫一号在太空实现交会对接时缓冲系统的稳定与否直接影响着交会对接的成功。因此，对弹簧-质量-阻尼系统的研究有着非常深的现实意义。
2 弹簧-质量-阻尼模型
数学模型是定量地描述系统的动态特性，揭示系统的结构、参数与动态特性之间关系的数学表达式。其中，微分方程是基本的数学模型，不论是机械的、液压的、电气的或热力学的系统等都可以用微分方程来描述。微分方程的解就是系统在输入作用下的输出响应。所以，建立数学模型是研究系统、预测其动态响应的前提 。通常情况下，列写机械振动系统的微分方程都是应用力学中的牛顿定律、质量守恒定律等。
弹簧-质量-阻尼系统是最常见的机械振动系统。机械系统如图2.1所示，

图2-1弹簧-质量-阻尼系统机械结构简图
其中、表示小车的质量，表示缓冲器的粘滞摩擦系数，表示弹簧的弹性系数，表示小车所受的外力，是系统的输入即,表示小车的位移，是系统的输出，即，i=1,2。设缓冲器的摩擦力与活塞的速度成正比，其中，，，，，。
2.1 系统的建立
由图2.1，根据牛顿第二定律，分别分析两个小车的受力情况，建立系统的动力学模型如下：
对有：
对有：
联立得到：
对：
对：
令，，，，，;
，
得出状态空间表达式：

所以，状态空间表达式为：
+

由此可以得出

已知：，，，，，
代入数据得：

2.1.1 系统传递函数的计算
在Matlab中，函数ss2tf给出了状态空间模型所描述系统的传递函数，其一般形式是[num,den]=ss2tf(A,B,C,D,iu)，其中iu是输入值。
用Matlab将状态空间表达式表示为传递函数：
在输入1单独作用的情况下
A=[0 0 1 0;0 0 0 1; -400 300 -9 6;150 -200 3 -4.5];
B=[0 0;0 0;1 0;0 0.5];
C=[1 0 0 0;0 1 0 0];
D=[0 0;0 0];
[num,den]=ss2tf(A,B,C,D,1)
运行程序，得到：
num =

0 -0.0000 1.0000 4.5000 200.0000
0 -0.0000 -0.0000 3.0000 150.0000

den =

1.0e+004 *

0.0001 0.0014 0.0623 0.1800 3.5000
在输入2单独作用的情况下：
A=[0 0 1 0;0 0 0 1; -400 300 -9 6;150 -200 3 -4.5];
B=[0 0;0 0;1 0;0 0.5];
C=[1 0 0 0;0 1 0 0];
D=[0 0;0 0];
[num,den]=ss2tf(A,B,C,D,2)
运行程序，得到：
num =

0 -0.0000 -0.0000 3.0000 150.0000
0 -0.0000 0.5000 4.5000 200.0000

den =

1.0e+004 *

0.0001 0.0014 0.0623 0.1800 3.5000
由此可知：
位移对外力的传递函数是：

位移对外力的传递函数是：

位移对外力的传递函数是：

位移对外力的传递函数是：

2.2 系统的能控能观性分析
在反馈控制理论中只讨论输入量对输出量的控制。而这两个量的关系唯一地由系统的传递函数所确定。一个稳定的系统，一定能控。同时，系统的输出量本身就是我们想要控制的量，对于一个实际的系统来说，输出量当然是可以被观测到的，因此在反馈控制理论中没有必要设立能控和能观这两个概念。
然而在现代控制理论中，能控和能观是两个重要的基本概念。我们把反映系统内部运动状态的状态向量作为被控量，而且它们不一定是实际上可观测到的物理量，至于输出量则是状态向量的线性组合，这就产生了从输入量到状态量的能控性问题和从输出量到状态量的能观测性问题。
在现代控制中，分析和设计一个控制系统，必须研究这个系统的能控性和能观性。状态方程描述了输入U(t)引起状态X（t）的变化过程；输出方程则描述了由状态变化引起的输出Y（t）的变化。能控性和能观性正是分别分析U(t)对状态X（t）的控制能力以及Y（t）对X（t）的反应能力。
2.2.1 系统能控性分析
设线性定常系统的状态方程为

式中 A——n×n矩阵
B——n×r矩阵
C——m×n矩阵
D——m×r矩阵
系统能控的充分必要条件为：能控判别阵的秩R()=n，
用Matlab计算能控矩阵的秩，从而对该系统的能控性进行判别，程序为：
A=[0 0 1 0;0 0 0 1; -400 300 -9 6;150 -200 3 -4.5];
B=[0 0;0 0;1 0;0 0.5];
C=[1 0 0 0;0 1 0 0];
D=[0 0;0 0];
Qc=ctrb(A,B)
R1=rank(Qc)
运行程序，得到：

R1 =

4
等于矩阵行数，由此可以判断，系统是完全能控的。
2.2.2 系统能观性分析
设线性定常系统的状态方程为：

式中 A——n×n矩阵
B——n×r矩阵
C——m×n矩阵
D——m×r矩阵
能观的充分必要条件为：能观判别阵的秩R()=n，
下面，用Matlab计算能控矩阵的秩，从而对该系统的能控性进行判断：
A=[0 0 1 0;0 0 0 1; -400 300 -9 6;150 -200 3 -4.5];
B=[0 0;0 0;1 0;0 0.5];
C=[1 0 0 0;0 1 0 0];
D=[0 0;0 0];
Qo=obsv(A,C)
R2=rank(Qo)
运行程序，得到：
R2 =

4
满秩，因此可以判断，该系统是完全能观的。
综上所述，这是一个既能控又能观的系统。

2.3 系统的稳定性分析
2.3.1 反馈控制理论中的稳定性分析方法
稳定性是一个系统可以被采用的最基本的条件，是系统的固有属性。
稳定系统的定义如下：设控制系统处于某一起始的平衡状态，在外力的作用下，它离开了平衡状态，当外作用消失后，如果经过足够长的时间它能够恢复到起始的平衡状态，则称这样的系统为稳定的系统，否则称为不稳定的系统。由稳定性的定义可见，稳定性是系统去掉外力作用后自身的一种恢复能力，所以是系统的一种固有特性。对于线性定常系统，它取决于系统本身的结构和参数，而与初始条件和外界作用无关。
线性定常系统稳定的充分必要条件是：闭环系统特征方程的所有特征根为负实数或具有负实部的共轭复数，即所有特征根位于复平面的左半平面。只要有一个闭环特征根分布在右半平面上，系统就是不稳定的；如果没有右半平面的根，但有纯虚根，则系统是临界稳定的；在工程上，处于不稳定和临界稳定的线性定常系统是不能采用的[1]。
在古典控制系统中，我们判断系统的稳定性经常用劳斯-赫尔维茨代数判据、时域分析法、根轨迹法、频域分析法等方法，但那只针对低阶系统。实际的工业生产中，经常会遇见一些特别复杂的系统。这时古典控制理论中的方法就有点捉襟见肘了。
1892年俄国学者李雅普诺夫提出的稳定性理论是确定系统稳定性的更一般性理论，它采用了状态向量描述，不仅适用于单变量、线性、定常的系统，而且适用于多变量，非线性、时变的系统。李雅普诺夫理论在建立一系列关于稳定性概念的基础上，提出了判断系统稳定性的两种方法：一种方法是利用线性系统微分方程的解来判断系统稳定性，称为李雅普诺夫第一法或间接法；另一种方法是首先利用经验和技巧来构造李雅普诺夫函数，进而利用李雅普诺夫函数来判断系统稳定性，称为李雅普诺夫第二法或直接法。
2.3.2 利用Matlab分析系统稳定性
随着计算机技术的发展，在现代控制理论中，我们经常采用Matlab判断系统的稳定性。对于线性定常系统，典型的系统输入信号类型有脉冲、阶跃、斜坡、加速度、正弦信号。系统的稳定性是对任何输入信号而言，即若一个系统是稳定的，则其在任何输入信号情况下对应的输出曲线是收敛的。然而，阶跃信号包含了另外几种常见输入信号的特性，所以我们常通过观察系统的单位阶跃响应曲线判断判断系统的稳定性。若系统的单位阶跃响应是收敛的，则系统一般是收敛的；否则，是发散的。
在Matlab中输入相应系统的状态空间表达式矩阵来求取系统的特征值：
A=[0 0 1 0;0 0 0 1; -400 300 -9 6;150 -200 3 -4.5];
B=[0 0;0 0;1 0;0 0.5];
C=[1 0 0 0;0 1 0 0];
D=[0 0;0 0];
eig(A)
运行程序，得到：
ans =

-5.7735 +22.3859i
-5.7735 -22.3859i
-0.9765 + 8.0332i
-0.9765 - 8.0332i
由此可以知道，经计算得出A阵的所有特征根均在复平面的左半平面，因此得出该系统是稳定的。
给系统加起阶跃信号：
A=[0 0 1 0;0 0 0 1; -400 300 -9 6;150 -200 3 -4.5];
B=[0 0;0 0;1 0;0 0.5];
C=[1 0 0 0;0 1 0 0];
D=[0 0;0 0];
step(A,B,C,D)
结果如下

图2-2 阶跃响应曲线
由图可以看出，在阶跃响应下，系统在一定时间内收敛于某一固定值，因此可以判断系统是稳定的，但同时我们也可以看出，系统的调节时间比较长，如果想要减少调节时间，那么需要重新配置极点，对系统进行改进。下面的章节将对系统进行极点的配置。
2.3.3 Simulink仿真结果
根据上述原理，用Matlab中的Simulink组件进行仿真。
根据状态空间表达式，搭建系统模型如下图所示：
我们分别对只有输入1作用下和只有输入2作用下的系统使用Simulink进行仿真，让其与Matlab图像进行对比

图2-3 Simulink模型图
（1）仅有作用时，系统的输出如下图所示

图2-4 u1作用时响应曲线
图中，绿色为输出1的曲线，蓝色为输出2的曲线。经分析：此曲线与对应Matlab曲线一致，系统稳定，但是超调量较大，调节时间较长。
（2）仅有作用，系统的输入如下所示：

图2-5 u2作用时响应曲线
图中，绿色为输出1的曲线，蓝色为输出2的曲线。经分析：同样，此曲线与对应的Matlab曲线一致，系统稳定，但是超调量较大，调节时间较长。
在共同作用下，系统的输出如下图所示：

图2-6 u1、u2同时作用时响应曲线
图中绿色为输出1的曲线，蓝色为输出2的曲线。经分析：此曲线与Matlab曲线一致，系统稳定，但是超调量较大，调节时间较长。需要进行极点配置，使系统得到更好的性能。
2.4 系统的极点配置
控制系统的性能主要取决于系统极点在根平面上的分布。因此，在系统设计中，通常是根据对系统的品质要求，规定闭环极点应有的分布情况。所谓的极点配置就是，就是通过选择反馈矩阵K,将闭环系统的极点恰好配置在根平面上所期望的位置，以获得所希望的动态性能。
2.4.1 状态反馈法
极点问题首先解决是否能通过状态反馈来实现给定的极点配置，即在什么条件下才有可能按照规定的要求来配置极点。其次是，这样的反馈阵如何确定的问题。

图2-7 状态反馈示意图
（1）采用状态反馈配置系统极点条件：
系统采用状态反馈，任意配置其闭环系统极点的充要条件为：系统完全能控。若系统不是完全能控的，就必须按能控性分解，只能任意配置可控的极点。
（2）极点配置的方法：
若原系统可控，则采用状态反馈阵，有可控。
设原系统的特征方程为。
其中，则有：
，
配置后的闭环特征方程为：
；
假设闭环系统希望的极点为，得到：
。
为使系统达到希望性能，对比式（1）和式（2）中系数，使之相等，即可求得状态反馈阵。采用状态反馈配置系统极点不改变系统可控性，它不能影响系统中不可控部分模块。
2.4.2 输出反馈法

图2-8 输出反馈示意图
对于完全能控的单变量系统，不能采用输出线性反馈来实现闭环系统极点任意配置。不能任意配置极点，正是输出线性反馈的基本弱点。为了克服这个弱点，在经典控制理论中，往往采取引入附加校正网络，通过增加开环零极点的方法改变根轨迹走向，从而使其落在指定的期望位置上。
对于完全能控的单变量系统，通过带动态补偿器的输出反馈时限极点任意配置的充要条件是：1. 系统完全能观测；2. 动态补偿器的阶数为n-1。
2.4.2 系统极点配置
在现代控制理论中是用系统内部的状态来描述系统的，所以经常从系统的状态引出信号作为反馈量。利用状态反馈只能改变系统能控部分的极点，而不能改变不能控部分的极点，因此利用状态反馈进行极点配置的充分必要条件是系统必须是完全能控的。
对一个可控系统，在采用状态反馈后，可以实现闭环极点的任意配置，即通过状态反馈的方法，使闭环系统的极点位于任意期望的位置上。对于

其中x是状态变量（n维），u是控制信号，这里选取控制信号为

因此，
系统的稳态响应和瞬态响应特性由矩阵的特征决定
虽然理论上系统的闭环极点离S左半平面越远越好，但是在工业生产实践中，系统极点离左半平面越远，系统的运动状态就变化的越快，这就要求执行机构快速运作，即使再好的执行元件也会短时间内被损坏掉。所以新的极点的绝对值大约是原系统极点绝对值的3至4倍左右。取P1= -15+40i；P2= -15-40i；P3= -3+10i；P4= -3-10i；
利用Matlab进行极点配置，希望可以减小超调量，缩短稳定时间以优化系统。Matlab程序如下：
A=[0 0 1 0;0 0 0 1; -400 300 -9 6;150 -200 3 -4.5];
B=[0 0;0 0;1 0;0 0.5];
C=[1 0 0 0;0 1 0 0];
D=[0 0;0 0];
p=[-15+40i,-15-40i,-3+10i,-3-10i];
k=place(A,B,p)
step(A-B*k,B,C,D)
运行程序，得到：
k =

-234.6522 131.8512 14.4561 6.3957
643.3762 -89.9765 6.7658 36.0878
图2-9 稳态响应曲线
由响应曲线可以看出该系统重新配置极点后，具有较快的调节时间，而且也减少了超调量，改善了系统的动态性能与稳态性能。
2.5系统的状态观测器

图2-10 状态观测器示意图
通过状态观测器可以任意配置系统的极点，从而使闭环系统具有期望的稳态和动态性能。但在工业生产中，系统的状态变量并非都是物理量，或者是难以测得的量。这样一来，系统的所有状态变量未必都可以直接测量得到，因此，状态反馈这种控制方式在许多实际控制问题中往往难以直接应用和实现。状态观测器就是利用系统的外部输入输出信息来确定系统内部的状态，进而，在系统的极点配置状态反馈中，用观测器得到的状态估计值代替系统的真实状态。下图为状态观测器的结构图：

图2-11 状态观测器示意图
使用MATLAB为本系统设置状态观测器，选用极点配置时的极点，程序如下图所示：
A=[0 0 1 0;0 0 0 1; -400 300 -9 6;150 -200 3 -4.5];
B=[0 0;0 0;1 0;0 0.5];
C=[1 0 0 0;0 1 0 0];
D=[0 0;0 0];
p=[-15+40i,-15-40i,-3+10i,-3-10i];
K1=place(A,B,p)
A1=A-B*K1
L1=(place(A',C',p))'
A2=A-L1*C
L2=(place(A1',C',p))'
A3=A1-L2*C
sys2=ss(A2,B,C,D)
sys2=ss(A3,B,C,D)
运行上面程序，得到：
L1 =

7.0833 30.0895
-30.5796 15.4167
-41.6552 -96.5401
168.1877 200.0790
A2 =

-7.0833 -30.0895 1.0000 0
30.5796 -15.4167 0 1.0000
-358.3448 396.5401 -9.0000 6.0000
-18.1877 -400.0790 3.0000 -4.5000
L2 =

3.7432 -7.1200
-21.4563 -3.7432
190.9894 93.5822
115.5037 -24.2083
A3 =

-3.7432 7.1200 1.0000 0
21.4563 3.7432 0 1.0000
-655.5795 -119.9176 -18.2856 30.9515
-81.9216 -402.0612 -29.0527 -17.7144
其中L1代表没进行状态反馈时的状态观测反馈矩阵，L2代表进行了状态反馈的状态观测矩阵。
2.6 利用离散的方法研究系统的特性
2.6.1 离散化定义和方法
利用数字计算机对线性定常连续系统求数值解是现代科学技术研究中常用的一种方法，它不但方便，而且精确。由于实际工业生产中线性定常连续系统被控对象需要在线控制等，必须将连续时间系统的状态方程转化为离散系统的状态方程，即将矩阵微分方程化成矩阵差分方程，这就是连续系统的离散化。
根据离散系统的构成设备不同可以将离散系统分为采样控制系统和数字控制系统：
a.采样控制系统：控制系统的构成中选择了采样开关（或者含有开关特性的设备）。
b.数字控制系统：控制系统的控制器选择了专用数字计算机。
通常，把系统中的离散信号是脉冲序列形式的离散系统，称为采样控制系统或脉冲控制系统；而把数字序列形式的离散系统，称为数字控制系统或计算机控制系统。
采样控制系统：采样控制系统是对来自传感器的连续信息在某些规定上的时间瞬时值上取值。例如，控制器系统中的误差信号可以是断续连续的脉冲信号，而相邻两个脉冲之间的误差信息，系统并没有收到。如果在有规律的间隔上，系统取得了离散信息，则这种采样称为周期采样；反之，如果信息之间的间隔是时变的，或随机的，则称为非周期采样，或随机采样。
在采样控制系统中，把连续信号转变为脉冲序列的过程称为采样过程，简称采样。实现采样的装置称为采样器，或采样开关。用表示采样周期，单位为；，表示采样频率，单位为；表示采样角频率，单位为。在采样控制系统中，把脉冲序列转变为连续信号的过程称为信号复现过程。实现复现过程的装置称为保持器。
采样周期的选择满足香农采样定理。采样周期太大会使信号失真，采样周期太小则容易造成计算过程的累积偏差或失去采样系统的特性。香农采样定理是在设计离散系统时必须要遵循的准则，它给出了自采样的离散信号不失真地恢复原连续信号所必需的理论上的最低采样频率。采样频率应该满足

即是采样角频率，应使其对连续信号中的最高频率分量，在一个周期内被采样2次以上（上半周与下半周都至少采样一次），则采样后的脉冲序列中将包含了连续信号的全部信息。但是，在仿真中所遇到的大多数被再现信号是没有频带限的，所以一般取采样频率再现信号主要频带中的最高频率的5~10倍。在离散控制系统的设计过程中，采样周期的确定依据的是现场检测的被调量信号的频率，对于频率较高的信号，采样周期的设定就小，而对于变化过程较慢的低频信号，采样周期的设定可以大一些。有关概念在工程上的实际应用会有专门的内容介绍。
线性连续系统状态方程离散化的实质是将矩阵微分方程化为矩阵差分方程，它是描述多输入多输出离散系统的一种方便的数学模型。
在推导离散化系统的方程时，假定系统是周期性采样，并且采样脉冲宽度远小于采样周期，采样周期T的选择满足香农采样定理，还假设系统具有零阶保持特性，即在两个采样瞬间之间，采样值不变，并等于前一个采样时刻的值。
通常离散化的方法有很多，例如欧拉法，梯形法，龙哥-库塔（Runge-Kutta）法，阿达姆斯（Adams）法等等。下面主要运用三种方法来对系统进行离散化并运用计算机进行模拟系统的特性，分析不同采样周期对系统的影响效果。
2.6.2 零阶保持器
零阶保持器可以将脉冲序列变成连续的方波信号，即将前一个采样周期的数值保留到下一个采样点到来的时候。在Matlab中输入函数如下：
A=[0 0 1 0;0 0 0 1; -400 300 -9 6;150 -200 3 -4.5];
B=[0 0;0 0;1 0;0 0.5];
C=[1 0 0 0;0 1 0 0];
D=[0 0;0 0];
p=[-15+40i,-15-40i,-3+10i,-3-10i];
k=place(A,B,p);
[H,I,J,K]=c2dm(A-B*k,B,C,D,0.1,'zoh')
dstep(H,I,J,K)
分别设置采样时间为0.1s，0.05s，0.01s，运行程序，得到下图：

t=0.1s

t=0.05s
t=0.01s
图2-12 零阶保持器离散化
2.6.3 一阶保持器
采用一阶保持器进行离散化，程序如下:
A=[0 0 1 0;0 0 0 1; -400 300 -9 6;150 -200 3 -4.5];
B=[0 0;0 0;1 0;0 0.5];
C=[1 0 0 0;0 1 0 0];
D=[0 0;0 0];
p=[-15+40i,-15-40i,-3+10i,-3-10i];
k=place(A,B,p);
[H,I,J,K]=c2dm(A-B*k,B,C,D,0.1,'foh')
dstep(H,I,J,K)
同样，分别设置采样时间为0.1s，0.05s，0.01s，运行程序，得到下图：

t=0.1s

t=0.05s

t=0.01s
图2-13 一阶保持器离散化

2.6.4 双线性变换法
采用双线性变换法进行离散，程序如下：
A=[0 0 1 0;0 0 0 1; -400 300 -9 6;150 -200 3 -4.5];
B=[0 0;0 0;1 0;0 0.5];
C=[1 0 0 0;0 1 0 0];
D=[0 0;0 0];
p=[-15+40i,-15-40i,-3+10i,-3-10i];
k=place(A,B,p);
[H,I,J,K]=c2dm(A-B*k,B,C,D,0.1,'tustin');
dstep(H,I,J,K)
同样，分别设置采样时间为0.1s，0.05s，0.01s，运行程序，得到下图：

t=0.1s

t=0.05s

t=0.01s
图2-14 双线性变换法离散化
3.总结
本文选用一典型的弹簧-质量-阻尼模型，首先判别它的能控性、能观性，分析它的稳定性，并对其进行状态反馈和极点配置，又设计了状态观测器，最后对其进行了离散化，得到了离散化后取不同采样时间的仿真图形。
在做这次作业的过程中，我碰到了许多困难。由于在本科学校时，Matlab是一门选修课，所以这个软件对我来说很陌生，基本不太会。于是从图书馆借了几本书学习如何使用Matlab。其中对我帮助尤其大的是黄忠霖先生编写的《控制系统MATLAB计算及仿真》。专门把控制系统在Matlab中用到的相关内容作了详细的介绍。学习这本书的过程，我受益匪浅。通过这次作业，我不仅加深了对现代控制理论这门课的理解，还学习了一个非常有用的软件，这是意外收获。另外就是，模型到公式的转化让枯燥死板的课本知识与实际紧密的结合了起来。这样在以后的学习生活中，我会更加有动力，同时也为以后的研究学习还有工作打下了坚实的基础。
4.参考文献
[1] 于希宁．自动控制理论[M]．北京：中国电力出版社，2008．
[2] 黄忠霖．控制系统MATLAB计算及仿真[M]．北京：国防工业出版社，2001．
[3] 徐宝云．计算机建模与仿真技术[M]．北京：北京理工大学出版社，2009．
[4] 韩璞，刘长良．控制系统数字仿真技术[M]．北京：中国电力出版社，2007．
[5] 白艳艳，张晓俊．建立弹簧-质量-阻尼系统模型的数轴法[J]．噪声与振动控制，2012,6．

29