飞行力学基础：刚体运动方程、坐标系定义、欧拉角与四元数、力和力矩建模
飞行控制系统设计的核心任务，是让飞行器按照期望的姿态、轨迹和速度运动。而这一切的前提，是能够用数学语言精确描述飞行器的运动状态及其变化规律。本章作为《飞行控制系统建模与调参手册》的力学基石，将系统梳理刚体运动方程、坐标系定义、姿态描述方法以及力和力矩的建模路径。这些内容不仅是后续气动导数辨识、控制律设计和仿真验证的基础，也是工程调试中判断“模型是否可信”的关键依据。理解本章，意味着你能够从物理本质出发，看懂飞控代码中每一个状态量的来源与去向。
2.1 坐标系定义与转换关系
飞行力学中，坐标系的选取直接决定了运动方程的形式和复杂度。工程实践中，最常用的四个坐标系为：地面惯性坐标系（E系）、机体坐标系（B系）、风轴坐标系（W系）和航迹坐标系（K系）。地面系固定于地球表面某点，忽略地球自转影响时近似为惯性系，是导航和位置解算的基准；机体系固连于飞行器质心，x轴指向机头，y轴指向右翼，z轴按右手定则指向下方，是力和力矩分解的自然参考；风轴系以空速向量为x轴，用于气动力的分解；航迹系则以地速向量为x轴，用于分析飞行轨迹。
坐标系之间的转换通过旋转矩阵实现。从地面系到机体系的转换通常采用Z-Y-X顺序的欧拉角旋转，即先偏航角ψ绕z轴旋转，再俯仰角θ绕新y轴旋转，最后滚转角φ绕新x轴旋转。对应的方向余弦矩阵为：
C_EB = R_x(φ)·R_y(θ)·R_z(ψ)
其中：
R_z(ψ) = [cosψ sinψ 0; -sinψ cosψ 0; 0 0 1]
R_y(θ) = [cosθ 0 -sinθ; 0 1 0; sinθ 0 cosθ]
R_x(φ) = [1 0 0; 0 cosφ sinφ; 0 -sinφ cosφ]风轴系与机体系之间仅相差两个角度：迎角α（速度向量在机体对称面内的投影与x轴的夹角）和侧滑角β（速度向量与对称面的夹角）。气动力的建模通常在风轴系中完成，再通过转换矩阵投影到机体系中参与运动方程求解。
2.2 欧拉角与四元数：姿态描述与转换
欧拉角以三个角度直观描述姿态，物理意义明确，适合人机交互和可视化。但其存在两个固有缺陷：一是万向锁问题，当俯仰角接近±90°时，偏航和滚转的旋转轴重合，自由度丢失；二是三角函数运算量大，且微分方程中存在奇异点。四元数则通过四个参数（一个标量加三个矢量分量）描述刚体旋转，避免了奇异问题，且运算效率高，是当前飞控系统姿态解算的主流方法。
四元数定义为 q = q0 + q1·i + q2·j + q3·k，其中 q0 为标量部分，满足归一化约束 q0² + q1² + q2² + q3² = 1。从四元数到欧拉角的转换关系为：
φ = atan2(2(q0·q1 + q2·q3), 1 - 2(q1² + q2²))
θ = asin(2(q0·q2 - q3·q1))
ψ = atan2(2(q0·q3 + q1·q2), 1 - 2(q2² + q3²))姿态运动学方程在四元数形式下为线性微分方程，避免了欧拉角方程中的三角运算和奇异点。工程中，陀螺仪输出的角速度通过四元数积分更新姿态，再定期用加速度计和磁力计数据通过互补滤波或卡尔曼滤波进行修正。
| 对比维度 | 欧拉角 | 四元数 |
|---|---|---|
| 参数数量 | 3个 | 4个（1个约束） |
| 物理直观性 | 高，直接对应姿态角 | 低，需转换理解 |
| 奇异性 | 俯仰±90°时万向锁 | 无奇异 |
| 计算效率 | 三角函数多，运算量大 | 仅乘加运算，效率高 |
| 插值平滑性 | 大角度插值可能跳变 | 球面线性插值平滑 |
| 工程应用 | 显示、地面站、初始对准 | 姿态解算、控制律内部计算 |
2.3 刚体六自由度运动方程
将飞行器视为理想刚体，其运动可分解为质心平动和绕质心转动，共六个自由度。在机体坐标系下，牛顿-欧拉方程展开为十二阶状态方程。力方程描述线速度变化，力矩方程描述角速度变化，运动学方程描述姿态变化率，导航方程描述位置变化率。
力方程的完整形式为：
m(du/dt + qw - rv) = X_A + X_T + X_G
m(dv/dt + ru - pw) = Y_A + Y_T + Y_G
m(dw/dt + pv - qu) = Z_A + Z_T + Z_G其中下标A表示气动力，T表示推力，G表示重力分量。注意重力在机体轴上的投影取决于当前姿态角，这体现了姿态与线运动之间的耦合。力矩方程中，惯性张量I通常假设为对角阵（忽略惯性积），但若飞行器质量分布不对称（如挂载不对称），则需考虑惯性积项。
运动学方程将角速度与姿态变化率联系起来。在四元数形式下，姿态更新方程为：
dq/dt = 0.5 · Ω(ω) · q
其中 Ω(ω) = [0 -p -q -r;
p 0 r -q;
q -r 0 p;
r q -p 0]导航方程则将机体速度转换到地面系，积分得到位置。完整的十二阶状态向量为 [u, v, w, p, q, r, φ, θ, ψ, x_E, y_E, h]（或使用四元数替代欧拉角）。
2.4 力和力矩的建模方法
力和力矩的建模是飞行力学中最具挑战性的部分，直接决定了模型的保真度。主要力源包括：气动力、推力、重力和地面接触力（起降阶段）。气动力通常采用导数模型（稳定性导数法），将力和力矩展开为状态量和控制量的函数。
以纵向通道为例，升力系数和俯仰力矩系数的线性化表达式为：
C_L = C_L0 + C_Lα·α + C_Lq·(q·c̄/2V) + C_Lδe·δe
C_m = C_m0 + C_mα·α + C_mq·(q·c̄/2V) + C_mδe·δe其中c̄为平均气动弦长，V为空速，δe为升降舵偏角。这些气动导数可通过风洞试验、CFD计算或系统辨识获得。对于飞控设计，通常需要线性化模型，即在配平点附近将非线性方程展开为状态空间形式：
ẋ = A·x + B·u
y = C·x + D·u推力建模需考虑发动机的动态响应（一阶惯性环节）、推力随高度和速度的变化。重力建模相对简单，只需将重力加速度投影到机体轴。对于小型无人机，还需考虑螺旋桨滑流对气动面的影响。
背景与约束：某翼展2.4m、起飞重量8kg的固定翼无人机，需设计俯仰姿态控制器。要求控制器带宽不低于2rad/s，相位裕度大于45°。项目周期紧张，无法进行大量风洞试验。
初始现象：初步采用经验气动导数搭建的模型，在仿真中俯仰响应振荡剧烈，阻尼比仅0.3，与试飞数据对比发现模型预测的短周期频率偏高约20%。
关键数据：试飞采集的俯仰角阶跃响应显示，超调量28%，调节时间3.2s；模型仿真超调量35%，调节时间2.1s。两者差异明显。
分析判断：对比模型与试飞数据，发现C_mα（俯仰静稳定性导数）估计偏大，导致模型过于稳定，短周期频率偏高。同时C_mq（俯仰阻尼导数）偏小，导致阻尼比偏低。
实际操作：采用频域系统辨识方法，利用扫频激励信号（频率范围0.1-10Hz）激励升降舵，记录俯仰角速率响应。通过FFT计算频率响应，拟合传递函数，辨识出修正后的气动导数：C_mα从-0.45修正为-0.38，C_mq从-8.5修正为-10.2。
结果对比：修正后模型仿真与试飞数据的俯仰阶跃响应超调量误差小于3%，调节时间误差小于0.3s，短周期频率误差小于5%。基于修正模型设计的控制器，在试飞中达到带宽2.3rad/s，相位裕度52°。
风险与代价：系统辨识需要额外的试飞架次，增加了约2个飞行起落的成本。但避免了因模型不准导致的控制器反复迭代，总体节省了约2周开发时间。
复盘经验：气动导数经验公式在小尺寸无人机上误差较大，建议在首飞前安排专门的辨识试飞科目。同时，模型验证不能只看稳态响应，必须对比动态响应指标（频率、阻尼）。
2.5 常见误区与工程总结
在飞行力学建模与调参过程中，以下误区值得特别注意：
- 坐标系混淆：将地面系下的速度分量误当作机体系分量使用，导致力和力矩计算错误。建议在代码中为每个变量明确标注坐标系后缀（如_EB表示从E到B的转换）。
- 忽略惯性积：对于非对称布局或挂载不对称的飞行器，惯性积不可忽略，否则力矩方程将产生显著误差。
- 四元数未归一化：长时间仿真或飞行后，四元数模长漂移导致姿态解算发散。务必在每个控制周期执行归一化。
- 线性化点选择不当：线性模型仅在配平点附近有效，若飞行包线宽广，需在不同飞行状态点分别建立线性模型，并设计增益调度。
- 气动导数外推风险：大迎角或高马赫数下，气动导数非线性显著，线性模型不再适用。需引入非线性修正项或查表模型。