☰
整车动力学仿真中的魔术轮胎模型:Simulink建模与参数辨识
2026/10/15 23:42:35 网站建设 项目流程

做整车动力学仿真的人,基本都绕不开轮胎模型。轮胎是车辆和地面唯一的接触体,纵向的驱动、制动,侧向的转弯、变道,最后都要靠那几块橡胶和地面之间的力来传递。而在工程里见得最多的,正是魔术轮胎模型,也就是常说的PAC轮胎模型。它用一条半经验的三角函数公式,就能把实验测出来的轮胎力曲线拟合得足够像,而且直接在MATLAB/Simulink 2021a里写成代码、搭成模型都很快。ABS/ESP控制算法开发、操纵稳定性仿真、底盘域控的MIL验证,这些场景里它都是默认选项之一。这篇文章就围绕这套模型在Simulink里的落地过程,把原理、代码、参数辨识和集成时的坑一次讲透。

1. 为什么车辆仿真绕不开一张“魔术曲线”:模型原理拆解

1.1 轮胎力计算的核心难点

轮胎模型难,不是因为它复杂,而是因为它所描述的物理过程高度非线性。侧偏角只有几度时,侧向力近似线性增长;再往大推,轮胎进入非线性区,力增长变缓,最后到达饱和;滑移率增加时,纵向力同样是先上升、过了峰值后反而下降。这种“先线性、后饱和、可能还会掉头向下”的曲线,用纯机理模型很难从头推导出来,因为橡胶材料的复杂特性、胎体结构、路面微观纹理全都参与其中。

所以工程上普遍采用半经验模型的思路:先确定一个带参数的经验公式,再用实验数据去拟合参数。魔术轮胎模型就是在这样的背景下出现的。它的核心做法是,用形如

y = D · sin(C · arctan(B·x - E·(B·x - arctan(B·x))))

的公式,把滑移率、侧偏角这类输入映射成轮胎力。这个式子乍看有点吓人,但本质就是正弦函数套了一个加强版的反正切函数。反正切负责让曲线从零点出发、逐渐趋于某个极限,正弦函数再把这个趋势收拢成峰值和饱和的形状。

正因为表达形式紧凑、参数含义相对清晰,这套模型在工程里落地非常方便。它不像有限元轮胎模型那样需要大量材料参数和网格,也不像查表模型那样只能依赖原始实验数据插值。你只要拥有几条不同载荷下的实验曲线,就能拟合出一套可用的参数,然后放心地丢给仿真环境使用。

1.2 B、C、D、E四个系数各自的职责

公式里的B、C、D、E是四个最关键的形状系数,很多刚接触的人一开始会混。我习惯用“曲线整容”来理解它:

  • D是峰值因子,直接决定曲线能拉多高,相当于曲线的“天花板”。对于侧向力曲线,D的大小基本对应路面附着系数乘垂直载荷。
  • C是形状因子,决定曲线整体是偏“方”还是偏“圆”,它控制着曲线在峰值附近以及峰值之后的走向形态。
  • B是刚度因子,主要影响零点附近的斜率。B越大,曲线起跳越快,小侧偏角下侧向力就越大,这对应轮胎的侧偏刚度。
  • E是曲率因子,调整峰值过后曲线的下降程度。E值不一样,峰值区域是平坦些还是急转直下,差别很明显,尤其对纵向力在打滑阶段的下降趋势影响很大。

另外公式里还有两个平移项Sh和Sv,分别对应横坐标偏置和纵坐标偏置。实际轮胎因为制造误差、胎压偏差等因素,零滑移/零侧偏角时并不完全出力为零,Sh和Sv就是用来补偿这种偏移的。

四个系数不能孤立地调。B、C、D三个共同决定零点斜率,近似关系是B·C·D。我一般在初估参数时,就用这个乘积去对照实验曲线的初始斜率,快速把B的范围定下来,再交给优化算法去精调。

1.3 纵向力、侧向力和回正力矩的公式差异

魔术轮胎模型并不仅仅算一个方向上的力,纵向力Fx、侧向力Fy、回正力矩Mz都有各自的公式。平时用得最多的是前两个。

纵向力的输入是纵向滑移率κ,工程里常用定义是κ = (ω·R - vx) / vx。注意不同文献可能分子分母有不同约定,驱动和制动时符号也不完全一致,这一点非常容易埋雷。带入公式后得到Fx,曲线形状是:滑移率从0开始增加,纵向力快速上升,到达峰值后缓慢下降,最终趋近于纯滑动状态下的摩擦力。

