弹簧质量阻尼试验参考指导书.doc

质量-弹簧-阻尼系统试验教学指导书

北京理工大学 机械和车辆学院
.3

试验一：单自由度系统数学建模及仿真
1 试验目标
（1）熟悉单自由度质量-弹簧-阻尼系统并进行数学建模；
（2）了解MATLAB软件编程，学习编写系统仿真代码；
（3）进行单自由度系统仿真动态响应分析。
2 试验原理
f(t)
x
k
c
m

单自由度质量-弹簧-阻尼系统，如上图所表示。由一个质量为m滑块、一个刚度系数为k弹簧和一个阻尼系数为c阻尼器组成。系统输入：作用在滑块上力f(t)。系统输出：滑块位移x(t)。
建立力学平衡方程：

改变为二阶系统标准形式：

其中：ω是固有频率，ζ是阻尼比。

2.1 欠阻尼(ζ<1)情况下，输入f(t)和非零初始状态响应：

2.2 欠阻尼(ζ<1)情况下，输入f(t)=f0*cos(ω0*t) 和非零初始状态响应：

输出振幅和输入振幅比值：
3 动力学仿真
依据数学模型，使用龙格库塔方法ODE45求解，任意输入下响应结果。
仿真代码见附件
4 试验
4.1 固有频率和阻尼试验
（1）将试验台设置为单自由度质量-弹簧-阻尼系统。
（2）关闭电控箱开关。点击setup菜单，选择Control Algorithm，设置选择Continuous Time Control，Ts=0.0042，然后OK。
（3）点击Command菜单，选择Trajectory，选择step,进入set-up,选择Open Loop Step设置(0)counts, dwell time=3000ms,(1)rep, 然后OK。此步是为了使控制器得到一段时间数据，并不会驱动电机运动。
（4）点击Data菜单，选择Data Acquisition,设置选择Encoder#1 ，然后OK离开；从Utility菜单中选择Zero Position使编码器归零。
（5）从Command菜单中选择Execute，用手将质量块1移动到2.5cm左右位置（注意不要使质量块碰触移动限位开关），点击Run, 大约1秒后，放开手使其自由震荡，在数据上传后点击OK。
（6）点击Plotting菜单，选择Setup Plot，选择Encoder #1 Position；然后点击Plotting菜单，选择Plot Data，则将显示质量块1自由振动响应曲线。
（7）在得到自由振动响应曲线图上，选择n个连续振幅显著振动周期，计算出这段振动时间t，由n/t即可得到系统频率，将Hz转化为rad/sec即为系统振动频率ω。
（8）在自由振动响应曲线图上，测量步骤7选择时间段内初始振动周期振幅X0和末尾振动周期振幅Xn。由对数衰减规律即可求得系统阻尼比。
（9）试验数据统计
序号
第1次试验
第2次试验
第3次试验
试验测试频率

试验测试阻尼比

滑块质量m

弹簧刚度k

阻尼系数c

频率理论值

阻尼比理论值

频率估量误差

阻尼比估量误差

（10）在仿真代码基础上，计算出试验结果对应理论结果。对比分析理论和试验结果差异。完成试验汇报。
4.2 幅频特征试验
（1）点击Command菜单，选择Trajectory，选择Sinuscidal，进入set-up,选择Open LoopStep设置Amplitude(0.5V), Frequency(2Hz),Repetition(8)，然后OK。
（2）从Utility菜单中选择Zero Position使编码器归零。从Command菜单中选择Execute，点击Run,在数据上传后点击OK。
（3）然后点击Plotting菜单，选择Plot Data，则将显示滑块受迫振动响应曲线。在响应曲线图上，测量出振动振幅，计算出振动频率并于输入正弦曲线频率比较。
（4）依据试验情况，改变输入正弦曲线频率大小，反复上述，纪录试验数据。
输入频率
滑块试验幅值
滑块仿真幅值
0.1Hz

0.2Hz

（5）在仿真代码基础上，实现正弦激励代码，计算出试验结果对应理论结果。对比分析理论和试验结果差异。完成试验汇报。

