行驶动力学建模、仿真及主动悬架控制器设计

目录1. 计算机仿真系统模型的建立 .......................................................... - 1 -2. LOG控制器设计 .............................................................................. - 2 -3. 计算实例........................................................................................... - 3 -4. MATLAB仿真过程.......................................................................... - 4 -5. 半车模型建模及仿真 ...................................................................... - 8 -5.1随机线性最优控制 ................................................................... - 9 -5.2预瞄控制 ................................................................................. - 11 -5.3结果比较 ................................................................................. - 12 -以单轮车辆模型为例，介绍行驶动力学计算机建模、仿真分析以及利用线性二次最优控制理论进行主动悬架LQG 控制器设计过程。

1. 计算机仿真系统模型的建立根据图7所示的主动悬架单轮车辆模型，运用牛顿运动定律，建立系统的运动方程，即： ()b b a s b w m x U K x x =--&& （4）()()w w a s b w t w g m x U K x x K x X =-+---&& （5）这里，采用一个滤波白噪声作为路面输入模型，即：00()2()2()g g x t f x t G uw t ππ=-+&& （6）式中，xg 为路面垂向位移（m ）；G0为路面不平度系数（m3/cycle ）；u 为车辆前进速度（m/s ）；w 为数字期望为零的高斯白噪声；f0为下截止频率（Hz ）。

图7 单轮车辆模型 结合式（4）、式（5）和式（6），将系统运动方程和路面输入方程写成矩阵形式，即得出系统的空间状态方程：X AX BU FW =++&（7）式中，()Tb w b w g X x x x x x =&&，为系统状态矢量；W=（w （t ）），为高斯白噪声输入矩阵；U=（Ua （t ）），为输入控制矩阵；0000001000001000002ss b b s t sa w ww K K m m K K K K m m m A f π⎡⎤-⎢⎥⎢⎥⎢--⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥-⎣⎦；11000b w m m B ⎡⎤⎢⎥⎢⎥⎢⎥-⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦；00002F ⎡⎤⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎢⎣2. LOG 控制器设计车辆悬架设计中的主要指标包括：①代表轮胎接地性的轮胎动载荷；②代表轮胎舒适性的车身垂向振动加速度；③影响车身姿态且与轮胎布置有关的悬架动行程。

因此，LQG 控制器设计中的性能指标J 即为轮胎动位移、悬架动行程和车身垂向振动加速度的加权平方和在时域T 内的积分值，其表达式为：22212301lim {[()()][()()]()}Tw g b w b tT J q x t x t q x t x t q x t d T →∞=-+-+⎰&& （8）式中，q1、q2和q3分别为轮胎动位移、悬架动行程和车身垂向振动加速度的加权系数。

加权系数的选取决定了设计者对悬架性能的倾向，如对车身垂向振动加速度项选择较大的权值，则考虑更多的是提高车辆操纵稳定性。

为方便起见，这里取车身垂向振动加速度的加权值q3=1。

将性能指标J 的表达式（8）改写成矩阵形式，即：1lim (2)TT T T tT J X QX U RU X NU d T →∞=++⎰ （9）式中，222222222121221100000000000000000s s bb s s b bK K q q m m Q K K q q q q m m q q ⎡⎤⎢⎥⎢⎥⎢⎥+--⎢⎥=⎢⎥⎢⎥⎢⎥--++-⎢⎥⎢⎥-⎣⎦；21b R m =；000a s N K K ⎡⎤⎢⎥⎢⎥⎢⎥=-⎢⎥⎢⎥⎢⎥⎣⎦ 当车辆参数值和加权系数值确定后，最优控制反馈增益矩阵可有黎卡提（Riccati ）方程求出，其形式如下：1()()0T T TPA A P PB N R B P N Q -+-+++= （10） 最优反馈控制增益矩阵T TK B P N =+，由车辆参数和加权系数决定。

根据任意时刻的反馈状态变量X （t ），就可得到t 时刻作动器的最优控制力Ua ，即：()()a U t KX t =- （11）3. 计算实例这里，以某轿车的后悬架为例，给出一个完整的计算实例，包括车辆模型参数、仿真路面输入参数、控制器的设计参数以及计算结果。

此例中车辆以20m/s 的速度在某典型路面上行驶，仿真时间T=50s 。

计算中输入的各参数及数值详见表2。

仿真计算中以式（6）所示的滤波白噪声作为路面输入模型。

白噪声的生成可直接调用MATLAB 函数WGN （M ，N ，P ）（此函数需要安装信号处理工具箱Communications toolbox ），其中M 为生成矩阵的行数，N 为列数，P 为白噪声的功率（单位为dB ）。

本例中取M=10001，N=1，P=20。

这意味着仿真计算中去一条白噪声，共10001个采集点，噪声强度为20dB 。

设定采样时间为0.005s 、车速为20m/s 时，相当于仿真路面长度为1000m ，仿真时间为50s 。

根据建立的系统状态方程式（7）及最优化性能指标函数式（9），利用已知的矩阵A 、B 、Q 、R 、N ，调用MATLAB 中的线性二次最优控制器设计函数[K ，S ，E]=LQR （A ，B ，Q ，R ，N ），即可完成最优主动悬架控制器的设计。

输出的结果中，K 为最优控制反馈增益矩阵，S 为黎卡提方程的解，E 为系统闭环特征根。

