4、动力学建模基础:牛顿-欧拉方程、刚体动力学、四旋翼的力和力矩来源

四旋翼飞行器的动力学模型是姿态控制、轨迹规划与参数辨识的数学根基。工程实践中,无论是PID调参、LQR设计还是基于模型的无迹卡尔曼滤波(UKF)辨识,都需要一个能够准确描述机体在力和力矩作用下运动规律的微分方程组。本章从牛顿-欧拉方程出发,建立刚体六自由度动力学模型,并深入剖析螺旋桨拉力、反扭矩与陀螺效应这三个核心力和力矩来源的物理本质与数学表达。

四旋翼动力学建模 — 力/力矩来源与牛顿-欧拉方程 机体坐标系 (B系) X_B (前) Z_B (下) Y_B (右) O_B 惯性坐标系 (E系) X_E (北) Z_E (地) 旋转矩阵 R 螺旋桨力/力矩 电机1 (前左) ω₁ 顺时针 电机2 (前右) ω₂ 逆时针 电机3 (后左) ω₃ 逆时针 电机4 (后右) ω₄ 顺时针 Tᵢ = C_T · ρ · A · (ωᵢ·R)² 力矩来源 反扭矩: Qᵢ = C_Q · ρ · A · ωᵢ²·R³ 方向与桨旋转方向相反 陀螺力矩: M_gyro = J_p · (ω × Ω) 高速旋转桨对机体角运动的耦合 牛顿-欧拉方程 平动 (牛顿第二定律) m·V̇ + ω×(m·V) = F F = [0, 0, -mg] + R·F_rotor F_rotor = [0, 0, ΣTᵢ] 转动 (欧拉方程) J·ω̇ + ω×(J·ω) = M M = M_roll + M_pitch + M_yaw + M_gyro + M_disturb J = diag(Ixx, Iyy, Izz) 状态向量 x = [x, y, z, φ, θ, ψ, u, v, w, p, q, r] 12维状态空间 C_T: 拉力系数 | C_Q: 扭矩系数 | ρ: 空气密度 | A: 桨盘面积 | R: 桨半径 | J_p: 螺旋桨转动惯量

4.1 刚体动力学与坐标系定义

四旋翼的动力学建模首先需要明确两个坐标系:惯性坐标系(E系)和机体坐标系(B系)。惯性系固定于地面,通常取北-东-地(NED)或东北-天(ENU)方向;机体系固连于飞行器质心,其X轴指向机头方向,Y轴指向右侧,Z轴垂直向下(右手定则)。两个坐标系之间的转换通过旋转矩阵R(由欧拉角φ、θ、ψ构成)完成,该矩阵在姿态解算和控制器设计中反复使用。

刚体动力学模型基于两个基本假设:机体为刚性体(忽略弹性形变)和质量分布恒定(忽略燃料消耗或载荷变化)。在此假设下,飞行器的运动可分解为质心的平动和绕质心的转动,分别由牛顿第二定律和欧拉方程描述。平动方程中,外力包括重力、螺旋桨总拉力以及空气阻力;转动方程中,外力矩包括滚转力矩、俯仰力矩、偏航力矩以及陀螺效应力矩。

核心公式:牛顿-欧拉方程的矢量形式为:

平动: m·(V̇ + ω×V) = F_gravity + F_rotor + F_aero
转动: J·ω̇ + ω×(J·ω) = M_rotor + M_gyro + M_aero

其中V为线速度向量,ω为角速度向量,J为3×3惯性张量矩阵。对于对称布局的四旋翼,J可近似为对角阵diag(Ixx, Iyy, Izz)。

4.2 螺旋桨拉力与反扭矩模型

螺旋桨是四旋翼唯一的升力与推力来源。根据动量理论和叶素理论,单个螺旋桨产生的拉力T与转速ω的平方成正比:T = C_T · ρ · A · (ω·R)²,其中C_T为拉力系数(无量纲),ρ为空气密度,A为桨盘面积(πR²),R为桨半径。在工程简化模型中,常将常数项合并为拉力系数k_T,即T = k_T · ω²。同理,螺旋桨旋转时受到的空气阻力产生反扭矩Q = k_Q · ω²,方向与旋转方向相反。