试验二：双自由度系统数学建模及仿真
1 试验目标
（1）熟悉双自由度质量-弹簧-阻尼系统并进行数学建模；
（2）了解MATLAB软件编程，学习编写系统仿真代码；
（3）进行双自由度系统仿真动态响应分析。
2 试验原理
f(t)
x1
k1
c1
m1
x2
k2
c2
m2
2.1 数学建模

双自由度质量-弹簧-阻尼系统，如上图所表示。由两个质量为m1和m2滑块、两个刚度系数为k1和k2弹簧和两个阻尼系数为c1和c2阻尼器组成。系统输入：作用在滑块上力f(t)。系统输出：滑块位移x1(t)和x2(t)。
建立力学平衡方程：

2.2 固有频率
将动力学方程写成矩阵形式：

得到系统质量矩阵M和刚度矩阵K。
解行列式可得固有频率方程：
可计算出固有频率方程：

两个振动模态，两个固有频率：高模和低模。
2.3 解耦
经过数学变换将微分方程改变为以下形式：

注意：y1和y2不是滑块位移。滑块位移x1(t)和x2(t)是y1和y2函数。
3 动力学仿真
依据数学模型，使用龙格库塔方法ODE45求解，任意输入下响应结果。
仿真代码见附件
4 试验
4.1 固有频率分析
（1）将试验台设置为双自由度质量-弹簧-阻尼系统，第一个滑块没有阻尼器能够不接，认为第一个阻尼为零。
（2）闭合控制器开关，点击setup菜单，选择Control Algorithm，设置选择Continuous Time Control，Ts=0.0042,然后OK。点击Command菜单，选择Trajectory，选择step,进入set-up,选择Open Loop Step设置(0)counts, dwell time=3000ms, (1)rep, 然后OK。此步是为了使控制器得到一段时间数据，并不会驱动电机运动。
（2）点击Data菜单，选择Data Acquisition,设置分别选择Encoder#1，Encoder#2，然后OK离开；从Utility菜单中选择Zero Position使编码器归零。
（4）从Command菜单中选择Execute，用手将质量块1移动到2.5cm左右位置（注意不要使质量块碰触移动限位开关），点击Run, 大约1秒后，放开手使其自由振荡，在数据上传后点击OK。
（5）点击Plotting菜单，选择Setup Plot，分别选择Encoder #1 Position，Encoder #2
Position；然后点击Plotting菜单，选择Plot Data，则将显示质量块1，2自由振动响应曲线。
（6）试验数据纪录：
试验条件：滑块质量m1和m2，弹簧刚度k1和k2，阻尼系数c1和c2。
试验数据：时间-滑块1位移数据；时间-滑块2位移数据。
问题1：两个滑块位移频率测量值是高模和低模频率么？
问题2：实际机械系统是多自由度，怎样经过试验法测试系统固有频率？
（7）试验汇报。关键点是理论和试验结果对比分析。
4.2幅频特征试验
（1）点击Command菜单，选择Trajectory，选择Sinuscidal,进入set-up,选择Open Loop Step设置(200counts)Amplitude, Frequency(2Hz),Repetition（8），然后OK。
（2）从Utility菜单中选择Zero Position使编码器归零。从Command菜单中选择Execute，点击Run,在数据上传后点击OK。
（3）然后点击Plotting菜单，选择Plot Data，则将显示质量块1，2受迫振动响应曲线。在响应曲线图上，即可测量出振动振幅。
问题1：单自由度和双自由度系统幅频特征有何差异？
问题2：高模贡献分析。

试验三：PID控制
1 试验目标
（1）学习PID闭环控制结构和系统闭环传输函数计算；
（2）PID控制器参数设计；
（3）控制性能分析。
2 试验原理

上图给出闭环控制系统原理框图。单自由度质量-弹簧-阻尼系统结构下，断开弹簧和阻尼，仅仅保留滑块质量m。

电控箱能够看做百分比增益khw。

其中：u是控制器输出。PID控制：

其中：e是比较器输出，参考输入和实际输出偏差值。

依据全部上式，可得闭环结构微分方程：

对应传输函数：

3 PID设计
PID控制器中设置积分因子ki为零，则为PD控制。传输函数变为：

闭环特征方程是分母：

设计频率ω=4Hz，三种阻尼（欠阻尼ζ=0.2，临界ζ=1.0，过阻尼ζ=2.0）控制器设计。
ω=4 Hz
khwkp
khwkd
ζ=0

