📚 课程目录第4章纵向长周期模态分析

4. 纵向长周期模态分析

任务定义

针对某固定翼无人机纵向小扰动线性模型，提取长周期（phugoid）模态的固有频率与阻尼比，分析其运动特征，并评估该模态对俯仰角控制律设计的影响。本案例以短周期模态已解耦为前提，仅处理长周期近似模型。

纵向长周期近似模型

在纵向小扰动方程中，长周期模态主要表现为速度与俯仰角的缓慢振荡，迎角变化较小。忽略迎角变化率项后，可得如下二阶近似系统（以速度偏差 ΔV 为输出，升降舵偏角 Δδ_e 为输入）：

ΔV̈ + 2ζ<sub>p</sub>ω<sub>p</sub>ΔV̇ + ω<sub>p</sub>²ΔV = K<sub>p</sub>Δδ<sub>e</sub>

其中特征参数由以下气动导数确定：

ω<sub>p</sub>² ≈ -2g·Z<sub>u</sub>/V₀ ，  2ζ<sub>p</sub>ω<sub>p</sub> ≈ -X<sub>u</sub>

式中 V₀ 为配平空速，g 为重力加速度，X_u、Z_u 为无量纲纵向力导数（已除以质量）。本例给定某无人机在高度 1000 m、V₀ = 50 m/s 配平状态下的导数值：X_u = -0.03 s⁻¹，Z_u = -0.01 s⁻¹，g = 9.81 m/s²。

模态参数计算

import numpy as np

# 给定配平状态与气动导数
V0 = 50.0          # 配平空速 m/s
g = 9.81           # 重力加速度 m/s²
Xu = -0.03         # 轴向力导数 s⁻¹
Zu = -0.01         # 法向力导数 s⁻¹

# 长周期自然频率与阻尼比（近似公式）
omega_p = np.sqrt(-2 * g * Zu / V0)
zeta_p = -Xu / (2 * omega_p)

print(f"长周期自然频率: {omega_p:.4f} rad/s")
print(f"长周期阻尼比:   {zeta_p:.4f}")
print(f"周期 T = {2*np.pi/omega_p:.2f} s")
print(f"半衰期 t_half = {np.log(2)/(zeta_p*omega_p):.2f} s")

预期输出（由给定公式直接推出）：

长周期自然频率: 0.0626 rad/s
长周期阻尼比:   0.2395
周期 T = 100.35 s
半衰期 t_half = 46.20 s

模态特性分析

计算结果显示：长周期自然频率仅 0.063 rad/s，对应周期约 100 s，远慢于短周期（通常 1~10 rad/s）。阻尼比约 0.24，属于欠阻尼，意味着速度扰动后将以缓慢衰减的振荡形式恢复。该模态下飞机保持近似等迎角、等升力系数，动能与势能相互转换，表现为俯仰角与速度的同相缓慢波动。

对俯仰角控制的影响评估

长周期模态对俯仰角控制律的影响主要体现在以下方面：

带宽分离：长周期频率远低于短周期，设计俯仰角内回路（角速率/姿态反馈）时可将长周期视为准静态扰动，不会与内回路产生耦合失稳。
速度反馈必要性：由于长周期欠阻尼且周期长，若仅用俯仰角反馈，速度漂移无法被有效抑制。需在速度回路引入空速反馈或油门补偿，否则长周期振荡将持续数十秒才衰减。
控制律增益约束：若在长周期频段（约 0.06 rad/s）提供过高增益，可能激发该模态，导致速度与俯仰角出现缓慢发散振荡。因此外回路（高度/速度）带宽应低于 0.3 rad/s，与内回路保持至少 10 倍频程间隔。

验证方法

将上述二阶模型转换为状态空间形式，在 Python 中仿真给定初始速度扰动 ΔV(0)=2 m/s 的零输入响应，观察 200 s 内的速度时间历程。判定标准：

响应曲线呈现衰减振荡，振荡周期约 100 s，与计算值偏差小于 5%；
峰值衰减符合阻尼比 0.24 对应的包络线（每周期幅值衰减约 e^{-2πζ/√(1-ζ²)} ≈ 0.22 倍）；
若在上述模型中加入俯仰角比例反馈增益 K=0.5，系统仍稳定，但速度响应衰减变慢，验证长周期对高增益的敏感性。

import numpy as np
from scipy.integrate import odeint
import matplotlib.pyplot as plt

# 状态空间实现：x = [ΔV, ΔVdot]
A = np.array([[0, 1],
              [-omega_p**2, -2*zeta_p*omega_p]])
B = np.array([0, 1.0])  # 单位增益，仅用于验证

def sys_deriv(x, t):
    return A @ x

t = np.linspace(0, 200, 2000)
x0 = [2.0, 0.0]  # 初始速度偏差 2 m/s
sol = odeint(sys_deriv, x0, t)

# 理论包络线
env = 2.0 * np.exp(-zeta_p*omega_p*t)

plt.figure(figsize=(8, 3))
plt.plot(t, sol[:, 0], label='ΔV response')
plt.plot(t, env, 'r--', label='damping envelope')
plt.plot(t, -env, 'r--')
plt.xlabel('Time (s)'); plt.ylabel('ΔV (m/s)')
plt.legend(); plt.grid(True)
plt.title('Phugoid Free Response (V0=50 m/s)')
plt.tight_layout()
plt.savefig('phugoid_response.png', dpi=100)
print("仿真完成，响应曲线已保存为 phugoid_response.png")

判定标准

完成上述仿真后，检查响应曲线是否满足：① 振荡周期在 95~105 s 范围内；② 第二个峰值（约 100 s 处）幅值约为初始幅值的 0.2~0.25 倍；③ 200 s 时响应幅值衰减至初始值的 5% 以下。满足以上三点即判定长周期模态分析正确，控制影响评估有效。若仿真结果与理论计算偏差超过 5%，需检查状态空间矩阵 A 的构建是否正确。