2023年五一杯（A题）无人机动力学方法建模！问题代码全解析！

背景中，随着科学技术的不断发展，无人机在许多领域都有着广泛的应用。对于空中执行定点投放任务的无人机，其投放精度不仅依赖于无人机的操作技术，而且还与无人机执行任务时所处状态和环境有关，例如在接近投放点时无人机的高度、速度，无人机所处位置的风速、投放点周围地理环境等因素。我们需要处理的是关于飞机投放的一系列问题，需要根据飞机的特性进行一定的。

qq_25834913

1465人浏览 · 2023-08-25 15:22:05

qq_25834913 · 2023-08-25 15:22:05 发布

问题重述

问题1：

（1）建立数学模型：

无人机投放模型在这个问题中的作用是建立数学模型来描述无人机投放爆炸物的过程，并且可以通过该模型来优化无人机投放的策略，从而提高命中率和效率。具体来说，该模型可以通过考虑无人机的飞行高度、飞行速度、俯冲角度、发射速度等因素来确定最佳的发射距离和发射时机，以确保物体能够准确地命中目标。此外，该模型还可以考虑外部因素，如风速和风向等，来调整无人机的飞行轨迹和姿态，以提高投放精度和稳定性。

对于本题的模型，有：

$\frac{1}{2}\rho v^{2}S C_{d}=(m+\frac{4}{3}\pi r^{3}\rho)g+F$

其中 ρ 为空气密度，S 为物资横截面积， $C_d$ 为物资的阻力系数， g 为重力加速度，F 为无人机与物资之间的牵引力。

$\frac{1}{2}\rho v^{2}S C_{d}=(m+\frac{4}{3}\pi r^{3}\rho)g$

当无人机投放物资时，物资与无人机之间断开连接，牵引力 F取0，上式可以化简为：

其中 $g=9.8m/s^2$ 为重力加速度。

（2）在无人机的飞行高度为 300m，飞行速度为 300km/h，风速为 5m/s，风向与水平面平行的情况下，使用代码解决：

import math
 
# 定义常量
v0 = 300      # 飞行速度，单位km/h
vw = 5        # 风速，单位m/s
h = 300       # 飞行高度，单位m
r = 0.2       # 球形物资半径，单位m
m = 50        # 球形物资质量，单位kg
g = 9.8       # 重力加速度，单位m/s^2
 
# 计算投放距离
d0 = v0**2/g * math.sin(0*2*math.pi/360) + vw*v0/g * math.cos(0*2*math.pi/360)
d180 = v0**2/g * math.sin(180*2*math.pi/360) + vw*v0/g * math.cos(180*2*math.pi/360)
d90 = v0**2/g * math.sin(90*2*math.pi/360) + vw*v0/g * math.cos(90*2*math.pi/360)
 
# 输出结果
print(f"无人机飞行方向与风向相同时，投放距离为：{d0:.1f}m")
print(f"无人机飞行方向与风向相反时，投放距离为：{d180:.1f}m")

问题2：

假设无人机在水平飞行过程中到达距离目标点的水平距离为 x ，飞行高度为 ℎ ，飞行速度为v ，俯冲角度为 α ，发射速度为 u 。则无人机发射炸弹的轨迹可以分解为水平方向和竖直方向两个分量。

在水平方向上，无人机在 $t_1=\frac{x}{v}$ 秒到达目标点，发射炸弹的时间为 $t_2=t_1+\frac{h}{v}tan\alpha$ 秒。发射炸弹时无人机的水平速度为vcos⁡α ，炸弹的水平初速度为 ucos⁡α。

在竖直方向上，炸弹自由落体运动，竖直初速度为 usin⁡α ，竖直加速度为g 。设炸弹飞行的时间为 $t_3$ 秒，则有：

$h=ut_3sin\alpha -\frac{1}{2}gt_3^2$

将 t3 的值代入到水平方向上的运动中，则可以求得无人机与目标点之间的距离 x1 ：

$d=\sqrt{x_1^2+h_1^2}$

假设无人机发射炸弹的距离为 d ，则需要满足 1000≤d≤3000 。为了使无人机在发射炸弹时仍能保持安全的飞行高度，假设无人机的飞行高度为 800m ，则有 ℎ≥300m 。

为了使发射策略可行，需要选择合适的俯冲角度 α 和发射速度u。假设 α 为定值，可以根据上述模型求出发射速度u与发射距离d的关系，并绘制出其图像，如下图所示。

由图像可知，当俯冲角度为 30∘ 时，发射速度最小，约为 426.8m/s ，此时发射距离为d≈1716.2m

