非线性规划

非线性规划

飞行管理问题

在约 10,000m 高空的某边长 160km 的正方形区域内,经常有若干架飞机作水平
飞行。区域内每架飞机的位置和速度向量均由计算机记录其数据,以便进行飞行管理。
当一架欲进入该区域的飞机到达区域边缘时,记录其数据后,要立即计算并判断是否会
与区域内的飞机发生碰撞。如果会碰撞,则应计算如何调整各架(包括新进入的)飞机
飞行的方向角,以避免碰撞。现假定条件如下:
1)不碰撞的标准为任意两架飞机的距离大于 8km;
2)飞机飞行方向角调整的幅度不应超过 30 度;
3)所有飞机飞行速度均为每小时 800km;
4)进入该区域的飞机在到达区域边缘时,与区域内飞机的距离应在 60km 以上;
5)最多需考虑 6 架飞机;
6)不必考虑飞机离开此区域后的状况。
请你对这个避免碰撞的飞行管理问题建立数学模型,列出计算步骤,对以下数据进
行计算(方向角误差不超过 0.01 度),要求飞机飞行方向角调整的幅度尽量小。
设该区域 4 个顶点的座标为(0,0),(160,0),(160,160),(0,160)。记录数据见表 1。

表1

飞机编号 横坐标\(\, x\) 纵坐标\(\, y\) 方向角(度)
1 150 140 243
2 85 85 223
3 150 155 220
4 145 50 159
5 130 150 230
新进入 0 0 52

注:方向角指飞行方向与\(x\)轴正向的夹角

为了方便以后的讨论,引进如下记号:

\(D\)为飞行管理区域的边长

\(\Omega\)为飞行管理区域,取直角坐标系使其为\([0,D]\times[0,D]\)

\(a\)为飞机飞行速度,\(a=800km/h\)

\((x_i^0,y_i^0)\)为第\(i\)架飞机的初始位置

\((x_i(t),y_i(t))\)为第\(i\)架飞机在\(t\)时刻的位置

\(\theta_i^0\)为第\(i\)架飞机的的原飞行角度,即飞行方向与\(x\)轴夹角,\(0\le\theta_i^0<2\pi\)

\(\Delta\theta_i\)为第\(i\)架飞机的方向角调整,\(-\frac\pi6\le\Delta\theta_i\le\frac\pi6\)

\(\theta_i=\theta_i^0+\Delta\theta_i\)为第\(i\)架飞机调整后的飞行方向角

根据相对运动观点在考察两架飞机\(i\)和\(j\)的飞行时,可以将飞机\(i\)视为不动,而飞机\(j\)以相对速度

\[v_{ij}=v_j-v_i=(a\cos\theta_j=a\cos\theta_i,a\sin\theta_j-a\sin\theta_i) \\=a(-2\sin\frac{\theta_j+\theta_j}2\sin\frac{\theta_j-\theta_i}2,2\cos\frac{\theta_j+\theta_j}2\sin\frac{\theta_j-\theta_i}2) \\=2a\sin\frac{\theta_j-\theta_i}2(\cos(\frac{\theta_j+\theta_j}2+\frac\pi2),\sin(\frac{\theta_j+\theta_j}2+\frac\pi2)) \]

不妨设\(\theta_j\ge\theta_i\),此时相对飞行方向角为\(\beta_{ij}=\frac\pi2+\frac{\theta_i+\theta_j}2\)

由于两架飞机的初始距离为

图1

非线性规划 $$ r_{ij}(0)=\sqrt{(x_i^0-x_j^0)^2+(y_i^0-y_j^0)^2}\\ \alpha_{ij}^0=\arcsin\frac8{r_{ij}(0)} $$

只要相对飞行方向角满足

\[\alpha_{ij}^0<\beta_{ij}<2\pi-\alpha_{ij}^0 \]

时,两架飞机不能碰撞

记\(\beta_{ij}^0\)为调整前第\(j\)架飞机相对于第\(i\)架飞机的相对速度与这两架飞机连线(从\(j\)指向\(i\)的矢量)的夹角(以连线矢量为基准,逆时针方向为正,顺时针方向为负),两架飞机不碰撞的条件为

\[|\beta_{ij}^0+\frac12(\Delta\theta_i+\Delta\theta_j)|>\alpha_{ij}^0 \]

其中

\[\beta_{mn}^0=相对速度v_{mn}的辐角-从n指向m的连线矢量的辐角\\ =\arg\frac{e^{i\theta_n}-e^{i\theta_m}}{(x_m+iy_m)-(x_n+iy_n)} \]