根据表2给出的仿真输入参数，本例中求得的最优反馈增益矩阵K 为： K=（711.88 -1241.5 -19284 -2038.5 20864） 同时，还得到了黎卡提方程的解：2.45590.0289 2.47458.66077.3090.02890.48860.02987.52627.23642.47450.0298 4.97448.6754 5.10338.66077.52628.67542710.12700.47.3097.2364 5.10332700.42693.7S -⎡⎤⎢⎥-⎢⎥⎢⎥=-⎢⎥---⎢⎥⎢⎥--⎣⎦在相同的仿真条件下，可将所设计的主动悬架系统与一个被动系统进行对比分析。

在被动悬架系统中，取悬架刚度Ks=22000N/m ，阻尼系数Cs=1000NS/m 。

除此之外，其他输入参数值均与主动悬架系统完全相同。

4.MATLAB仿真过程1）生成路面输入模型代码如下：a=wgn(10001,1,20);t=0:0.005:50;road_file(:,1)=t';road_file(:,2)=a;save road_file road_file2）参数输入代码如下：load road_file.mat%载入路面数据模型Ks=22000;mb=320;Kt=200000;mw=40;f0=0.1;G0=0.000005;u=20; Kb=20000;Ks1=22000;Cs=1000; %输入仿真有关参数A=[0,0,-Ks/mb,Ks/mb,0; %建立主动悬架的状态矩阵0,0,Ks/mw,(-Kt-Ks)/mw,Kt/mw;1,0,0,0,0;0,1,0,0,0;0,0,0,0,-2*pi*f0];A1=[-Cs/mb,Cs/mb,-Ks1/mb,Ks1/mb,0; %建立被动悬架的状态矩阵Cs/mw,-Cs/mw,Ks1/mw,(-Kt-Ks1)/mw,Kt/mw;1,0,0,0,0;0,1,0,0,0;0,0,0,0,-2*pi*f0];B=[1/mb,0;-1/mw,0;0,0;0,0;0,2*pi*sqrt(G0*u)];B1=[0,0;0,0;0,0;0,0;0,2*pi*sqrt(G0*u)];C=[1,0,0,0,0;0,1,0,0,0;0,0,1,0,0;0,0,0,1,0;0,0,0,0,1];D=[0,0;0,0;0,0; 0,0; 0,0];K=[711.88,-1241.5,-19284,-2038.5,20864]; K1=[0,0,0,0,0];3）用Simulink 创建仿真框图 ● 状态变量],,,,[g w b w b x x x x x X ••= ● 输入与系统模块，如下图：● 输出模块，如下图：gg w w b b x x x x x x dtd =-=-==•路面输入轮胎动位移悬架动行程车身加速度)(整体程序框图如下：4）结果分析可以直接通过双击scope查看输出的波形图，为更好比较主动悬架与被动悬架的差别，下面通过输出到workspace的状态变量编程绘图并计算均方根值。

代码如下：%% 绘制车身加速度曲线，并计算均方根值%ba-主动悬架车身加速度%ba1-被动悬架车身加速度ba=diff(X.data(:,1))./diff(X.time);ba1=diff(X1.data(:,1))./diff(X1.time);subplot(2,1,1)plot(X.time(1:end-1),ba)subplot(2,1,2)plot(X1.time(1:end-1),ba1)BA=norm(ba,2)./(length(ba).^0.5);BA1=norm(ba1,2)./(length(ba1).^0.5);%% 绘制悬架动行程曲线，并计算其均方根值%sws-主动悬架动行程%sws1-被动悬架动行程figure()sws=X.data(:,3)-X.data(:,4);sws1=X1.data(:,3)-X1.data(:,4); subplot(2,1,1) plot(X.time,sws) subplot(2,1,2) plot(X.time,sws1)SWS=norm(1000*sws,2)./(length(sws).^0.5); SWS1=norm(1000*sws1,2)./(length(sws1).^0.5); %% 绘制轮胎动位移曲线，并计算其均方根值 %dtd-主动悬架动位移 %dtd1-被动悬架动位移 figure()dtd=X.data(:,4)-X.data(:,5); dtd1=X1.data(:,4)-X1.data(:,5); subplot(2,1,1) plot(X.time,dtd) subplot(2,1,2) plot(X.time,dtd1)DTD=norm(1000*dtd,2)./(length(dtd).^0.5); DTD1=norm(1000*dtd1,2)./(length(dtd1).^0.5);结果如下：① 车身加速度曲线时间(s)车身加速度(m /s 2)主动悬架车身加速度时间(s)车身加速度(m /s 2)被动悬架车身加速度② 悬架动行程曲线时间(s)位移(m )主动悬架动行程时间(s)位移(m )被动悬架动行程③ 轮胎动位移时间(s)位移(m )主动悬架轮胎动位移时间(s)位移(m )被动悬架轮胎动位移● 主动悬架与被动悬架性能指标均方根值比较5. 半车模型建模及仿真● 半车模型建模悬架的半车模型有四个自由度，可选取前后轮垂直位移、前后悬架与车身连接处垂直位移四个自由度，写出状态空间方程：以],,,,,,,,,[201012341234x x z z z z z z z z X ••••=作为系统状态变量，状态空间方程如下：白噪声输入。