侧向力的输入是侧偏角α,公式形式完全一样,只是自变量从κ换成了α。侧偏角增大时,侧向力先近似线性增长,然后进入非线性区,最后饱和。滑动状态下的侧向力会明显下降,这也是车辆在极限工况下失去稳定性的根源。

回正力矩Mz可以独立用一套公式去拟合,也可以偷懒用“侧向力 × 气胎拖距”来近似。前者精度高但参数多,后者简单,适合只关心车辆运动而不过多关注转向手感的场景。我在整车级动稳定控制仿真里,大部分时候会先忽略Mz,等到做转向手感或方向盘力矩相关课题再单独加上。

这三个量虽然公式长得像,但参数完全是独立标定的,不能把纵向力那套B、C、D直接搬到侧向力上。拟合实验数据时要分开处理,分别保存独立的参数集。

2. 在Simulink里把公式变成可跑模型:接口、代码与封装

2.1 输入输出接口怎么定

很多人一上来就写代码,结果到最后发现接口不统一,换个工况就要重新接线。我建议在写第一行公式之前,先把输入输出定义清楚。

对于纯纵向/纯侧向分开计算的基础版本,输入一般有四个:侧偏角α、纵向滑移率κ、垂直载荷Fz、路面附着系数μ。输出就是Fx和Fy。其中α单位用rad,κ无量纲,Fz单位用N,μ根据工况取值比如干燥沥青0.9、湿滑路面0.4。全部统一到国际单位制里,后面会少很多麻烦。

Fz要不要作为输入,取决于你需不需要模拟载荷转移。如果只是做单胎分析,可以固定Fz,直接作为参数常量写在参数结构体里;如果做整车模型,车辆加减速和侧倾会引起四个车轮载荷变化,Fz就必须动态输入。我强烈建议从开始就把Fz设计成输入,这样模型扩展性会好很多,后面接入整车环境不用返工。

参数统一放在一个struct结构体里,比如P.Bx、P.Cx、P.Dx这种命名方式。不要散落成十几个单独变量,尤其是当你在Simulink里做批量参数扫描时,结构体的优势非常明显。

2.2 一份可直接套用的核心MATLAB函数

在Simulink 2021a的MATLAB Function块里,核心代码大致可以写成下面这样:

function [Fx, Fy] = pacTyre(alpha, kappa, Fz, P) % P是参数结构体,包含纵向和侧向两套系数 % 纵向力 xB = kappa + P.Shx; Fx = P.Dx * sin(P.Cx * atan(P.Bx * xB - P.Ex * (P.Bx * xB - atan(P.Bx * xB)))) + P.Svx; % 侧向力 yB = alpha + P.Shy; Fy = P.Dy * sin(P.Cy * atan(P.By * yB - P.Ey * (P.By * yB - atan(P.By * yB)))) + P.Svy; % 如果需要按垂直载荷对峰值进行缩放 Fx = Fx * (Fz / P.FzRef); Fy = Fy * (Fz / P.FzRef); end

这里默认Fy和Fx分别按独立公式计算,也就是纯工况。最后两行是我自己习惯加的简单载荷缩放。注意,这只是一个工程近似。严格来说,魔术轮胎模型的参数本身就会随Fz变化,最准确的做法是把B、D、C都拟合成Fz的相关函数,或者在不同Fz下插值。但如果你刚开始接触这个模型,用D随Fz线性缩放,其它参数固定,往往也够在基础仿真里用了。

以及,公式里的自变量都是弧度制的角度,输入alpha千万不能直接填度数,否则曲线会完全变形。Simulink里如果是角度传感器给出的信号,一定记得先乘pi/180。

2.3 Simulink子系统搭建步骤

MATLAB Function块只是核心,要让它真正在Simulink模型里跑起来,还需要完整接线和配置。我习惯的搭建步骤是:

  1. 打开MATLAB 2021a,在命令行输入simulink,新建空白模型。
  2. 在Library Browser中找到User-Defined Functions,拖入MATLAB Function块。
  3. 双击MATLAB Function块,进入编辑器,把上面的函数粘贴进去。
  4. 再拖入若干Inport(对应alpha、kappa、Fz)和Outport(对应Fx、Fy),按顺序连线。
  5. 在MATLAB Function块的编辑界面里,用Edit Data按钮把P的参数维度、类型定义好。P是结构体时,建议直接在工作区里用脚本定义一个P实例,并勾选参数类型继承,让代码生成器自动推断。
  6. 单独配置参数结构体:写一个MATLAB脚本,在仿真前运行,把P写入工作区。