\(i\)是虚数单位,这里利用复数的辐角,来计算角度\(\beta_{mn}^0\,(m,n=1,2,\cdots,6)\)

本问题中的优化目标可以有不同的形式:如是所有飞机的最大调节量最小;所有飞机的调整量绝对值之和最小等。这里以所有飞机的调整量绝对值之和为目标函数,可以得到如下的数学规划模型:

\[\min\sum_{i=1}^6|\Delta\theta_i|\\ s.t\quad|\beta_{ij}^0+\frac12(\Delta\theta_i+\Delta\theta_j)|>\alpha_{ij}^0\,i,j=1,2\cdots,6,i\neq j\\ |\Delta\theta_i|\le30°,\qquad i=1,2,\cdots,6 \]

利用如下的程序,使用python:

import numpy as np

x0 = np.array([150, 85, 150, 145, 130, 0]) # 初始横坐标
y0 = np.array([140, 85, 155, 50, 150, 0]) # 初始纵坐标
theta_i = np.array([243, 236, 220.5, 159, 230, 52])
xy0 = np.concatenate([[x0], [y0]], axis=0)
d0 = np.zeros((6, 6))
for i in range(xy0.shape[1]):  # 计算矩阵各个列向量之间的距离
    for j in range(xy0.shape[1]):
        d0[i, j] = np.linalg.norm(xy0[:, i] - xy0[:, j])
d0[d0 == 0] = np.inf  # 设置点到自身的距离为无穷大
a0 = np.arcsin(8 / d0) / np.pi * 180  # arcsin()返回弧度制,换算成角度,获得不碰撞的最小角度
xy1 = np.array([])
for i, j in zip(x0, y0):  # 计算x+yi
    xy1 = np.append(xy1, np.complex(i, j))
xy2 = np.array([])
for i in theta_i:# 计算e_iθ
    xy2 = np.append(xy2, np.exp(np.complex(0, i / 180 * np.pi)))
b0 = np.zeros((6, 6))
for i in range(6):
    for j in range(6):
        if i != j: # 不计算对角线
            b0[i, j] = np.angle((xy2[j] - xy2[i]) / (xy1[j] - xy1[i]), deg=True)  # 范围角度值,范围在-180°到180°
            # 获得两点之间的相对角度
np.savetxt('./a0.csv', a0, delimiter=',', fmt='%f')  # 写入数据
np.savetxt('./b0.csv', b0, delimiter=',', fmt='%f')

求得\(\alpha_{ij}^0\)的值如表2所示

表2

0 5.39119 32.230953 5.091816 20.963361 2.234507
5.39119 0 4.804024 6.61346 5.807866 3.815925
32.230953 4.804024 0 4.364672 22.833654 2.125539
5.091816 6.61346 4.364672 0 4.537692 2.989819
20.963361 5.807866 22.833654 4.537692 0 2.309841
2.234507 3.815925 2.125539 2.989819 2.309841 0

求得\(\beta_{ij}^0\)如表3所示

表3

0 -70.736358 51.75 -155.82017 -6.934949 -165.525066
-70.736358 0 91.128904 137.756437 87.695154 -171
51.75 91.128904 0 -167.523689 121.213757 -179.689191
-155.82017 137.756437 -167.523689 0 -174.030766 176.474394
-6.934949 87.695154 121.213757 -174.030766 0 -178.085617
-165.525066 -171 -179.689191 176.474394 -178.085617 0

上述飞行管理数学规划模型可如下输入LINGO求解:

model:
sets: ! 变量设置;
plane/1..6/:delta;
link(plane,plane):alpha,beta;
endsets
data:
alpha=@file('C:\Users\reion\OneDrive\Documents\机器学习实战\a0.csv'); ! 导入alpha_ij的数据;
beta=@file('C:\Users\reion\OneDrive\Documents\机器学习实战\b0.csv');! 导入beta_ij的数据;
enddata
min=@sum(plane:@abs(delta)); ! 最小化角度的绝对值之和;
@for(plane:@bnd(-30,delta,30)); ! 角度的调整范围在-30°到30°之间;
@for(link(i,j)|i#ne#j:@abs(beta(i,j)+0.5*delta(i)+0.5*delta(j))>alpha(i,j)); ! 调整后的角度的绝对值要小于a_ij;
end

求得的最优解为\(\Delta\theta_5=-28.05682 °\),其余角度为0°

上一篇:Linux-权限管理(你听过777、755、644吗)


下一篇:厄米特矩阵(Hermittan Matrix)