用Python和SymPy搞定汽车二自由度模型：从理论方程到仿真代码（保姆级推导）

用Python和SymPy搞定汽车二自由度模型：从理论方程到仿真代码（保姆级推导）

在车辆动力学研究中，二自由度模型是理解汽车侧向和横摆运动的基础。许多教科书会给出微分方程的最终形式，却很少展示如何从原始公式推导出可执行的代码。这正是本文要解决的问题——我们将用Python的SymPy库，像解数学题一样拆解每个推导步骤，最终生成可直接用于仿真的状态空间方程。

1. 理解二自由度模型的物理意义

汽车运动本质上包含六个自由度，但在特定假设下可简化为两个关键自由度：

• 侧向运动（y轴）：车辆横向移动的动力学特性
• 横摆运动（z轴旋转）：车辆绕垂直轴的旋转运动

为什么选择这两个自由度？通过以下假设实现简化：

1. 忽略悬架影响：假设车身是刚性体，排除垂直运动和侧倾/俯仰
2. 恒定前进速度：x轴方向加速度为零，无需考虑纵向动力学
3. 线性轮胎特性：侧向加速度≤0.4g时，轮胎力与滑移角呈线性关系

# 导入基础库 import numpy as np import sympy as sp from IPython.display import display, Math # 定义符号变量 v, u, r = sp.symbols('v u r') # 侧向速度、纵向速度、横摆角速度 delta, a, b = sp.symbols('delta a b') # 前轮转角、前轴距、后轴距 m, I_z = sp.symbols('m I_z') # 质量、绕z轴转动惯量 C_f, C_r = sp.symbols('C_f C_r') # 前/后轮侧偏刚度

2. 运动学方程的符号推导

运动学分析的核心是建立速度-加速度关系。我们采用车辆坐标系（车身固定坐标系），考虑非惯性系下的加速度分量：

$$ a_y = \dot{v} + u r $$

这个看似简单的公式实际包含两个物理效应：

• $\dot{v}$：侧向速度变化率
• $u r$：由于横摆运动产生的向心加速度

# 运动学方程推导 v_dot = sp.symbols('\dot{v}') a_y = v_dot + u*r # 侧向加速度 display(Math(f"a_y = {sp.latex(a_y)}"))

3. 动力学方程的建立与展开

根据牛顿第二定律和力矩平衡原理，建立动力学方程：

侧向力平衡：$$ m a_y = F_{yf} \cos\delta + F_{yr} $$

横摆力矩平衡：$$ I_z \dot{r} = a F_{yf} \cos\delta - b F_{yr} $$

在小角度假设下（$\cos\delta \approx 1$），轮胎侧向力可表示为： $$ F_{yf} = -C_f \alpha_f, \quad F_{yr} = -C_r \alpha_r $$

其中侧偏角计算公式为： $$ \alpha_f = \delta - \frac{v + a r}{u}, \quad \alpha_r = -\frac{v - b r}{u} $$

# 定义轮胎力和侧偏角 alpha_f = delta - (v + a*r)/u alpha_r = -(v - b*r)/u F_yf = -C_f * alpha_f F_yr = -C_r * alpha_r # 动力学方程 r_dot = sp.symbols('\dot{r}') eq1 = sp.Eq(m*a_y, F_yf + F_yr) # 简化后的侧向力平衡 eq2 = sp.Eq(I_z*r_dot, a*F_yf - b*F_yr) # 横摆力矩平衡 display(Math(f"{sp.latex(eq1)} \\\\ {sp.latex(eq2)}"))

4. 状态空间方程的矩阵形式推导

将方程整理为状态空间形式 $\dot{X} = A X + B U$：

1. 状态变量 $X = [v, r]^T$
2. 输入 $U = \delta$
3. 输出通常选择为横摆角速度 $r$ 和侧向加速度 $a_y$

通过SymPy求解代数方程：

# 解方程组得到状态导数 solutions = sp.solve([eq1.subs(a_y, v_dot + u*r), eq2], [v_dot, r_dot]) # 提取结果 v_dot_expr = solutions[v_dot].simplify() r_dot_expr = solutions[r_dot].simplify() # 转换为矩阵形式 A = sp.Matrix([ [v_dot_expr.coeff(v), v_dot_expr.coeff(r)], [r_dot_expr.coeff(v), r_dot_expr.coeff(r)] ]) B = sp.Matrix([ [v_dot_expr.coeff(delta)], [r_dot_expr.coeff(delta)] ]) display(Math(f"A = {sp.latex(A)} \\\\ B = {sp.latex(B)}"))

典型参数下的数值矩阵示例：

# 数值化矩阵示例 params = { m: 1500, I_z: 2500, a: 1.2, b: 1.5, C_f: 80000, C_r: 100000, u: 20 } A_num = A.subs(params).evalf() B_num = B.subs(params).evalf() print("数值化A矩阵:\n", np.array(A_num, dtype=float)) print("数值化B矩阵:\n", np.array(B_num, dtype=float))

5. 时域仿真与结果可视化

使用Scipy的数值积分器进行时域仿真：

from scipy.integrate import solve_ivp import matplotlib.pyplot as plt def vehicle_model(t, x, A, B, delta_func): delta = delta_func(t) return A @ x + B * delta # 阶跃转向输入 def delta_step(t): return 0.1 if t >= 1 else 0 # 仿真参数 t_span = (0, 5) x0 = [0, 0] # 初始状态 # 运行仿真 sol = solve_ivp(vehicle_model, t_span, x0, args=(A_num, B_num, delta_step), dense_output=True) # 可视化结果 t = np.linspace(0, 5, 500) v, r = sol.sol(t) plt.figure(figsize=(12, 4)) plt.subplot(121) plt.plot(t, v, label='侧向速度 [m/s]') plt.legend() plt.subplot(122) plt.plot(t, np.rad2deg(r), label='横摆角速度 [deg/s]') plt.legend() plt.tight_layout()

常见问题排查指南：

1. 符号正负混淆：

   • 检查轮胎力公式中的负号是否与坐标系定义一致
   • 确认侧偏角计算时分子项的加减组合

2. 数值不稳定：

   • 确保纵向速度u不为零
   • 检查参数单位是否统一（如角度用弧度）

3. 物理合理性验证：

   • 阶跃响应中横摆角速度应趋于稳态值
   • 侧向速度变化率应与理论预期相符

6. 模型扩展与实际应用

基础模型可通过以下方式增强实用性：

速度自适应模型：

def create_speed_dependent_model(u_val): params = {m: 1500, I_z: 2500, a: 1.2, b: 1.5, C_f: 80000, C_r: 100000, u: u_val} return A.subs(params).evalf(), B.subs(params).evalf()

典型应用场景：

• 自动驾驶控制算法开发
• 车辆稳定性控制策略验证
• 转向特性灵敏度分析

参数影响规律总结：

在完成所有推导后，建议将完整代码封装为类：

class TwoDOFVehicle: def __init__(self, params): self.params = params self._compute_system_matrices() def _compute_system_matrices(self): # 实现矩阵计算逻辑 pass def simulate(self, delta_func, t_span, x0): # 实现仿真逻辑 pass

实际项目中，这种模型通常作为更复杂仿真系统的子系统。例如在自动驾驶仿真中，二自由度模型可以快速验证控制算法的基础性能，然后再使用更高阶模型进行精细验证。