4. 多体动力学建模理论
齿轮传动系统的多体动力学仿真，本质上不是“把齿轮画出来让它转”，而是把一对或多对齿轮、轴、轴承、箱体等离散为若干刚体或柔性体，用统一的数学语言描述它们的位姿、约束与受力，再通过数值积分得到随时间演化的运动与载荷。本节围绕多刚体系统建模的基本理论展开，重点说明广义坐标如何选取、约束方程如何建立、拉格朗日方程如何组装，以及从模型到求解的完整流程。
一、多刚体系统的基本描述
多刚体系统由若干刚体通过运动副连接而成。每个刚体在空间中的位置和姿态需要一组独立参数描述。对于平面问题，一个刚体常用三个参数：质心坐标 (x, y) 和转角 θ；对于空间问题，则常用质心坐标 (x, y, z) 加姿态参数（如欧拉角或四元数）。在齿轮传动中，若只关心传动误差和动态啮合力，常将轴系简化为“绕定轴转动”的集中惯量，此时每个齿轮只需一个转角参数即可；若需评估箱体振动或轴承座受力，则必须保留空间六自由度。
建模时首先要明确“哪些构件必须保留为独立刚体”。以一级直齿轮副为例，主动轮、从动轮、两根轴、四个轴承座通常需要分别建模；若箱体刚度远大于轴系，可将其视为刚体或直接作为参考系。每个刚体的质量、质心位置、转动惯量张量是后续动力学方程的基本参数，必须与三维模型或图纸一致。
二、广义坐标与约束方程
广义坐标是描述系统位形的最少一组独立参数，记作 q = [q₁, q₂, …, qₙ]ᵀ。对于含 n 个刚体的平面系统，若每个刚体取 3 个笛卡尔坐标，则共有 3n 个坐标，但其中许多并不独立，需要通过约束方程消去。约束方程的一般形式为：
Φ(q, t) = 0 （完整约束）
Φ(q, q̇, t) = 0 （非完整约束，齿轮传动中较少用）齿轮传动中最典型的约束是齿轮副约束。设主动轮转角为 θ₁、从动轮转角为 θ₂，传动比为 i = z₂/z₁，则理想啮合约束可写为 θ₂ − θ₁/i = 0（或 r₁θ₁ + r₂θ₂ = 0，取决于正方向定义）。实际仿真中，为模拟齿侧间隙和啮合刚度，常把该约束“软化”为弹簧-阻尼力元，而不是硬约束。此外，轴承约束通常将轴颈中心与轴承座中心用径向弹簧-阻尼连接，允许微小相对位移；若采用理想铰接，则约束为轴颈中心与轴承中心重合。
约束方程建立后，系统自由度从 3n 降为 3n − m（m 为独立约束数）。这一步直接决定后续方程规模，也影响求解效率。实际建模中应尽量用“力元”替代“硬约束”，因为硬约束会引入拉格朗日乘子，增加方程刚性。
三、拉格朗日方程与动力学组装
多刚体系统的动力学方程可由拉格朗日方程导出。对于含约束的系统，常用第一类拉格朗日方程：
d/dt(∂T/∂q̇) − ∂T/∂q + Φ_qᵀλ = Q
Φ(q, t) = 0其中 T 为系统动能，Q 为广义力，Φ_q 为约束雅可比矩阵，λ 为拉格朗日乘子。将动能写成 T = ½ q̇ᵀM(q)q̇，代入后可得矩阵形式：
M(q)q̈ + Φ_qᵀλ = Q(q, q̇, t)
Φ(q, t) = 0这就是多刚体系统动力学求解的标准微分-代数方程组（DAE）。在齿轮传动中，M 通常为对角或分块对角矩阵（各刚体惯量），Q 包含啮合力、轴承力、驱动力矩和负载力矩。啮合力本身又依赖于相对位移和相对速度，因此 Q 是 q 和 q̇ 的函数，使方程具有非线性。
实际组装时，可按“体—铰—力元”三层结构组织：先为每个刚体写出其惯量矩阵，再根据运动副写出约束方程和雅可比，最后把各力元按作用点映射到广义力向量。对于齿轮副，若采用硬约束，则 Φ 中增加一行 θ₂ − θ₁/i；若采用弹簧-阻尼力元，则 Q 中增加 −kₘ(r₁θ₁ + r₂θ₂) 和 −cₘ(r₁θ̇₁ + r₂θ̇₂) 对应的广义力项。
四、动力学求解流程与数值实现
从模型到结果，典型流程如下：
- 定义坐标系与单位：统一采用国际单位制，明确各转角正方向（如逆时针为正），并在全文和模型中保持一致。
- 确定刚体参数：从三维模型或图纸获取质量、质心、转动惯量。示例：主动轮 m₁ = 2.5 kg，J₁ = 0.012 kg·m²；从动轮 m₂ = 6.8 kg，J₂ = 0.095 kg·m²。
- 建立约束与力元：齿轮副取传动比 i = 2.5，啮合刚度取时变值 kₘ(t) = k₀ + k₁cos(ωₘt)，其中 k₀ = 2×10⁸ N/m，k₁ = 0.3k₀，ωₘ 为啮合频率。
- 组装方程：形成 M、Φ_q、Q 的数值矩阵。对于 n 个刚体、m 个约束，方程规模为 (n+m) 个未知量。
- 选择积分器：刚性较强的齿轮系统推荐 Newmark-β 法或隐式 Runge-Kutta 法。时间步长建议取啮合周期的 1/20～1/50。示例：啮合频率 500 Hz，则步长取 1×10⁻⁴ s 左右。
- 施加初始条件：通常从静态平衡位置出发，给主动轮一个初始转速 ω₁ = 100 rad/s，从动轮按传动比匹配。
- 求解与后处理：提取 θ₁(t)、θ₂(t)、啮合力 Fₙ(t)、轴承力 F_b(t)，计算传动误差 TE = θ₂/i − θ₁（或按定义取相反数）。
五、建模中的常见问题与处理
第一，约束过刚导致求解困难。理想齿轮副硬约束会使方程指数增大，实际中建议用大刚度弹簧近似，刚度取 10⁸～10⁹ N/m 量级，并配合适当阻尼（如 cₘ = 0.01kₘ/ωₘ）。第二，时变啮合刚度处理。直齿轮可取矩形波或傅里叶级数近似，斜齿轮需考虑接触线长度变化。第三，轴承建模。滚动轴承可用径向弹簧-阻尼，滑动轴承需引入油膜力，后者与转速和间隙强相关。第四，初始平衡。建议先做静力学分析，得到各刚体初始位移，再以此为初值进行动力学积分，避免初始瞬态掩盖真实响应。
多刚体系统建模理论为齿轮传动动力学提供了统一的数学框架。掌握广义坐标选取、约束方程建立、拉格朗日方程组装和数值求解流程，是后续进行齿轮系统动态特性分析、故障机理研究和参数优化的基础。实际建模中应始终关注“模型规模与精度”的平衡，用合理的简化换取可接受的求解效率，同时通过能量检查和步长收敛性验证保证结果可信。