ζ=0.2

ζ=1

ζ=2.0

4 试验
4.1 频率
（1）在控制器断开情况下，拆除和质量块1连接弹簧，使其它元件远离质量块1运动范围，为其安装4个500g铜块，加上小车本身质量，标定总质量m=2.6kg。
（2）试验标定khw值：依据估量出khw值，设置控制器ki=0和kd=0，调整kp来估量系统系统频率ω=4Hz。注意： kp不能大于0.08。
（3）闭合控制器开关；点击Data菜单，选择Data Acquisition,设置选择Encoder#1和Commanded Position information；点击Command菜单，选择Trajectory，选择step，设置(0)counts,dwell time=3000ms,(1)rep。
（4）点击setup菜单，选择Control Algorithm，设置选择Continuous Time Control，Ts=0.0042,选择PID，进入Setup Algorithm，输入kp值（ki=0和kd=0）（输入值不能大于kp=0.08），然后OK。移动质量块1到-0.5cm位置（要求，朝电机方向为负）选择Implement Algorithm，然后Ok。
注意：以后步开始每一步，要进行下一步之前，全部要和运动装置保持一定安全距离；选择Implement Algorithm后，控制器将会立即加载，若出现不稳定或很大控制信号时，运动装置可能反应很猛烈；若加载后，系统看上去稳定，要先用一轻质不尖锐物体轻轻碰触质量块以验证其稳定性。
（5）点击Command菜单，进入Execute，用手将质量块移动到2cm左右位置，点击Run,移动质量块大约到3cm位置，然后释放（不要拿着质量块多于一秒，以免电机过热而断开控制）。
（6）点击Plotting菜单，选择Setup Plot，选择Encoder #1；然后点击Plotting菜单，选择Plot Data，则将显示质量块1时间响应曲线。
试想一下，若将百分比增益系数kp增加一倍，则系统响应频率将有什么改变？
4.2 阻尼
（1）确定kd值（不能大于0.04），使得khwkd=50N/(m/s)，反复步骤4，除了输入kd值和ki=0和kp=0。
（2）先用尺子检验系统稳定性，然后用手往返移动质量块来感受系数kd带来粘性阻尼影响（注意不要极度迫使质量块运动，以免电机过热而断开控制）。
（3）增大kd值（kd<0.04），反复以上步骤，看能否感受到阻尼增大？
4.3 位置控制
（1）控制电机驱动滑块1移动到某一个特定位置，获取时间-位移数据后，计算闭环控制特征，包含：超调量、上升时间、调整时间、稳态误差。
（2）使用仿真方法计算得到一样条件下动态响应结果，对比分析理论和试验结果差异。
4.4 正弦激励
（1）使用正弦输入驱动滑块作正弦运动，获取时间-位移数据后，计算闭环控制系统幅频特征。
（2）使用仿真方法计算得到一样条件下幅频响应结果，对比分析理论和试验结果差异。
4.5 柔性结构控制
（1）试验台只有质量块，为刚性结构机械系统。假如安装上一根弹簧，则组成简单柔性结构系统（含有柔性步骤）。此时机械系统因为能量交互可能产生振动。怎样设计控制器来动态控制机械振动。

试验四：扰动抑制
1 试验目标
（1）外扰动设置；
（2）PID控制器设计；
（3）控制性能分析。
2 试验原理
d
x
r

＋
＋

上图给出闭环控制系统原理框图。单自由度质量-弹簧-阻尼系统结构下，断开弹簧和阻尼，仅仅保留滑块质量m。位移输出受到外扰动d干扰。

电控箱能够看做百分比增益khw。

其中：u是控制器输出。PID控制：

其中：e是比较器输出，参考输入和实际输出偏差值。

依据全部上式，可得闭环结构微分方程：

对应传输函数：

3 PID设计
PID控制器中设置积分因子ki为零，则为PD控制。传输函数变为：

闭环特征方程是分母：

附件1: 试验1仿真代码
% exp 1
clc
clear
global SK SC SM SF;
sampling=1/1000;
SK=200;
SC=0;
SM=2.6;
State=zeros(1,2);
State(1)=1;
State(2)=0.0;
sss=zeros(1,1);
for k=1:fix(5/sampling)
t=k*sampling;
sss(k,1)=t; %time