四旋翼的电机布局采用X型或十字型配置,相邻电机旋转方向相反:1号和4号电机顺时针旋转,2号和3号电机逆时针旋转(或反之)。这种设计使得四个反扭矩在偏航轴上相互抵消,当四个电机转速相等时,偏航力矩为零。通过差动调节转速,可以产生三个轴上的控制力矩:

控制通道 电机转速关系 产生的力矩 物理机制
滚转 (Roll) Δω = ω₁ - ω₂ + ω₃ - ω₄ M_x = L·k_T·(ω₂² - ω₄²) 左右拉力差产生绕X轴力矩
俯仰 (Pitch) Δω = ω₁ + ω₂ - ω₃ - ω₄ M_y = L·k_T·(ω₁² - ω₃²) 前后拉力差产生绕Y轴力矩
偏航 (Yaw) Δω = ω₁ - ω₂ - ω₃ + ω₄ M_z = k_Q·(ω₁² - ω₂² + ω₃² - ω₄²) 反扭矩差产生绕Z轴力矩

其中L为电机轴到机体中心的距离(力臂)。值得注意的是,滚转和俯仰力矩由拉力差产生,而偏航力矩完全由反扭矩差产生,这解释了为什么四旋翼必须采用正反桨交替布局。

4.3 陀螺效应与耦合动力学

陀螺效应是四旋翼动力学中不可忽略的耦合项。每个高速旋转的螺旋桨相当于一个陀螺仪,当机体绕某一轴转动时,螺旋桨的角动量方向改变,会产生垂直于转动方向的陀螺力矩。该力矩的数学表达式为:M_gyro = J_p · (ω_body × Ω_prop),其中J_p为螺旋桨绕其旋转轴的转动惯量,ω_body为机体角速度,Ω_prop为螺旋桨旋转角速度向量。

具体而言,当飞行器绕X轴(滚转)以角速度p旋转时,四个螺旋桨的角动量方向随之改变,产生绕Y轴(俯仰)的陀螺力矩;反之亦然。这种交叉耦合在高速机动或大角速度飞行时尤为显著,可能导致控制器性能下降甚至发散。在建模时必须将陀螺力矩项显式加入欧拉方程中:

M_gyro = J_p · [q·(Ω₁ - Ω₂ + Ω₃ - Ω₄), 
                 -p·(Ω₁ - Ω₂ + Ω₃ - Ω₄), 
                 0]

其中p为滚转角速度,q为俯仰角速度。由于正反桨的角速度方向相反,陀螺力矩在偏航轴上的分量相互抵消,因此只影响滚转和俯仰两个通道。

工程风险提示:在高速偏航机动或大风扰条件下,陀螺效应力矩可能达到控制力矩的10%-20%。若在参数辨识中忽略陀螺项,会导致辨识出的惯性参数(Ixx、Iyy)偏大,进而影响控制增益的整定。建议在辨识实验中采用小幅度激励信号,并优先使用频率响应法分离陀螺效应的影响。

4.4 完整动力学方程与参数辨识准备

综合上述分析,四旋翼的完整动力学方程可写为:

平动方程(惯性系下):
m·ẍ = (cosφ·sinθ·cosψ + sinφ·sinψ)·ΣTᵢ
m·ÿ = (cosφ·sinθ·sinψ - sinφ·cosψ)·ΣTᵢ
m·z̈ = (cosφ·cosθ)·ΣTᵢ - m·g

转动方程(机体系下):
Ixx·ṗ = (Iyy - Izz)·q·r + L·k_T·(ω₂² - ω₄²) + J_p·q·(Ω₁ - Ω₂ + Ω₃ - Ω₄)
Iyy·q̇ = (Izz - Ixx)·p·r + L·k_T·(ω₁² - ω₃²) - J_p·p·(Ω₁ - Ω₂ + Ω₃ - Ω₄)
Izz·ṙ = (Ixx - Iyy)·p·q + k_Q·(ω₁² - ω₂² + ω₃² - ω₄²)

