12.1一个飞行管理问题

第十二章数学建模竞赛中的部分优化问题429在中国大学生数学建模竞赛（Chinaundergraduatemathematicalcontestinmodeling,CUMCM）中，曾经出现过大量的优化建模赛题.本章从中选择了部分典型赛题，举例分析其优化建模过程，说明如何应用LINDO/LINGO软件包求解这些赛题.12.1一个飞行管理问题12.1.1问题描述1995年全国大学生数学建模竞赛中的A题（“一个飞行管理问题”）.在约10000m高空的某边长为160km的正方形区域内，经常有若干架飞机作水平飞行.区域内每架飞机的位置和速度向量均由计算机记录其数据，以便进行飞行管理.当一架欲进入该区域的飞机到达区域边缘时，记录其数据后，要立即计算并判断是否会与区域内的飞机发生碰撞.如果会碰撞，则应计算如何调整各架（包括新进入的）飞机飞行的方向角，以避免碰撞.现假定条件如下：（1）不碰撞的标准为任意两架飞机的距离大于8km；（2）飞机飞行方向角调整的幅度不应超过30°；（3）所有飞机飞行速度均为800km/h；（4）进入该区域的飞机在到达该区域边缘时，与该区域内的飞机的距离应在60km以上；（5）最多需考虑6架飞机；（6）不必考虑飞机离开此区域后的状况.请你对这个避免碰撞的飞行管理问题建立数学模型，列出计算步骤，对以下数据进行计算（方向角误差不超过0.01°），要求飞机飞行方向角调整的幅度尽量小.设该区域4个顶点的坐标为（0，0），（160，0），（160，160），（0，160）.记录数据见表12-1.表12-1飞机位置和方向角记录数据飞机编号横坐标x纵坐标y方向角)/(飞机编号横坐标x纵坐标y方向角)/(11501402434145501592858523651301502303150155220.5新进入0052说明：方向角指飞行方向与ｘ轴正向的夹角.试根据实际应用背景对你的模型进行评价和推广.12.1.2模型1及求解模型建立这个问题显然是一个优化问题.设第ｉ架飞机在调整时的方向角为（题目中已经给出），调整后的方向角为=+（=1,2,…，6）.题目中就是要求飞机飞行方向角调整的幅度尽量小，因此优化的目标函数可以是.（1）第十二章数学建模竞赛中的部分优化问题429为了建立这个问题的优化模型，只须要明确约束条件就可以了.一个简单的约束是飞机飞行方向角调整的不应超过30°，即||30°.(2)题中要求进入该区域的飞机在到达该区域边缘时，与该区域内飞机的距离应在60km以上，这个条件是个初始条件，很容易验证目前所给数据是满足的，因此本模型中可以不予考虑.剩下的关键是要满足题目中描述的任意两架飞机不碰撞的要求，即任意两架位于该区域内的飞机的距离应大于8km.但这个问题的难点在于飞机是动态的，这个约束不好直接描述，为此我们首先需要描述每架飞机的飞行轨迹.记飞机飞行速率为v（=800km/h）,以当前时刻为0时刻.设第ｉ架飞机在调整时的位置坐标为（，）（已知条件），t时刻的位置坐标为（，），则=+，=+.如果要严格表示两架位于该区域内的飞机的距离应大于8km，则需考虑每架飞机在该区域内的飞行时间的长度.记为第ｉ架飞机飞出去与的时刻，即=argmin{t0:+=0或160，（4）或者+=0或160}.（5）记t时刻第ｉ架飞机与第ｊ架飞机的距离为（t），并记（t）=-64，这时在该区域内飞机不相撞的约束条件就变成了（t）=-640（0t）.其中=min{，}.(6)此外，经过计算，可以得到（t）=+-64=++，（7）其中第十二章数学建模竞赛中的部分优化问题429=2vt，（8）=2[-(-)+,（9）=+-64.（10）所以，（t）是一个关于的二次函数，表示的是一条开口向上的抛物线.当=-/2，即t=-/4v（记为）时，函数（t）取最小值-/4+.注意到（0）=0（初始时刻不相撞）,如果-/20（即0），则此时约束条件（5）一定成立，所以对约束条件（5）只需考虑以下两种可能情况：如果0且，只需要（t）在右端点的函数值非负即可，即（）0；（11）如果0且0，只需要（t）的最小值（）=-/4+0即可，即40.(12)实际上，约束（11）表示的是（t）在右端点的函数值非负，这个约束在（12）的条件下也是自然成立的，所以可以对约束（11）不再附加0且的条件.于是，我们的模型就是;（13）S．t．||30°，（14）（）0(15)40.（当0且0）.(16)第十二章数学建模竞赛中的部分优化问题429模型求解上面这是一个非线性规划模型，虽然是严格满足题目要求的模型，但得到的模型逻辑关系比较复杂，约束（16）是在一定条件下才成立的约束，而且其中的计算式（4）也是含有相当复杂的关系式，使用LINGO软件不太容易将模型很方便地输入，因为逻辑处理不是LINGO的优势所在.即使能想办法把把这个模型输入到LINGO,也不一定能求出好的解(笔者尝试过,但LINGO运行时有时会出现系统错误,可能是系统有漏洞,无法继续求解).而且,在实时飞行调度中显然需要快速求解,所以下面我们想办法简化模型.这个模型麻烦之处就在于,要求严格表示两架位于该区域内的飞机距离应大于8km,所以需要考虑每架飞机在该区域内的飞行时间,比较繁琐.注意到区域对角线的长度只有160km，任何一架飞机在所考虑的范围内停留的时间不会超过=160/800=0.20.283（h），因此这里我们简化一下问题：不再单独考虑每架飞机在该区域停留的时间，而以最大时间（这里已经是一个常数）代替之，此时所有=.这实际上强化了问题的要求,即考虑了有些飞机可能已经飞出该区域,但仍不允许两架飞机的距离小于8km.这个简化的模型可以如下输入LINGO软件：MODEL:TITLE飞行管理问题的非线性规划模型;SETS:Plane/1..6/:x0,y0,cita0,cita1,d_cita;!cita0表示初始角度，cita1为调整后的角度,d_cita为调整的角度;link(plane,plane)|&1#LT#&2:b,c;ENDSETSDATA:x0y0cita0=1501402438585236150155220.5145501591301502300052;max_cita=30;T_max=0.283;V=800;ENDDATAINIT:d_cita=000000;第十二章数学建模竞赛中的部分优化问题429ENDINIT@for(plane:cita1-cita0=d_cita);@for(link(i,j):b(i,j)=-2*(x0(i)-x0(j))*@sin((cita1(i)+cita1(j))*3.14159265/360)+2*(y0(i)-y0(j))*@cos((cita1(i)+cita1(j))*3.14159265/360);c(i,j)=(x0(i)-x0(j))^2+(y0(i)-y0(j))^2-64;);!避免碰撞的条件;!右端点非负;@for(link(i,j):[Right](2*V*T_max*@sin((cita1(i)-cita1(j))*3.14159265/360))^2+b(i,j)*(2*V*T_max*@sin((cita1(i)-cita1(j))*3.14159265/360))+c(i,j)0);!最小点非负;@for(link(i,j):[Minimum]@if(b(i,j)#lt#0#and#-b(i,j)/4/V/@sin((cita1(i)-cita1(j))*3.14159265/360)#gt#0#and#-b(i,j)/4/V/@sin((cita1(i)-cita1(j))*3.14159265/360)#lt#T_max,b(i,j)^2-4*c(i,j),-1)0);!@for(link(i,j):@if(b(i,j)#lt#0,b(i,j)^2-4*c(i,j),-1)0);@for(link:@free(b));!调整角度上下限，单位为角度;@for(plane:@bnd(-max_cita,d_cita,max_cita));[obj]MIN=@SUM(plane:(d_cita)^2);![obj]MIN=@SUM(plane:@abs(d_cita));END注意：上面的模型中的方向角单位一律用角度，而LINGO只接受弧度，所以程序中一律进行了转换.求解这个模型,得到Localoptimalsolutionfound.Objectivevalue：295.4937ModelTitle:飞行管理问题的非线性规划模型VariableValueReducedCostD_CITA(1)-10.5.9800.000000D_CITA(2)0.0000000.000000D_CITA(3)0.0000000.000000D_CITA(4)6.5154250.000000D_CITA(5)10.006810.000000D_CITA(6)6.5054250.000000这个结果得到的是一个局部极小点，调整角度较大.能找到更好的解吗?如果不用全局求解程序,通常很难得到稍大规模的非线性规划问题的全局最优解.所以我们启动LINGO全局求解程序求解这个模型,可以得到全局最优解如下:Globaloptimalsolutionfoundatiteration:93第十二章数学建模竞赛中的部分优化问题429Objectivevalue:6.953944ModelTitle:飞行管理问题的非线性规划模型VariableValueD_CITA(1)0.2719480E-02D_CITA(2)0.5613433E-02D_CITA(3)2.059140D_CITA(4)-0.4985421D_CITA(5)-0.5407837E-03D_CITA(6)1.570129(这里只给出的值)可以看到，在0.01°的误差要求下，需要调整第3，4，6三架飞机的角度，分别调整2.06°，-0.50°，1.57°.调整量的平方和为6.95.其实,使用全局求解程序,通常也不一定要等到得到全局最优解,而是观察求解状态窗口,看到一个较好的当前解(或当前最好解在较长时间内不发生变化)时,就可以终止程序,用当前最好的局部最优解作为最后的恶结果.列如,对于本列,LINGO求出全局最优解大约需要1min,而实际上5s内LINGO就得到了与全局最优解类似的解.此外,上面的模型还可以进一步简化,列如可以假设要求飞机永远不相撞,即认为为无限大，这时显然约束（15）也是多余的,而且约束(16)中只需要0的条件就可以了.也就是说，上面的程序中的对应部分（约束[Right]和[Minimum]）可以改写为更简单的形式：！有端点非负，不再需要；！最小点非负，简化为以下形式；@for(link(i,j):@if(b（i，j）#lt#0，b（i，j）^2-4*c（i，j），-1）0）；实际计算显示，此时得到的结果与前面计算的结果几乎没有差别.备注优化的目标函数除了外，也可以设定为或等，用LINGO求解的过程是完全类似的，计算结果略有差异，这里就不再对两个目标函数具体计算了.甚至可以考虑让参与调整的飞机的数量尽量小，这种想法在实际中也不能说没有道理，但与题目中的要求不符，而且解题难度并没有减小，意义似乎不大.实际调度中，由于计算上面的调度方案需要时间，将调度信息告诉飞机驾驶员并作出调整方向角的操作也需要时间，因此如果考虑一定的反应滞后时间，应该是比较合理的.也就是说，如果反应时间是10s，则计算式中应采用飞机沿当前的方向角飞行10s以后的位置作为计算的基础.12.1.3模型2及求解从12.1.2节可以看出，求解模型1的非线性规划模型是比较困难的，输入后也很可能找不到好的解甚至出现错误.此外，演示版软件还会受到求解规模的限制，尤其可能无权使用全局求解程序.因此，如果能把这个问题简化成比较简第十二章数学建模竞赛中的部分优化问题429单的规划模型，将是非常有价值的.模型建立如图12-1，把两架飞机