SF=0;
sss(k,2)=SF;

[TimeOdeArray,VarOdeArray] = ode45( @mdlDerivatives1, [t t+sampling],State);
[m,n]=size(TimeOdeArray);
TimeAtEndOfArray = TimeOdeArray(m,1);
if( abs(TimeAtEndOfArray - (t+sampling) ) >= abs(0.001*sampling) ) warning('numerical integration failed'); break; end
VAR = VarOdeArray(m, : );%Ò»¸ö²½³¤½áÊøµÄ״̬
State=VAR;%³õʼ״̬¸üУ¬½øÐÐÏÂÒ»´Îµü´ú
sss(k,3)=VAR(1); % response
sss(k,4)=VAR(2); % velocity
end

plot(sss(:,1),sss(:,2),'r',sss(:,1),sss(:,3),'k',sss(:,1),sss(:,4),'b');

% mdlDerivatives1
function dx=mdlDerivatives1(T,x)
global SK SC SM SF;
% y(1)=x;
% y(2)=xd;
dx=zeros(2,1);
%sloshing dynamics
dx(1)=x(2);
dx(2)=SF/SM-SK/SM*x(1)-SC/SM*x(2);
附件2: 试验2仿真代码
% exp 2
clc
clear
global SK1 SK2 SC1 SC2 SM1 SM2 SF;
sampling=1/1000;
SK1=200;
SK2=200;
SC1=0;
SC2=0;
SM1=2.6;
SM2=2.6;
State=zeros(1,4);
State(1)=1.0;
State(3)=0.0;
sss=zeros(1,1);
for k=1:fix(5/sampling)
t=k*sampling;
sss(k,1)=t; %time

SF=0;
sss(k,2)=SF;

[TimeOdeArray,VarOdeArray] = ode45( @mdlDerivatives2, [t t+sampling],State);
[m,n]=size(TimeOdeArray);
TimeAtEndOfArray = TimeOdeArray(m,1);
if( abs(TimeAtEndOfArray - (t+sampling) ) >= abs(0.001*sampling) ) warning('numerical integration failed'); break; end
VAR = VarOdeArray(m, : );%Ò»¸ö²½³¤½áÊøµÄ״̬
State=VAR;%³õʼ״̬¸üУ¬½øÐÐÏÂÒ»´Îµü´ú
sss(k,3)=VAR(1); % response 1
sss(k,4)=VAR(2); % velocity 1
sss(k,5)=VAR(3); % response 2
sss(k,6)=VAR(4); % velocity 2
end
plot(sss(:,1),sss(:,2),'r',sss(:,1),sss(:,3),'k',sss(:,1),sss(:,5)-sss(:,3),'b');

% mdlDerivatives2
function dx=mdlDerivatives2(T,x)
global SK1 SK2 SC1 SC2 SM1 SM2 SF;
% x(1)=x1;
% x(2)=x1d;
% x(3)=x2;
% x(4)=x2d;
dx=zeros(4,1);
%sloshing dynamics
dx(1)=x(2);
dx(2)=(SF-SC1*x(2)+SC1*x(4)-SK1*x(1)+SK1*x(3))/SM1;
dx(3)=x(4);
dx(4)=(-1*SC2*x(4)-SC1*x(4)+SC1*x(2)-SK1*x(3)-SK2*x(3)+SK1*x(1))/SM2;

附件3 快速傅里叶变换FFT计算频率代码
%frequency response
sampling=1/1000;
Fs=fix(1/sampling);
data=sss(1:fix(5/sampling),3);
m=length(data);
nfft=2^nextpow2(m);
y=fft(data,nfft);%¸µÀïÒ¶±ä»»
Ayy=abs(y)*2/nfft;
f0=(0:nfft/2-1)*Fs/nfft; % Hz , 1/s
f1=f0*2*pi; % rad/s
zzz(:,1)=f0;
zzz(:,2)=f1;
zzz(:,3)=Ayy(1:nfft/2);
plot(zzz(:,1),zzz(:,3));
xlabel('frequency/Hz');
ylabel('magnitude');
title('freqency response');

1 Command/Trajectory
输入命令
2 Command/Execute
实施
3 plotting/Setupplot
绘图运动轨迹

4 Setup改变控制器