2. 纯比例导引律与真实比例导引律对比

任务定义

本案例实现二维平面内拦截弹对匀速运动目标的纯比例导引(PPN)与真实比例导引(TPN)制导律仿真。给定拦截弹初始位置、速度大小、目标初始位置与常值速度,分别用两种导引律生成加速度指令,数值积分弹目运动方程,对比脱靶量、拦截时间与需用过载。

制导律模型

设拦截弹位置为(xm, ym),速度为Vm,弹道偏角为γm;目标位置为(xt, yt),速度为Vt,弹道偏角为γt。视线角λ = atan2(yt−ym, xt−xm),视线角速率λ̇由相对运动计算。两种导引律的加速度指令均垂直于视线方向:

  • PPN:aPPN = N·Vm·λ̇,加速度方向为拦截弹速度方向旋转90°指向视线旋转方向;
  • TPN:aTPN = N·Vr·λ̇,其中Vr为弹目相对速度沿视线方向的分量(接近速度),加速度方向垂直于视线。

导航常数取N=3,拦截弹速度恒定,加速度指令仅改变弹道偏角。仿真终止条件为弹目距离小于等于脱靶量判定半径(取0.5 m)或达到最大仿真时间。输出脱靶量、拦截时间、最大需用过载,并绘制弹道轨迹与需用过载曲线。

仿真代码

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

# 物理常量与初始条件
N = 3.0                # 导航常数
Vm = 300.0             # 拦截弹速度 m/s
Vt = 150.0             # 目标速度 m/s
missile_r0 = np.array([0.0, 0.0])     # 拦截弹初始位置
target_r0 = np.array([3000.0, 0.0])   # 目标初始位置
gamma_m0 = np.deg2rad(30.0)           # 拦截弹初始弹道偏角
gamma_t0 = np.deg2rad(180.0)          # 目标初始弹道偏角(朝-x方向飞行)
t_max = 30.0                          # 最大仿真时间
R_miss = 0.5                          # 脱靶量判定半径

def los_rate(r_rel, v_rel):
    """计算视线角速率 lambda_dot"""
    rx, ry = r_rel
    vx, vy = v_rel
    cross = rx * vy - ry * vx
    dot = rx * vx + ry * vy
    R2 = rx**2 + ry**2
    if R2 < 1e-6:
        return 0.0
    return cross / R2

def guidance_ppn(r_rel, v_rel, Vm):
    """纯比例导引加速度矢量(垂直于视线)"""
    lam_dot = los_rate(r_rel, v_rel)
    a_mag = N * Vm * abs(lam_dot)
    # 视线方向单位向量
    R = np.linalg.norm(r_rel)
    e_los = r_rel / R
    # 视线法向单位向量(指向视线旋转方向)
    e_perp = np.array([-e_los[1], e_los[0]])
    # 根据视线角速率符号确定加速度方向
    if lam_dot >= 0:
        return a_mag * e_perp
    else:
        return -a_mag * e_perp

def guidance_tpn(r_rel, v_rel, Vm):
    """真实比例导引加速度矢量(垂直于视线,用接近速度)"""
    lam_dot = los_rate(r_rel, v_rel)
    R = np.linalg.norm(r_rel)
    e_los = r_rel / R
    Vr = -np.dot(v_rel, e_los)   # 接近速度(正值为接近)
    if Vr < 0:
        Vr = 0.0
    a_mag = N * Vr * abs(lam_dot)
    e_perp = np.array([-e_los[1], e_los[0]])
    if lam_dot >= 0:
        return a_mag * e_perp
    else:
        return -a_mag * e_perp

def dynamics(t, state, guidance_func):
    """状态: [xm, ym, gamma_m, xt, yt]"""
    xm, ym, gamma_m, xt, yt = state
    r_rel = np.array([xt - xm, yt - ym])
    v_m = np.array([Vm * np.cos(gamma_m), Vm * np.sin(gamma_m)])
    v_t = np.array([Vt * np.cos(gamma_t0), Vt * np.sin(gamma_t0)])
    v_rel = v_t - v_m

    a_vec = guidance_func(r_rel, v_rel, Vm)
    # 加速度垂直于速度方向,改变弹道偏角
    gamma_m_dot = np.cross(v_m, a_vec) / Vm**2

    return [Vm * np.cos(gamma_m), Vm * np.sin(gamma_m), gamma_m_dot,
            Vt * np.cos(gamma_t0), Vt * np.sin(gamma_t0)]

