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 m(在判定半径0.5 m附近触发终止,实际脱靶量应接近0.5 m量级);
- PPN与TPN的拦截时间应接近,差异通常在1 s以内;
- 轨迹图中两条弹道均平滑收敛至目标轨迹终点附近,无发散振荡。
若需对比不同初始条件,可修改missile_r0、target_r0、gamma_m0、gamma_t0等参数后重新运行。例如将目标初始弹道偏角改为90°(朝+y方向飞行),可观察两种导引律在迎击与尾追场景下的性能差异。