P.Bx = 2.5; P.Cx = 1.8; P.Dx = 9200; P.Ex = 0.97; P.Shx = 0; P.Svx = 0; P.By = 5.5; P.Cy = 1.6; P.Dy = 7600; P.Ey = 0.95; P.Shy = 0; P.Svy = 0; P.FzRef = 4200;
  1. 配置模型求解器:在Model Settings里选择固定步长,求解器选ode4,步长我常用1e-4秒。轮胎模型的非线性很强,步长太大会在集成到整车模型时出现数值振荡。
  2. 保存模型,先用Simulink自带的Signal Editor或者Ramp信号做开环扫掠测试,看输出曲线是否符合预期。

模型封装时,我建议把参数结构体P定义成一个Simulink.Parameter对象,或者直接放在基础工作区。这样后面做批处理仿真、参数扫描都方便,不用每次打开模型手动改数值。

2.4 示例参数与开环曲线验证

为了验证模型代码有没有写对,我们可以在MATLAB脚本里先把曲线画出来,不急着进Simulink。这样定位问题更快速。

kappa = -0.3:0.005:0.3; alpha = -0.2:0.005:0.2; [Fx, Fy] = pacTyre(0, kappa, 4200, P); figure; plot(kappa, Fx); xlabel('kappa'); ylabel('Fx'); figure; plot(alpha, Fy); xlabel('alpha'); ylabel('Fy');

用上面那套示例参数,你会发现纵向力在滑移率大约0.1~0.18之间达到峰值,之后缓慢下降;侧向力在小侧偏角范围内基本线性,到0.12~0.15 rad左右开始饱和。这个形态和真实轮胎特性是对得上的。

需要强调,这套参数是我随手给出的示意值,不是某个真胎的实验参数。实际工程中参数必须从实验数据或者目标车型的轮胎手册里辨识得到,不能指望网上找个参数表直接开跑。但用这种示例参数做模型搭建、流程验证,已经完全足够了。

3. 参数从哪来:辨识流程与验收标准

3.1 从纯工况到组合工况的取舍

魔术轮胎模型在工程里分成几个层次。第一层是纯工况,也就是单纯改变滑移率、固定侧偏角为0,去拟合纵向力;或者单纯改变侧偏角、固定滑移率为0,去拟合侧向力。这个层次参数少,实现简单,适合先跑通流程。

第二层是组合工况,也就是侧偏角不为0的同时还有滑移,此时纵向力和侧向力会互相耦合、互相削弱。一个常见的处理方式是用权重函数合成:

Fx = Fx0 · Gxa(α) Fy = Fy0 · Gyk(κ)

其中Fx0和Fy0分别是纯工况下的纵向力和侧向力,Gxa和Gyk是取值范围0到1的加权函数。这样实现起来不复杂,又能明显提升极限工况下的仿真可信度。

具体G函数的形式,建议从实验数据里拟合,而不是强行套固定公式。真实轮胎在较大侧偏角时还能提供一定纵向力,但如果侧偏角已经接近饱和,纵向力的响应会显著变弱,这个趋势必须体现在权重函数里。

如果目标是做标准道路工况的操纵稳定性仿真,纯工况往往也够;但如果要模拟极限避障、冰雪路面起步这类大滑移+大侧偏同时出现的情况,组合工况基本是必需品。

3.2 最小二乘辨识的实操步骤

参数辨识听起来很学术,做起来其实就是调整参数,让模型曲线尽可能贴近实验曲线。我用MATLAB完成这个任务已经很多次,流程基本稳定:

  1. 把实验数据整理成三列:输入x(滑移率或侧偏角)、垂直载荷Fz、输出力F。
  2. 按Fz分组建档案,每组数据单独拟合,避免让不同载荷下的特性互相干扰。
  3. 设初始参数。这个初始值很关键,随便给可能让优化器陷入局部极小。经验做法:
    • D先用该组数据里的最大力测量值来初始化;
    • C先给1.5~2之间的常见值;
    • B用零点附近的斜率除以C和D的乘积来反算;
    • E先给0.9~1左右。
  4. 调用lsqcurvefit或lsqnonlin,目标函数是计算曲线与实验曲线的残差平方和。设置边界约束,比如C保持在0.8~2.5,防止优化器给出不合理的形状。
  5. 拟合完成后,把Sh、Sv一起放开再优化一遍,使曲线在零输入附近也能对齐。