def simulate(guidance_func, label):
    state0 = [missile_r0[0], missile_r0[1], gamma_m0,
              target_r0[0], target_r0[1]]
    # 事件:弹目距离小于R_miss时终止
    def hit_event(t, state):
        xm, ym, _, xt, yt = state
        return np.hypot(xt - xm, yt - ym) - R_miss
    hit_event.terminal = True
    hit_event.direction = -1

    sol = solve_ivp(dynamics, [0, t_max], state0, method='RK45',
                    events=hit_event, max_step=0.05,
                    args=(guidance_func,), rtol=1e-8, atol=1e-8)

    # 提取轨迹
    xm = sol.y[0]; ym = sol.y[1]; gamma_m = sol.y[2]
    xt = sol.y[3]; yt = sol.y[4]
    t = sol.t

    # 计算脱靶量(最终距离)
    miss_dist = np.hypot(xt[-1] - xm[-1], yt[-1] - ym[-1])
    intercept_time = t[-1]

    # 计算需用过载(加速度大小/g,取g=9.8)
    a_hist = []
    for i in range(len(t)):
        r_rel = np.array([xt[i] - xm[i], yt[i] - ym[i]])
        v_m = np.array([Vm * np.cos(gamma_m[i]), Vm * np.sin(gamma_m[i])])
        v_t = np.array([Vt * np.cos(gamma_t0), Vt * np.sin(gamma_t0)])
        v_rel = v_t - v_m
        a_vec = guidance_func(r_rel, v_rel, Vm)
        a_hist.append(np.linalg.norm(a_vec) / 9.8)
    a_hist = np.array(a_hist)
    max_a = np.max(a_hist) if len(a_hist) > 0 else 0.0

    print(f"--- {label} ---")
    print(f"脱靶量: {miss_dist:.3f} m")
    print(f"拦截时间: {intercept_time:.3f} s")
    print(f"最大需用过载: {max_a:.3f} g")
    print()

    return xm, ym, xt, yt, t, a_hist

# 运行两种导引律仿真
xm_ppn, ym_ppn, xt_ppn, yt_ppn, t_ppn, a_ppn = simulate(guidance_ppn, "PPN")
xm_tpn, ym_tpn, xt_tpn, yt_tpn, t_tpn, a_tpn = simulate(guidance_tpn, "TPN")

# 绘制弹道轨迹
plt.figure(figsize=(10, 5))
plt.subplot(1, 2, 1)
plt.plot(xm_ppn, ym_ppn, 'b-', linewidth=1.5, label='PPN导弹')
plt.plot(xt_ppn, yt_ppn, 'r--', linewidth=1.5, label='目标')
plt.plot(xm_tpn, ym_tpn, 'g-', linewidth=1.5, label='TPN导弹')
plt.plot(xt_tpn, yt_tpn, 'k--', linewidth=1.0, alpha=0.5)
plt.scatter([missile_r0[0]], [missile_r0[1]], color='black', marker='o', s=30, label='拦截弹起点')
plt.scatter([target_r0[0]], [target_r0[1]], color='purple', marker='s', s=30, label='目标起点')
plt.xlabel('x (m)'); plt.ylabel('y (m)')
plt.title('拦截弹道轨迹对比')
plt.legend(); plt.grid(True); plt.axis('equal')

# 绘制需用过载
plt.subplot(1, 2, 2)
plt.plot(t_ppn, a_ppn, 'b-', linewidth=1.5, label='PPN')
plt.plot(t_tpn, a_tpn, 'g-', linewidth=1.5, label='TPN')
plt.xlabel('时间 (s)'); plt.ylabel('需用过载 (g)')
plt.title('需用过载对比')
plt.legend(); plt.grid(True)

plt.tight_layout()
plt.savefig('guidance_comparison.png', dpi=150)
plt.show()

运行与验证

运行环境需安装numpy、scipy、matplotlib。执行代码后,终端输出两种导引律的脱靶量、拦截时间与最大需用过载,同时生成guidance_comparison.png图像文件。验证判定标准:

  1. 两种导引律的脱靶量均应小于1 m(在判定半径0.5 m附近触发终止,实际脱靶量应接近0.5 m量级);
  2. PPN与TPN的拦截时间应接近,差异通常在1 s以内;
  3. 轨迹图中两条弹道均平滑收敛至目标轨迹终点附近,无发散振荡。

若需对比不同初始条件,可修改missile_r0、target_r0、gamma_m0、gamma_t0等参数后重新运行。例如将目标初始弹道偏角改为90°(朝+y方向飞行),可观察两种导引律在迎击与尾追场景下的性能差异。

PPN与TPN制导律仿真验证流程 初始条件 Vm=300 m/s, Vt=150 m/s 导弹(0,0), γm=30° 目标(3000,0), γt=180° 制导律 PPN a = N·Vm·λ̇ TPN a = N·Vr·λ̇ 数值积分 (RK45) solve_ivp 求解运动方程 终止条件: R ≤ 0.5 m 最大时间: 30 s max_step=0.05 s 输出指标 脱靶量 (m) 拦截时间 (s) 最大需用过载 (g) 验证判定标准 ✓ 脱靶量 < 1 m ✓ 拦截时间差异 < 1 s ✓ 弹道平滑收敛无发散 对比分析 可视化输出 guidance_comparison.png — 弹道轨迹对比 + 需用过载对比