问题3：

无人机的飞行稳定性可以用无人机的俯仰角和偏航角的变化率来描述，即：

$Stability=\frac{\mathrm{d}\theta }{\mathrm{d} t}+\frac{\mathrm{d}\psi }{\mathrm{d} t}$

其中， θ表示俯仰角， $\psi$ 表示偏航角。这个数值越小，说明无人机的飞行越稳定。

无人机的命中精度可以用命中目标的距离来描述，与无人机的飞行稳定性呈反比关系，即：

$Accuracy\propto \frac{1}{Stability}$

在实际应用中，可以通过无人机的传感器数据来计算无人机的俯仰角和偏航角的变化率，并根据上述公式来评估无人机的飞行稳定性和命中精度。

import numpy as np
import matplotlib.pyplot as plt


def dynamic_equation(x, u):
    # 状态方程
    A = np.array([[1, 0, 0, dt, 0, 0],
                  [0, 1, 0, 0, dt, 0],
                  [0, 0, 1, 0, 0, dt],
                  [0, 0, 0, 1-0.5*rho*Cd*S/m*dt, 0, 0],
                  [0, 0, 0, 0, 1-0.5*rho*Cd*S/m*dt, 0],
                  [0, 0, 0, 0, 0, 1-0.5*rho*Cd*S/m*dt]])

    # 输入方程
    B = np.array([[0, 0, 0, 0],
                  [0, 0, 0, 0],
                  [0, 0, 0, 0],
                  [0, 0, 0, -0.5*rho*S*vx0**2/m*dt],
                  [0, 0, 0, -0.5*rho*S*vx0**2/m*dt],
                  [0, 0, 0, -0.5*rho*S*vx0**2/m*dt]])

    # 状态更新
    x_new = np.dot(A, x) + np.dot(B, u)

    return x_new

# 定义无人机飞行过程的仿真函数
def simulate_flight(x0, u, t):
    # 初始化状态和控制输入
    x = x0
    u = u.reshape(-1, 1)

    # 初始化状态列表和控制输入列表
    x_list = [x]
    u_list = [u]

    # 循环仿真
    for i in range(len(t)):
        # 计算下一个状态
        x = dynamic_equation(x, u)

        # 记录状态和控制输入
        x_list.append(x)、、模糊处理 完整版看文章下面~
        u_list.append(u)

    # 将列表转换为数组
    x_array = np.array(x_list)
    u_array = np.array(u_list)

    return x_array, u_array

# 无人机和环境参数设置
h = 800  # 飞行高度，单位：m
v0 = 300  # 无人机飞行速度，单位：km/h
v = np.linspace(300, 400, 101) / 3.6  # 无人机相对地面速度，单位：m/s
vw = np.array([6, 0])  # 风速，单位：m/s
gamma = np.deg2rad(45)  # 俯冲角，单位：rad
g = 9.8  # 重力加速度，单位：m/s^2

# 计算无人机稳定性
S = 2 * np.pi * (0.5 ** 2)  # 球形爆炸物的参考面积
Cd = 0.5  # 球形爆炸物的阻力系数
m = 50  # 球形爆炸物的质量，单位：kg
rho = 1.2  # 空气密度，单位：kg/m^3
K = 0.5 * rho * S * Cd / m  # 阻力系数
u = np.sqrt(v ** 2 + (v0 * np.sin(gamma)) ** 2)  # 爆炸物相对空气速度
D = K * u ** 2  # 阻力大小
H = h - np.sqrt((h ** 2) / (np.tan(gamma) ** 2 + 1))  # 爆炸物发射高度
t = (H - 300) / (v0 * np.cos(gamma))  # 爆炸物发射时间
x0 = v0 * t  # 无人机前进距离
x = x0 + (v + vw[0]) * t  # 爆炸物水平位移距离
y = H - (v0 * np.sin(gamma) + (g + vw[1]) * t) * t / 2  # 爆炸物垂直位移距离
delta_x = 20 / 2  # 爆炸物命中误差，单位：cm
sigma = delta_x / 3  # 标准差
P = 1 / (sigma * np.sqrt(2 * np.pi)) * np.exp(-(x - delta_x) ** 2 / (2 * sigma ** 2))  # 命中概率密度函数
hit_rate = np.trapz(P, x)  # 命中率

# 可视化结果
plt.plot(v, P)
plt.xlabel('Horizontal displacement (m)')
plt.ylabel('Probability density')
plt.title('Hit probability density')
plt.show()

print('The hit rate is %.2f%%.' % (hit_rate * 100))