在这个过程里,我最常遇到的问题就是参数强相关。比如B和C同时变化,可能在某个范围内产生相似的曲线形状,导致拟合结果不稳定。解决办法是固定其中一个,只优化另一个;或者用多组Fz下的数据做联合约束,减少自由度。

3.3 验收指标与参数外推经验

拟合结束不代表工作完成,还要有明确的验收标准。我常用两个指标:最大绝对误差和均方根误差RMSE,并且统一折算成相对垂直载荷的百分比。对于一个还算靠谱的拟合,峰值附近的误差控制在5%以内是比较理想的,曲线过零点附近的误差允许稍微大一点,因为那里的力值本来就小,相对误差容易偏大。

还要做交叉验证。拿一组不参与拟合的数据丢进模型里跑一遍,看预测效果如何。如果只拟合了一组载荷下的数据,模型在其他载荷下完全失真,说明你手里数据太少,或者D随Fz的变化规律没有表达出来。

参数外推也是工程里绕不开的话题。轮胎实验数据往往只覆盖一定范围的Fz和侧偏角,模型用起来却可能遇到超出范围的情况。我的建议是:D和B尽量建立成Fz的线性或二次多项式函数;超出数据范围时,宁可让模型保守地趋于饱和,也不要外推出一个奇怪的大斜率。车辆动力学仿真里,一个不合理的轮胎力往往比一个偏保守的轮胎力更容易把结果带偏。

4. 集成到整车系统后出现的四个典型坑

4.1 零速附近的滑移率分母爆炸

滑移率的计算公式里带有车轮速度或车辆速度作为分母,车速接近0时,分母趋近于0,滑移率会瞬间变成很大的数,甚至正负震荡。如果不处理,轮胎模型会输出一个离谱的大力,直接把整车模型震飞。

我的习惯是在计算滑移率之前先做速度门槛限制:

vxSafe = max(abs(vx), 0.1); % 0.1 m/s以下按0.1处理 kappa = (omega * Rw - vx) / vxSafe;

当车速非常低时,再单独判断:如果vx小于某个阈值,比如0.5m/s,直接把轮胎力过渡到零或者按特殊低速模型处理。为什么不用eps?因为eps在数值上虽然不会除零,但轮胎力依然会异常放大,必须设置一个物理上合理的门槛值。

这个坑几乎每个刚开始做整车仿真的人都会踩,而且一旦不处理,现象非常明显:车辆静止起步那一瞬间,纵向力会跳到一个不可理喻的数值,整个车直接往后弹。

4.2 轮胎力与车辆方程之间会形成代数环

整车模型里,轮胎力影响整车加速度,加速度积分出速度,速度又反过来算滑移率,滑移率再决定轮胎力。这个回路在同一时刻内存在,Simulink会检测到代代数环,然后通过迭代求解来处理。迭代本身没问题,但它会拖慢仿真速度,严重情况下还会引起数值振荡。

解决思路通常是允许一层延迟:在轮胎力输出端或者滑移率计算端加一个Unit Delay块,让信号在时间上错开一个步长,从代数环变成递推关系。这样做会带来一个步长的相位滞后,但对绝大多数控制算法开发和车辆动力学分析来说,1e-4秒的步长滞后几乎可以忽略。

如果你用的是变步长求解器,代数环有时能迭代过去,但会在某些时刻突然卡得非常慢。所以我一直坚持用定步长加Unit Delay的组合,宁可多一个步长延迟,也要让仿真行为可控、可复现。

4.3 坐标符号反了,曲线变成“狗啃状”

魔术轮胎模型原始实验数据和公式约定大多是在特定坐标系下定义的,而整车模型里可能采用ISO坐标系或者SAE坐标系。这两种坐标系下的侧偏角符号恰好相反,如果搞反了,侧向力方向和实际转向方向会不对应。

最直观的现象是:给一个正的侧偏角,按理说应产生正的侧向力,结果出来的力却反向,整车模型会表现为“转向方向与横摆方向相反”,曲线画出来像左右不对称的狗啃状。

我建议在Simulink模型入口处加一个明确的坐标转换模块,把所有方向的定义统一成模型内部约定。不要依赖建模者“记着同一个方向”,那太脆弱了。接口注释里写清楚:alpha为正时,表示车轮前进方向偏左还是偏右,以及对应的Fy方向,让后来接手的人一眼就能看懂。

4.4 单位混乱导致的参数数量级错误

