汽车半悬挂系统建模与分析(现代控制理论大作业)

XX大学现代控制理论——汽车半主动悬架系统的建模与分析姓名：XXX学号：XXXX专业：XXXX一．课题背景汽车的振动控制是汽车设计的一个重要研究内容，涉及到汽车的平顺性和操纵稳定性。悬架系统是汽车振动系统的一个重要子系统，其振动传递特性对汽车性能有很大影响。因此设计性能良好的悬架系统以减少路面激励的振动传递，从而提高汽车的平顺性和操纵稳定性是汽车振动控制研究的重要课题。悬架系统是汽车车身与轮胎间的弹簧和避震器组成整个支撑系统，用于支撑车身，改善乘坐舒适度。而半主动悬架是悬架弹性元件的刚度和减振器的阻尼系数之一可以根据需要进行调节控制的悬架。目前，半主动悬架研究主要集中在调节减振器的阻尼系数方面，即将阻尼可调减振器作为执行机构，通过传感器检测到汽车行驶状况和道路条件的变化以及车身的加速度，由ECU根据控制策略发出脉冲控制信号实现对减振器阻尼系数的有级可调和无级可调。二．系统建模与分析1.1半主动悬架系统的力学模型以二自由度1/4半主动悬架模型为例，并对系统作如下假设：(1)悬挂质量与非悬挂质量均为刚体；(2)悬架系统具有线性刚度和阻尼；(3)悬架在工作过程中不与缓冲块碰撞；(4)轮胎具有线性刚度，且在汽车行驶过程中始终与地面接触。综上，我们将该系统等效为两个质量块M，m；两个弹簧系统Ks，Kt；一个可调阻尼器（包含一个常规阻尼器Cs和一个变化阻尼力F），如图1所示。图1系统力学模型1.2半主动悬架系统的数学模型由减振器的简化模型得:NSFCVF对m进行分析：211201122()tssdzdzdzmKzzKzzCFdtdtdt即：1011212()tssmzKzzKzzCzzF对M进行分析：2212122()ssdzdzdzMKzzCFdtdtdt即：21212ssMzKzzCzzF选取状态变量：1102213142xzzxzzxzxz，，，输入变量：uF输出变量：1122yxyx，综上可得，系统状态空间表达式为：11032214331234423411tssssssxzzxxzzxxKKCCxxxxxFmmmmmKCCxxxxFMMMM整理得：0010000110110tssssssKKCCmmmmmKCCMMMMxxu10000100yx三．数值化分析选取系统参数为：M=391kg，m=50.7kg，Ks=60KN/m，Kt=362KN/m，Cs取1KN·s/m。状态空间表达式变为：001000011071401183.4319.7219.72-0.020-153.452.56-2.560.0026xxu10000100yx四．能控性与能观性分析00100001101000,,71401183.4319.7219.72-0.0201000-153.452.56-2.560.0026Abc4.1能控性分析能控性矩阵：23(,,,)MbAbAbAb通过matlab计算得：Rank(M)=4，满秩，故系统可控。4.2能观性分析能观性矩阵：23,,,TNCCACACA通过matlab计算得：Rank(N)=4，满秩，故系统可观。五．稳定性分析存在唯一平衡点x=0，对矩阵A进行特征值计算：通过MATLAB计算，我们得到特征值为：-10.2018+90.5683i，-10.2018-90.5683i，-0.9382+11.4463i，-0.9382-11.4463i。由于矩阵A的特征值均有负实部，所以系统是大范围渐近稳定的。六．状态观测器设计因为系统完全能观，所以可以设计状态观测器。6.1全维观测器将系统极点配置为：-1，-2，-3，-4.MATLAB程序：A=[0,0,1,0;0,0,-1,1;-7140,1183.43,-19.72,19.72;0,-153.45,2.56,-2.56]；b=[0;0;-0.02;0.0026]；c=[1,0,0,0;0,1,0,0]；opt=[-1,-2,-3,-4];G=(place(A,c',opt))’;输出结果为：7877.00131.08140.09822.03395.101621.13873.789405.1G所以，全维观测器方程为：001001.940578.3873001101.162110.3395ˆˆˆ()71401183.4319.7219.72-0.020.98220.81400-153.452.56-2.560.00260.01310.7877yyxxu+6.2降维观测器由于rank（c）=2，n=4,所以将系统极点配置为-1，-2.构造变换阵作线性变换，设10010001000010001,1000100001000100TT。则，119.719.771401183.42.62.60153.410001100ATAT，10.020.00260010,000010BTBCCT。MATLAB程序：opt2=[-1,-2];T=[0,0,1,0;0,0,0,1;1,0,0,0;0,1,0,0;];Tni=inv(T);A_2=Tni*A*T;B_2=Tni*B;C_2=C*T;A_11=A_2(1:2,1:2);A_21=A_2(3:4,1:2);G2=(place(A_11',A_21',opt2))';输出结果为：56.0272.1912G。所以，降维观测器方程为：110.980.0271401183.40.02ˆˆ0.042.040153.40.0026119.72ˆˆ20.56wxyuxwy七．最优控制对于半主动悬架系统，最优控制器的设计目的就是寻找最优控制F，使实现控制所需的能量为最小：dtxxqxqJ024222211，其中，1q，2q分别为轮胎动变形加权系数，悬架动挠度加权系数，为车身加速度加权系数。将目标性能泛函改写成二次型性能指标形式：dtRuuQxxJTT0)(，这里，0000000000000021qqQ,为半正定常数矩阵；21MR，为正定常数矩阵。所以，最优控制存在，且唯一：)()(1tPxBRtuT式中，P为44维正定常数矩阵，满足黎卡提矩阵代数方程：01QPBPBRPAPATT采用试探法取三组不同权系数1q、2q，运用MATLAB进行计算分析：（1）q1=3.35e5，q2=40.5e5；（2）q1=3.35e8，q2=40.5e8；（3）q1=3.35e9，q2=40.5e9；Matlab程序：%最优控制clc;clear;M=391;A=[0,0,1,0;0,0,-1,1;-7140,1183.43,-19.72,19.72;0,-153.45,2.56,-2.56];B=[0;0;-0.02;0.0026];C=[1,0,0,0;0,1,0,0];D=0;R=1/M^2;%求不同Q、R下的状态反馈阵KQ1=3.35e5;Q2=40.5e5;Q=[Q1,0,0,0;0,Q2,0,0;0,0,0,0;0,0,0,0];[KPe]=lqr(A,B,Q,R)Ac=(A-B*K);Bc=B;Cc=C;Dc=D;T=0:0.05:5;U=0.2*ones(size(T));[Y,X1]=lsim(Ac,Bc,Cc,Dc,U,T);Q1=3.35e8;Q2=40.5e8;Q=[Q1,0,0,0;0,Q2,0,0;0,0,0,0;0,0,0,0];[KPe]=lqr(A,B,Q,R)Ac=(A-B*K);Bc=B;Cc=C;Dc=D;T=0:0.05:5;U=0.2*ones(size(T));[Y,X2]=lsim(Ac,Bc,Cc,Dc,U,T);Q1=3.35e9;Q2=40.5e9;Q=[Q1,0,0,0;0,Q2,0,0;0,0,0,0;0,0,0,0];[KPe]=lqr(A,B,Q,R)Ac=(A-B*K);Bc=B;Cc=C;Dc=D;T=0:0.05:5;U=0.2*ones(size(T));[Y,X3]=lsim(Ac,Bc,Cc,Dc,U,T);figure;holdon;plot(T,X1(:,1),'--','color','black');plot(T,X2(:,1),'-','color','green');plot(T,X3(:,1),'-.','color','red');xlabel('时间（s）');ylabel('轮胎动变形（m）');holdoff;legend('q1=3.35e5，q2=40.5e5','q1=3.35e8，q2=40.5e8','q1=3.35e9，q2=40.5e9');figure;holdon;plot(T,X1(:,2),'--','color','black');plot(T,X2(:,2),'-','color','green');plot(T,X3(:,2),'-.','color','red');xlabel('时间（s）');ylabel('悬架动挠度（m）');holdoff;legend('q1=3.35e5，q2=40.5e5','q1=3.35e8，q2=40.5e8','q1=3.35e9，q2=40.5e9');figure;holdon;plot(T,X1(:,3),'--','color','black');plot(T,X2(:,3),'-','color','green');plot(T,X3(:,3),'-.','color','red');xlabel('时间（s）');ylabel('悬架动载荷（N）');holdoff;legend('q1=3.35e5，q2=40.5e5','q1=3.35e8，q2=40.5e8','q1=3.35e9，q2=40.5e9');figure;holdon;plot(T,X1(:,4),'--','color','black');plot(T,X2(:,4),'-','color','green');plot(T,X3(:,4),'-.','color','red');xlabel('时间（s）');ylabel('车身加速度（m/s2）');holdoff;legend('q1=3.35e5，q2=40.5e5','q1=3.35e8，q2=40.5e8','q1=3.35e9，q2=40.5e9');matlab仿真结果如下：00.511.522.533.544.55-6-4-2024681012x10-8时间（s）轮胎动变形（m）q1=3.35e5，q2=40.5e5q1=3.35e8，q2=40.5e8q1=3.35e9，q2=40.5e9图2轮胎动变形变化趋势00.511.522.533.544.5500.511.522.53x10-7时间（s）悬架动挠度（m）q1=3.35e5，q2=40.5e5q1=3.35e8，q2=40.5e8q1=3.35e9，q2=40.5e9图3悬架动挠度的变化趋势00.511.522.533.544.55-3-2-1012345x10-6时间（s）悬架动载荷（N）q1=3.35e5，q2=40.5e5q1=3.35e8，q2