数学建模 飞机管理优化模型.docx

数学建模 飞机管理优化模型飞行管理优化模型 摘 要 本文建立了关于飞行管理问题的简洁数学模型。首先我们对新进入的飞机作出判断，通过模型给出了计算机模拟求解，看其是否会与其他六架飞机相碰，若会，则再次通过建型求出使各飞机安全通过区域应调整的方向角，对模型给出了非线性优化的具体算法。在模型改进中，我们对风速、人的反应时间及飞机的实际速度等对方向角的调整的影响做了简单的分析与评价，使得模型更易在实际应用中推广。 模型：针对问题一，建立了碰撞检测模型。首先，对已给数据进行分析，并利用VC编程，模拟在6驾飞机都不改变飞行方向的条件下的飞行情况，结果是会相撞的。 模型：针对问题二，建立了多元非线性动态优化模型。在确保互不相撞的前提下，要使得飞机调整的角度尽可能小，满足Min f=åi=16bi2，及最优解。运用MATLAB软件编程给出了具体算法。第i架飞机需调整角度分别为：0.0000 0.0000 2.0683 -0.4896 -0.0055 1.5611 关键词：飞机碰撞 方向角调整 非线性优化 模拟仿真 一、问题重述 在约10000米的高空某边长为160公里的正方形区域内，经常有若干架飞机作水平飞行。区域内每架飞机的位置和速度向量均由计算机记录其数据，以便进行飞行管理。当一驾欲进入该区域的飞机到达区域边缘时，记录其数据后，要立即计算并判断是否会与区域内的飞机发生碰撞。如果会碰撞，则应计算如何调整各架飞机飞行的方向角，以避免碰撞。现假定条件如下： 1）不碰撞的标准为任意两架飞机的距离大于8公里； 2）飞机飞行方向角调整的幅度不应超过30度； 3）所有飞机飞行速度均为每小时800公里； 4）进入该区域的飞机在到达区域边缘时，与区域内飞机的距离应在60公里以上； 5）最多需考虑6架飞机； 6）不必考虑飞机离开此区域后的状况。 1 请你对这个避免碰撞的飞机管理问题建立数学模型，列出计算步骤，对以下数据进行计算，要求飞机飞行方向角调整的幅度尽量小。 设该区域内4个顶点的坐标为(0,0),(160,0),(160,160),(0,160)。 记录数据为： 飞机编号 1 2 3 4 5 新进入 横坐标x 150 85 150 145 130 0 纵坐标y 40 85 155 50 150 0 方向角 243 236 220.5 159 230 52 注：方向角指飞行方向与x轴方向的夹角。 试根据实际应用背景对你的模型进行评价与推广。 问题要求 判断飞机是否会碰撞； 若会碰撞，请给出方向角最优调整方案。 二、模型假设 1、不碰撞的标准为任意两架飞机的距离大于8公里； 2、飞机飞行方向角调整的幅度不超过30度； 3、所有飞机在同一平面飞行，且速度均为每小时800公里； 4、当飞机刚进入区域时就开始调整方向角，全过程只调整一次，且调整时间忽略不计； 5、飞机在区域内飞行时不考虑其他一切外界因素的影响； 6、相对正方形区域而言，将飞机作为质点考虑； 7、最多需考虑6架飞机，不必考虑飞机离开此区域后的状况。 2 三、符号说明 : 第i架飞机的坐标单位：km； ：第j架飞机的坐标单位：km； v : 飞机飞行速度，其值为800 km/h； ti: 第i架飞机在区域内飞行时间 单位：h； t : 飞机在区域内飞行的任一时刻 单位：h； ai、aj : 第i、j架飞机的方向角 单位：°； bi、bj : 第i、j架飞机需调整的方向角 单位：°； dij: 表示第i架飞机与第j架飞机之间的距离，单位：km。 四、问题分析 发生碰撞 飞机安全通过正方形区域 调整方向角 不 碰判断是否碰撞 由上图可知：模型假设飞机在同一高度飞行，所以只需讨论飞机在二维平面的碰撞及调节管理问题。以实现进入该区域的飞机都能安全通过该区域的目标。 当新飞入的飞机刚进入该正方形区域时，首先应通过模型对其作出判断，看其在区域内飞行时是否会与区域内的其他飞机发生碰撞。并在二维坐标系中标出各飞机的位置坐标，结合其方向角借助计算机进行仿真模拟。 若通过模型验证不会碰撞，则使各飞机按原路线飞行即可安全通过正方形区域；否则，就需通过模型计算得出各飞机需要调整的角度，并且使各角度的3 改变量尽可能小，即最优解。在原方向角的基础之上调整计算得出的方向角，再用模型对其进行检验，提高结果的可信度，以确保各飞机能安全通过该正方形区域。 另外，为了使结果更加让人信服，在对时间的范围选取和最优解的处理问题上，我们进行了细致的分析。 时间最小：ti表示为第i架飞机在区域内飞行的时间，tj表示为第j架飞机在区域内飞行的时间，取t=min，即计算第i架飞机与第j架飞机间的距离时保证两架飞机都在区域内,且在t时间内第i架飞机不与第j架飞机相撞，即保证所讨论范围在正方形区域内。 平方和最小：由于所调角度bi，如果直接求和，正负将会抵消，所以目标函数定为：Min f=åi=16bi2。 五、模型的建立与求解 模型一：碰撞检测模型 在平面直角坐标系中设出各飞机的坐标分别为(xi,yi)，方向角分别为ai(i=1,2,3,4,5,6)，则第i架飞机t时刻的位置为： (xi+vtcosai,yi+vtsinai) 第i架飞机在区域内时满足： 0£xi+vtcosai£1600£yi+vtsinai£160 t³0求得0£t£ti,ti即为第i架飞机在区域内的飞行时间。 所以，第6架飞机在区域内飞行的时间为t6。 在t=min时间内，第6架飞机与第i架飞机都在正方形区域内， 第六架飞机与其他五架飞机的距离为 di6=xi+vtcosai-(x6+vtcosa6)2+yi+vtsinai-(y6+vtsina6)2>8 4 如果满足 di6即认为不会碰撞。 根据题目所给数据，标出各飞机在二维坐标中的位置。其中左上角为坐标原点，横轴水平向右为x轴正方向，纵轴竖直向下为y轴正方向。1、2、3、4、5、6分别为飞机编号，红色圆圈表示两架飞机互不相撞应满足的范围，即圆的半径为4km。 图5.1 飞机6进入该区域时各飞机位置 在VC模拟环境下对其进行判断，根据程序模拟结果的截图可直观的看出新进入的编号为6的飞机会与编号为5的飞机和编号为3的飞机碰撞： 图5.2 飞机5和飞机6碰撞图 5 表5.1：飞机5和飞机6碰撞时各架飞机间距 飞机编号 1 2 3 4 5 6 1 / / / / / / 2 89.7594 / / / / / 3 48.6257 98.5514 / / / / 4 60.343 83.7181 22.9526 / / / 5 44.2661 80.1122 18.4395 16.1179 / / 6 51.5642 77.686 22.4788 9.9093 7.8443 / 图5.3飞机3和飞机6碰撞图 表5.2：飞机3和飞机6碰撞时各架飞机间距 飞机编号 1 2 3 4 5 6 1 / / / / / / 2 90.1614 / / / / / 3 51.4065 99.0578 / / / / 4 68.3697 91.4379 23.0886 / / / 5 45.9234 80.2436 18.8310 23.2464 / / 6 52.1156 91.9763 7.9422 17.5740 12.3961 / 根据以上判断，新进入的飞机会与区域内的飞机发生碰撞，所以需对各飞机的方向角进行调整。 6 模型二：多元非线性优化模型 经过模型的判断，如果会碰撞，则需计算如何调整各架飞机的方向角，以避免碰撞，并且使各架飞机调整的角度尽可能小，即求一个最优解。设各飞机飞行方向的调整角度为约束变量，建立模型，由于bi，如果直接求和，正负将会抵消，所以目标函数定为： Min f=åi=16bi2约束条件： 按调整后的角度飞行，任意两架飞机在区域内的距离大于8公里，即 dij=xi+vtcos(ai+bi)-xj+vtcos(aj+bj)2+yi+vtsin(ai+bi)-yj+vtsin(aj+bj)2 dij>64 (i<j,i,j=1,2,3,4,5,6)； 根据模型知，第i、j架飞机在区域内飞行的时间分别为：ti、tj； 在t=min时间内，第i架飞机与第j架飞机都在正方形区域内。 飞机飞行最大调整角度bi、bj满足：故飞机调整角度最优解模型如下： Min f=åi=16-30o£bi£30o-30£bj£30oobi2dijf64 s.t -30o£bi£30o -30o£bj£30o 这是一个非线性优化问题，借助MATLAB的fimincon函数编程计算得出最优解。 表5.3：各架飞机方向角调整情况表 飞机编号i 调整前的方向角ai 调整的方向角bi 调整后方向角ai1 243 0 243 2 236 0 236 3 220.5 4 159 5 230 6 52 2.0683 -0.4896 -0.0055 1.5611 222.5683 7 +bi158.5104 229.9945 53.5611 由所求结果，利用模型一对结果进行检验，用VC模拟调整角度后的飞机的动态飞行路线 ，检验结果的截图见图5.4。 图5.4 调整角度后模拟飞行路线截图 通过图5.4与图5.3比较，各飞机均可安全飞出正方形区域，完全符合安全飞行标准。 六、模型的灵敏度分析 由于飞机的内部性能和飞行员的技术水平等一系列内部因素及不断变化的外部条件都可能对结果产生影响，因此讨论这些因素对最优解的影响是十分必要的。 Dt表示新飞如的飞机进入区域Dt后各飞机开始调整方向。则 d=ijìxíiî+Dtcosa+vtcosai+bii()ë-éxêëj+Dtcosa+vtcosaj+bjj(ü)ùýúûþ2+ìíîyi+Dtsina+vtsinai+bii()-éyêj+Dtsina+vtsinaj+bjj(ü)ùýúûþ2dijf64 s.t -30o£bi£30o -30o£bj£30o 由以上条件求得最优解。显然，Dt越大，模型无解的可能性越大，即各8 飞机在新飞入的飞机进入区域后越晚调整，飞机在正方形区域内发生碰撞的概率越大。 v0表示区域内平均风速，其方向角为a0。所以此时飞机飞行的水平速度为vcosai+v0cosa0，竖直速度为vsinai+v0sina0，即飞机此时沿曲线飞行。风速与飞机飞行方向的夹角越大，飞机的操纵性越差，即飞机在侧风中飞行时，发生安全事故的可能性越大。 七、模型的评价 1、优点 本模型成功解决了飞行管理问题，建立了较为满意的模型，并且应用了仿真模拟系统很好地检验了结果的正确性，可信度较高； 模型建立时在时间t的处理问题上，精确地将t限定在区域内，即只考虑飞机在区域内的碰撞情况，所以能够保证取得最优解； 本模型的约束条件较少，计算本应很复杂，但使用了Matlab中的fmincon函数使得计算大大简化，且计算结果够精确。 2、缺点 模型是在理想条件下进行的，不考虑天气、风向等一切因素的影响，且设定飞行员接受命令到执行的时间为零，飞行调整的角度为绝对精确，因此可能与实际情况有一定差距； 考虑飞机是在同一平面内飞行，未考虑飞机间的垂直距离； 由于在计算过程中题目给的飞机飞行速度800公里/小时，而现实中不可能每架飞机在每一时刻的速度都为800km/h，所以使得模型计算的结果与实际情况有所差别。 八、模型的改进与推广 本模型是建立在多个假设条件之下的，排除各类影响因素的情况下，是比较简洁和明晰的，为更深入地研究提供了借鉴基础。 1.在现实生活中，各架飞机的速度不完全一样，可以假设六架飞机的速度分别为vi(i=1,2,3,4,5,6)，那么六架飞机之间的相互距离为： 9 dij=ìxíiî+vtcosai+bii()-éxêëj+vtcosaj+bjj(ü)ùýúûþ22ì+íîyi+vtsinai+bii()-éyêëj+vtsinaj+j(übj)ùýúûþ(i<j,i=1,2,3,4,5,6) 则利用已建立模型的原理同样可以使问题得到解决。 2. 把飞行区域从二维引申到三维空间，更切合了生活实际。一定程度上增加了解题的复杂度，但中心思想应该是不变的，通过两点之间距离的限定作为重要的约束条件，增加垂直方向坐标，在三维空间里做进一步的讨论。 3.当把飞机限定在6驾时，其实讨论的是该区域的极限，当区域内飞机的数量少于6驾时，数据分析会更简单一些。根据情况可以适当缩减改变角的范围，提高运行的效率。 4.可以将风速、飞行员的反应时间及飞机调整方向角所需时间等因素考虑进去，然后做更进一步的讨论与分析。 本题的求最优解模型不仅仅适用于在正方形区域内六架飞机的碰撞判断及方向角调整的问题，对于其他的两点之间距离判断及如何做出优化调整方案都具有积极的指导意义。如在航海运输方面设计轮船的最优航行路线等。总之，凡是涉及到要求设计最优化路线等的问题都可以应用到本模型来解决。 八、参考文献 1张兆宁，张晓燕，交叉航路碰撞风险研究，航空计算技术报，第37卷，第二期，XX年三月，第14页。 2姜启源、谢金星、叶俊，数学模型 第三版, 北京；高等教育出版社，2003。 3黄忠裕，初等数学建模, 成都：四川大学出版社，2004.12。 4应爱玲、徐肖豪, 空域飞行侧向碰撞危险的Reich，模型方法研究中国民航学院学报，第20卷第四期，第610页。 5宫文娟，利用代数方法求解碰撞检测问题，山东大学硕士学位论文，绪论。 10 附录 飞机对应的结构体： typedef struct PLANE_TYPE char nIndex4;/飞机编号 bool bIsIn; /判断当前该飞机是否在该安全区域中 double x; /飞机横坐标 double y; /飞机纵坐标 double angle; /飞机相对于X轴方向角 PLANE; *本程序实现的主要函数* void CSimulationDlg:MoveToShow(PLANE &pe,IplImage* img) /计算经历时间Tstep后飞机的位置 pe.x = pe.x+V*Tstep*cos(pe.angle*PI/180)/3600; pe.y = pe.y+V*Tstep*sin(pe.angle*PI/180)/3600; if(pe.x <= 0 | pe.x >= WIDE | pe.y <= 0 | pe.y >= WIDE) /判断飞机是否在安全区中 pe.bIsIn = false; return; int x = pe.x*ROT; int y = pe.y*ROT; cvCircle(img,cvPoint(x,y),4*ROT,CV_RGB(255,0,0),1,8,0);/标出每架飞机的圆形区域 cvPutText(img,pe.nIndex,cvPoint(x-3*ROT,y+2*ROT), &m_font, cvScalar(0,255,255,0); /标出飞机编号 / bool CSimulationDlg:IsCrash(int &centerX,int &centerY,int& p1,int& p2) /判断安全区中是否有飞机否相撞 double s = 0.0; centerX = -1; centerY = -1; for(int i=0;i<NUMofPLANE;i+) if(!m_arrPlanei.bIsIn) continue; for(int j=i+1;j<NUMofPLANE;j+) if(!m_arrPlanej.bIsIn) continue; s = pow(m_arrPlanei.x-m_arrPlanej.x,2); s += pow(m_arrPlanei.y-m_arrPlanej.y,2); 11 if(s < 64) centerX = (m_arrPlanei.x+m_arrPlanej.x); centerY = (m_arrPlanei.y+m_arrPlanej.y); p1=i; p2=j; return true; return false; / IplImage* CSimulationDlg:NowToShow(bool &b,char* &info) /将所有飞机显示出来 IplImage* img; int m = NUMofPLANE; info = NULL; b = true; img = GetBackground("base.jpg"); /将移动后的每一架飞机显示在图片上 for(int i=0;i<NUMofPLANE;i+) if(m_arrPlanei.bIsIn) MoveToShow(m_arrPlanei,img); int centerX,centerY; int p1,p2; if(IsCrash(centerX,centerY,p1,p2)/判断飞机是否相撞 /飞机如果相撞则用蓝色圆将相撞的两架飞机圈出 cvCircle(img,cvPoint(centerX,centerY),8*ROT,CV_RGB(0,255,0),1,8,0); /保存相撞时的图片 cvSaveImage("crash.jpg",img); b = false; info = new char32; /保留信息 sprintf(info,"After %.3fh, %d and %d crash",m_dElapse*V/3600/800,p1,p2); KillTimer(0);/终止定时器 return img; for(i=0;i<NUMofPLANE;i+)/用于判断安全区域中还是否有飞机 if(!m_arrPlanei.bIsIn) m-; 12 if(m = 0) /如果安全区中没有飞机 info = new char32; sprintf(info,"After %.3fh,No plane in the area",m_dElapse*V/3600/800); KillTimer(0); return img; 注：本程序在vc 6.0 MFC中实现，借助opencv包完成飞机管理模型的动态模拟。以上三个函数为本程序的主要函数。本程序中，利用SetTimer设置定时器，时间步长为0.01秒 13 附录 function f=result %¸根据目标函数和约束函数，利用fmincon求解 deltaini=zeros(1,6); ¸ vlb=-pi*ones(1,6)/6;vub=pi*ones(1,6)/6; options=optimset('LargeScale','off','Algorithm' ,'active-set');%,'Tolx',1.0000e-032); dt,fval=fmincon('airfun',deltaini,vlb,vub,'airfunco',options); d1=dt*180/pi,d2=d1+243 236 220.5 159 230 52,fval=d1*d1' / function f=airfun(delta) %目标函数:mindelta2 f=delta*delta' / function c,ceq=airfun(delta) %约束函数，要求任意两架飞机在安全区中的的最短距离的平方都不小于64 x0=150 85 150 145 130 0;y0=140,85,155,50,150,0; alpha0=243 236 220.5 159 230 52*pi/180;v=800; co=cos(alpha0+delta);si=sin(alpha0+delta); for i=1:6 for j=i+1:6 d(i,j)=getmin(x0(i),v*co(i),x0(j),v*co(j),y0(i),v*si(i),y0(j),v*si(j); end end c=64-d(1,2:6),d(2,3:6),d(3,4:6),d(4,5:6),d(5,6); 14 ceq=; / function f=getmin(xi,coi,xj,coj,yi,sii,yj,sij) %求第i架飞机和第j架飞机在安全区中所相距的最短距离 %设第i架飞机在安全区中飞行时间为ti %设第j架飞机在安全区中飞行时间为tj %则tUp=min(ti,tj) %在tUp时间内，只要i和j两架飞机中有一架通过，则i和j不会发生碰撞 v=800; tUp=min(gettime(xi,coi/v,yi,sii/v),gettime(xj,coj/v,yj,sij/v); banana=(t)(xi+t*coi)-(xj+t*coj)2+(yi+t*sii)-(yj+t*sij)2; x,f = fminbnd(banana,0,tUp); / function f=gettime(x0,co,y0,si) %求飞机在安全区中飞行的时间 v=800; f=min(max(-x0/v/co,(160-x0)/v/co),max(-y0/v/si,(160-y0)/v/si); 15