该方程组包含12个状态变量(位置3个、姿态3个、线速度3个、角速度3个),以及待辨识参数:质量m、惯性矩Ixx/Iyy/Izz、拉力系数k_T、反扭矩系数k_Q、螺旋桨转动惯量J_p、力臂L。在实际辨识中,m和L可通过称重和测量直接获得,k_T和k_Q可通过台架试验标定,而Ixx/Iyy/Izz和J_p则需要通过飞行实验或地面激振试验辨识。

4.5 案例拆解:悬停状态下的动力学参数辨识实验

以某型四旋翼(质量1.5kg,轴距450mm)为例,说明动力学模型在参数辨识中的应用。实验目标是通过悬停状态下的激励信号,辨识Ixx和Iyy。实验条件为无风室内环境,飞行器悬停于1.5m高度,采用扫频信号激励滚转通道。

实验步骤:

  1. 台架预标定:将电机安装于测力台架,测量不同PWM输入下的拉力T和反扭矩Q,拟合得到k_T = 1.21e-5 N/(rad/s)²,k_Q = 2.34e-7 N·m/(rad/s)²。
  2. 激励设计:在滚转通道注入频率为0.5-5Hz的线性扫频信号,幅值限制在±5°姿态角以内,避免触发非线性效应。
  3. 数据采集:以500Hz采样率记录陀螺仪角速度p、q和电机转速指令ω₁-ω₄,同步记录IMU姿态数据。
  4. 频响分析:将采集数据代入转动方程,利用FFT计算从电机差速指令到滚转角速度的传递函数,在频域中拟合惯性参数。
  5. 结果验证:辨识得到Ixx = 0.021 kg·m²,Iyy = 0.019 kg·m²,与CAD模型计算值(0.022/0.020)误差小于5%。

复盘要点:实验中发现,若激励幅值超过8°,陀螺效应项导致频响曲线在2Hz附近出现明显谐振峰,影响拟合精度。因此将激励幅值限制在5°以内,并在拟合模型中显式包含陀螺力矩项,最终辨识结果的一致性显著提升。此外,电机转速指令到实际转速的延迟(约20ms)也会引入相位误差,需要在数据处理时进行补偿。

经验提醒:在动力学建模时,建议将空气阻力项F_aero = -k_d·V(k_d为阻力系数)纳入平动方程。虽然悬停状态下阻力影响较小,但在高速前飞或抗风场景下,忽略阻力会导致位置估计漂移。阻力系数可通过风洞实验或飞行数据回归获得,典型值在0.1-0.5 N·s/m范围内。

4.6 常见误区与建模复盘

在四旋翼动力学建模中,初学者常犯以下错误:

  • 忽略陀螺效应:在低速悬停时陀螺力矩较小,但在快速机动中不可忽略。建议在模型中始终保留陀螺项,并通过仿真对比验证其影响。
  • 惯性张量简化为标量:虽然四旋翼对称布局使非对角项较小,但若载荷偏心或结构不对称,非对角项不可忽略。此时应使用完整的3×3惯性张量。
  • 拉力系数视为常数:k_T实际受空气密度、桨叶磨损、电池电压影响。在长时间飞行或高海拔场景下,应实时修正k_T值。
  • 忽略电机动态延迟:电机从指令到转速响应存在一阶惯性延迟(时间常数约20-50ms),在高速控制回路中必须建模为:ω̇ = (ω_cmd - ω)/τ_m。
  • 坐标系混淆:平动方程在惯性系下表达,转动方程在机体系下表达,两者通过旋转矩阵耦合。若混用坐标系,会导致严重的模型错误。

建模完成后,建议通过以下步骤验证模型准确性:首先在仿真环境中对比模型输出与真实飞行数据(开环响应);其次进行残差分析,检查模型预测误差是否满足白噪声假设;最后通过交叉验证,使用不同飞行工况的数据检验模型的泛化能力。只有经过充分验证的动力学模型,才能作为参数辨识和控制设计的可靠基础。