轮胎实验数据经常以kN、km/h、kPa为单位,而模型内部如果用N、m/s、rad,就需要换算。最容易出问题的是B系数,它不是纯量,它的量纲取决于自变量的量纲和C、E的配合。当你把滑移率从“百分比”换成0到1的小数时,B的值也要跟着缩放,否则曲线会被压扁或者炸开。

解决方法是建立一套固定的内部单位制,并且写一个单位转换脚本。所有外部数据进模型之前,先过一遍转换函数,统一输出到SI单位。这个脚本要留在模型文件旁边,方便日后追溯。

我在实际使用中还发现,有人会把Fz单位写成kN,D还是按N来标定,结果峰值力瞬间放大1000倍。这种错误跑出来的曲线非常夸张,一看就知道不对,但有些人会先去怀疑代码,而不是查单位。

5. 从ABS到状态估计:真实工程里的扩展用法

5.1 ABS开发中如何用附着曲线

在ABS控制算法开发里,魔术轮胎模型的主要价值不只是“让车能跑”,而是提供一条完整的滑移率-附着系数曲线。控制策略要在这条曲线上找到峰值附着对应的滑移率,然后把车轮滑移率控制在峰值点左侧的安全区域。

用Simulink做ABS的模型在环仿真时,我会把PAC轮胎模型作为被控对象模型里的核心模块,路面附着系数设成可切换的参数,这样就能模拟高附着、低附着、对开路面的切换,观察ABS算法在不同路面上能否自适应调节压力。

这里面有个细节:ABS仿真相比普通操纵稳定性仿真,更关心滑移率在0到0.3之间那段曲线的形态。峰值点的位置和峰值后的下降斜率,直接决定了控制算法的动作时机和动作强度。因此拟合参数时,要重点保证峰值附近的拟合精度,尤其是低附着路面下峰值滑移率往往会更小,训练数据必须覆盖足够低的路面附着。

5.2 车辆状态估计中的非线性映射

另一个和我常打交道的场景是车辆状态估计。轮胎力计算公式可以作为卡尔曼滤波或扩展卡尔曼滤波的观测方程,通过测量到的加速度和轮速,反向估计车辆的质心侧偏角、横摆角速度甚至路面附着系数。

这个过程里,魔术轮胎模型的价值在于它给出了一种可导的非线性映射。你能够求出轮胎力对侧偏角、滑移率的偏导数,用于计算雅可比矩阵。换成纯查表模型,虽然也能工作,但梯度信息就要靠差分去近似,误差和计算量都不理想。

当然,这种基于模型的估计器对参数准确性要求很高。一个常见做法是把路面附着系数作为待估计状态,在滤波过程中实时调整D的取值。模型里的D不再是固定值,而是随估计器输出变化,相当于让模型自己去适应路面。

5.3 我建议的落地顺序与几个习惯

最后分享几点这些年积累的实操习惯,希望能帮你少走弯路。

第一,永远先画曲线再进整车系统。无论代码改了什么,打开脚本画一条Fx-kappa曲线、Fy-alpha曲线,看一眼形状再继续往下做。很多错误在单曲线阶段就能发现,不需要跑到整车模型里去看整车轨迹慢慢排查。

第二,参数文件和工作区分离。我会把一个模型的参数脚本放在固定文件里,每次仿真前运行。同一套代码,只要改脚本里的数值就能换成不同轮胎,接整车模型时也方便做多组参数对比。

第三,对参数初值保持敏感。辨识过程中如果初值给得离谱,优化结果经常充满“人工感”,曲线虽然拟合得不错,参数之间却互相补偿。这就好比用两根斜率不一样的线去逼近一条弯曲的线,中间虽然贴合,端点和外推区完全不可信。所以不要盲目点头交给算法,结果出来要再用人眼判断一遍。

第四,越简单越好。如果模型只需要体现ABS控制逻辑的正确性,那一个纵向力模型加一个简单垂向载荷分配就够了,完全没必要堆上完整组合工况和回正力矩。模型复杂度要和你的研究目标匹配,而不是和“学术完整度”匹配。

第五,代码生成部署前的检查。如果你准备把Simulink模型用Embedded Coder生成C代码嵌入到硬件里,注意三角函数在嵌入式芯片上的开销。必要时可以把预先计算好的轮胎力曲线生成查表数据,在代码里用插值替代实时三角函数,代价是精度略微下降,但实时性提升明显。这个优化要在功能验证完成之后再做,不要一开始就陷入性能优化。

需要专业的网站建设服务?

联系我们获取免费的网站建设咨询和方案报价,让我们帮助您实现业务目标

立即咨询