提示:文章写完后,目录可以自动生成,如何生成可参考右边的帮助文档
前言
提示:这里可以添加本文要记录的大概内容:
本篇博文主要用来记录汽车系统动力学课程的学习,仅供课程阶段参考,与实际工程应用还有较大距离。基于Simulink的魔术公式轮胎动力学仿真
本实例主要有以下任务。
1.建立魔术公式轮胎动力学仿真模型。
2.绘制不同载荷下的轮胎纵向力—滑动率关系曲线。
3.绘制不同载荷下的轮胎侧向力—侧偏角关系曲线。
4.绘制不同载荷下的轮胎回正力矩—侧偏角关系曲线。
提示:以下是本篇文章正文内容,下面案例可供参考
一、魔术公式
基于魔术公式的轮胎动力学仿真所需参数见表2-1。
纯工况Pacejka89轮胎模型计算轮胎纵向力、侧向力和回正力矩时共需要43个特性参数,纵向力需要确定11个特性参数;侧向力需要确定14个特性参数;回正力矩需要确定18个回正力矩特性参数。目前应用广泛的参数如表(2-1)所示。
表2-1基于魔术公式的轮胎动力学仿真所需参数
纯纵滑工况下,Pacejka89轮胎纵向力表达式为:
F
x
0
=
D
x
sin
(
C
x
arctan
(
B
x
s
x
0
−
E
(
B
x
s
x
0
−
arctan
(
B
x
s
x
0
)
)
)
)
+
s
v
x
(
1
)
{F_{x0}} = {D_x}\sin ({C_x}\arctan ({B_x}{s_{x0}} - E({B_x}{s_{x0}} - \arctan ({B_x}{s_{x0}})))) + {s_{vx}}(1)
Fx0=Dxsin(Cxarctan(Bxsx0−E(Bxsx0−arctan(Bxsx0))))+svx(1)
式中,
s
x
0
=
s
x
+
s
h
x
{s_{x0}} = {s_x} + {s_{hx}}
sx0=sx+shx 为纵向力组合自变量;
s
x
{s_x}
sx为轮胎纵向滑动率,当
s
x
>
0
{s_x} > 0
sx>0时,为驱动滑转率,当
s
x
<
0
{s_x} < 0
sx<0为制动滑移率,
C
x
{C_x}
Cx为纵向力曲线形状因子,表示纵向力曲线的形状;
D
x
{D_x}
Dx为纵向力峰值因子,表示最大纵向力值;
B
x
{B_x}
Bx为纵向力曲线刚度因子;
B
x
C
x
D
x
{B_x}{C_x}{D_x}
BxCxDx为纵向力零点处的纵向刚度;
E
x
{E_x}
Ex为纵向力曲率因子,表示曲线峰值附件的形状;
s
v
{s_v}
sv为纵向力曲线垂直偏移因子;
s
h
x
{s_{hx}}
shx为纵向力曲线水平偏移因子。
纵向力参数表达式为:
C
x
=
B
0
D
x
=
B
1
F
z
2
+
B
2
F
z
2
B
x
=
(
B
3
F
z
2
+
B
4
F
z
)
e
−
B
5
F
z
/
(
C
x
D
x
)
E
x
=
B
6
F
z
2
+
B
7
F
z
+
B
8
s
v
x
=
0
s
h
x
=
B
9
F
z
+
B
10
(
2
)
\begin{array}{c} {C_x} = {B_0}\\ {D_x} = {B_1}F_z^2 + {B_2}F_z^2\\ {B_x} = ({B_3}F_z^2 + {B_4}{F_z}){e^{ - {B_5}{F_z}}}/({C_x}{D_x})\\ {E_x} = {B_6}F_z^2 + {B_7}{F_z} + {B_8}\\ {s_{vx}} = 0\\ {s_{hx}} = {B_9}{F_z} + {B_{10}} \end{array}(2)
Cx=B0Dx=B1Fz2+B2Fz2Bx=(B3Fz2+B4Fz)e−B5Fz/(CxDx)Ex=B6Fz2+B7Fz+B8svx=0shx=B9Fz+B10(2)
纯侧偏工况下,Pacej]a89 轮胎侧向力表达式为:
F
y
0
=
D
y
sin
(
C
y
arctan
(
B
y
s
y
0
−
E
(
B
y
s
y
0
−
arctan
(
B
y
s
y
0
)
)
)
)
+
s
v
y
(
3
)
{F_{y0}} = {D_y}\sin ({C_y}\arctan ({B_y}{s_{y0}} - E({B_y}{s_{y0}} - \arctan ({B_y}{s_{y0}})))) + {s_{vy}}(3)
Fy0=Dysin(Cyarctan(Bysy0−E(Bysy0−arctan(Bysy0))))+svy(3)
式中,
α
=
α
+
s
h
y
\alpha = \alpha + {s_{h{\rm{y}}}}
α=α+shy为侧向力组合自变量: $\alpha $为轮胎侧偏角;
C
y
{C_y}
Cy为侧向力曲线形状因子,
表示侧向力曲线的形状;
B
y
{B_y}
By为侧向力曲线刚度因子:
D
y
{D_{\rm{y}}}
Dy为侧向力峰值因子,表示最大侧向
力值;
B
y
C
y
D
y
{B_y}{C_y}{D_y}
ByCyDy为侧向力零点处的侧向刚度;
E
y
{E_y}
Ey为侧向力曲率因子,表示曲线峰值附近的形
状;
s
v
y
{s_{vy}}
svy 为侧向力曲线垂直偏移因子;
s
h
y
{s_{h{\rm{y}}}}
shy为侧向力曲线水平偏移因子。
纵向力参数表达式为:
C
y
=
A
0
D
y
=
A
1
F
z
2
+
A
2
F
z
2
B
y
=
A
3
sin
(
2
arctan
F
z
/
A
4
)
(
1
−
A
5
∣
γ
∣
)
/
(
C
y
D
y
)
E
y
=
A
6
F
z
+
A
7
s
v
y
=
A
8
γ
+
A
9
F
z
+
A
10
s
h
y
=
A
11
F
z
γ
+
A
12
F
z
+
A
13
(
4
)
\begin{array}{c} {C_y} = {A_0}\\ {D_y} = {A_1}F_z^2 + {A_2}F_z^2\\ {B_y} = {A_3}\sin (2\arctan {F_z}/{A_4})(1 - {A_5}\left| \gamma \right|)/({C_{\rm{y}}}{D_y})\\ {E_y} = {A_6}{F_z} + {A_7}\\ {s_{vy}} = {A_8}\gamma + {A_9}{F_z} + {A_{10}}\\ {s_{hy}} = {A_{11}}{F_z}\gamma + {A_{12}}{F_z} + {A_{13}} \end{array}(4)
Cy=A0Dy=A1Fz2+A2Fz2By=A3sin(2arctanFz/A4)(1−A5∣γ∣)/(CyDy)Ey=A6Fz+A7svy=A8γ+A9Fz+A10shy=A11Fzγ+A12Fz+A13(4)
纯侧偏工况下,Pacejka89 轮胎回正力矩表达式为:
M
z
0
=
D
z
sin
(
C
z
arctan
(
B
z
α
z
−
E
z
(
B
z
α
z
−
arctan
(
B
z
α
z
)
)
)
)
+
s
v
z
(
5
)
{M_{z0}} = {D_z}\sin ({C_z}\arctan ({B_z}{\alpha _z} - {E_z}({B_z}{\alpha _z} - \arctan ({B_z}{\alpha _z})))) + {s_{vz}} (5)
Mz0=Dzsin(Czarctan(Bzαz−Ez(Bzαz−arctan(Bzαz))))+svz(5)
式中,
α
z
=
α
+
s
h
z
{\alpha _{\rm{z}}} = \alpha + {s_{hz}}
αz=α+shz为回正力矩组合自变量:
C
z
{C_z}
Cz为回正力矩曲线形状因子,表示回正力矩曲线的形状:
D
z
{D_z}
Dz 为回正力矩峰值因子,表示最大回正力矩值;
B
z
{B_z}
Bz 为回正力矩曲线刚度因子;
B
z
C
z
D
z
{B_z}{C_z}{D_z}
BzCzDz为回正力矩零点处的扭转刚度;
E
z
{E_z}
Ez为回正力矩曲率因子,表示曲线峰值附近的形状;
s
v
z
{s_{v{\rm{z}}}}
svz为回正力矩曲线垂直偏移因子;
s
h
z
{s_{hz}}
shz为回正力矩曲线水平偏移因子。
回正力矩参数表达式为:
C
z
=
C
0
D
z
=
C
1
F
z
2
+
C
2
F
z
B
z
=
(
C
3
F
z
2
+
C
4
F
z
)
(
1
−
C
6
∣
γ
∣
)
e
−
C
5
F
z
/
(
C
z
D
z
)
E
z
=
(
C
7
F
z
2
+
C
8
F
z
+
C
9
)
(
1
−
C
6
∣
γ
∣
)
s
v
z
=
C
11
γ
+
C
12
F
z
+
C
13
s
h
z
=
γ
(
C
14
F
z
2
+
C
15
F
)
+
C
16
F
z
+
C
17
(
6
)
\begin{array}{c} {C_z} = {C_0}\\ {D_{\rm{z}}} = {C_1}F_z^2 + {C_2}{F_z}\\ {B_z} = ({C_3}F_z^2 + {C_4}{F_z})(1 - {C_6}\left| \gamma \right|){e^{ - {C_5}{F_z}}}/({C_z}{D_z})\\ {E_z} = ({C_7}F_{\rm{z}}^2 + {C_8}{F_z} + {C_9})(1 - {C_6}\left| \gamma \right|)\\ {s_{vz}} = {C_{11}}\gamma + {C_{12}}{F_z} + {C_{13}}\\ {s_{hz}} = \gamma \left( {{C_{14}}F_z^2 + {C_{15}}F} \right) + {C_{16}}{F_z} + {C_{17}} \end{array} (6)
Cz=C0Dz=C1Fz2+C2FzBz=(C3Fz2+C4Fz)(1−C6∣γ∣)e−C5Fz/(CzDz)Ez=(C7Fz2+C8Fz+C9)(1−C6∣γ∣)svz=C11γ+C12Fz+C13shz=γ(C14Fz2+C15F)+C16Fz+C17(6)
二、仿真建模
1.simulink建模
根据上述公式,在simulink中建立模型如图(1)所示。
图(1)simulink模型
其中,子系统“纵向力参数”对应于公式(2),子系统“纵向力计算”对应公式(1);子系统“侧向力参数”对应于公式(4),子系统“侧向力计算”对应公式(3);子系统“回正力矩参数”对应于公式(6),子系统“回正力矩计算”对应公式(5)。
2.仿真结果
设置仿真时间为10s,滑移率和侧偏角输入模块采用斜坡函数,斜率为2,值域[-10,10]。采用ode45算法求解,将轮胎侧偏角α,轮胎侧向力Fz,轮胎回正力矩Mz等输入到工作空间后,采用附录中的matlab语句绘图,得到侧向力、回正力矩与侧偏角的关系曲线如图(2)(3)(4)所示:
图(2)轮胎纵向力-滑移率特性
从仿真结果来看,同一垂直载荷情况下,随着轮胎滑移率的增加,轮胎侧向力数值呈现线增大后减小,最后稳定在一定数值的趋势。
图(3)轮胎侧向力-侧偏角特性
根据仿真结果来看,随着垂直载荷的增加,同一轮胎侧偏角对应的侧向力也随之增加。从图(3)中可得当轮胎侧偏角增加大一定数值后,轮胎侧向力数值几乎不再增加
图(4)轮胎回正力矩-侧偏角特性
根据仿真结果来看,随着垂直载荷的增加,同一轮胎侧偏角对应的回正力矩也随之增加。
从图(4)中得到,同一垂直载荷的情况下,随着侧偏角的增大,回正力矩数值先增大后减小。
3.模型及绘图文件附录
simulink模型搭建的具体内容在本篇文章附件,可在matlab2015b及以上版本打开并运行,仅供参考。附录为绘图脚本。
legend('载荷3kN','载荷5kN','载荷8kN');
xlabel('轮胎滑移率(%)');
ylabel('轮胎纵向力(N)');
title('轮胎纵向力—滑动率关系曲线');
%%%%%%%%%%%%%%%%%%%%%
figure();plot(s,Fy,'Linewidth',1);
legend('载荷3kN','载荷5kN','载荷8kN');
xlabel('轮胎侧偏角(°)'); ylabel('轮胎侧向力(N)');
title('轮胎侧向力—侧偏角关系曲线');
%%%%%%%%%%%%%%%%%%%
figure();plot(s,Mz,'Linewidth',1);
xlabel('轮胎侧偏角(°)');
ylabel('轮胎回正力矩(N.m)');
title('轮胎回正力矩—侧偏角关系曲线');
legend('载荷3kN','载荷5kN','载荷8kN');
在这里插入代码片