数学建模论文嫦娥三号软着陆轨道设计与控制策略

1、.2014高教社杯全国大学生数学建模竞赛承 诺 书我们仔细阅读了全国大学生数学建模竞赛章程和全国大学生数学建模竞赛参赛规则（以下简称为“竞赛章程和参赛规则”，可从全国大学生数学建模竞赛网站下载）。我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。我们知道，抄袭别人的成果是违反竞赛章程和参赛规则的，如果引用别人的成果或其他公开的资料（包括网上查到的资料），必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出。我们郑重承诺，严格遵守竞赛章程和参赛规则，以保证竞赛的公正、公平性。如有违反竞赛章程和参

2、赛规则的行为，我们将受到严肃处理。我们授权全国大学生数学建模竞赛组委会，可将我们的论文以任何形式进行公开展示（包括进行网上公示，在书籍、期刊和其他媒体进行正式或非正式发表等）。我们参赛选择的题号是（从A/B/C/D中选择一项填写）： A 我们的报名参赛队号为（8位数字组成的编号）： 18004050 所属学校（请填写完整的全名）： 中南大学 参赛队员 (打印并签名) ：1. 董建洋 2. 郑振鹏 3. 蒋源 指导教师或指导教师组负责人 (打印并签名)： 侯木舟 （论文纸质版与电子版中的以上信息必须一致，只是电子版中无需签名。以上内容请仔细核对，提交后将不再允许做任何修改。如填写错误，论文可能被

3、取消评奖资格。） 日期： 2014 年 9月 14 日赛区评阅编号（由赛区组委会评阅前进行编号）：2014高教社杯全国大学生数学建模竞赛编 号 专 用 页赛区评阅编号（由赛区组委会评阅前进行编号）：赛区评阅记录（可供赛区评阅时使用）：评阅人评分备注全国统一编号（由赛区组委会送交全国前编号）：全国评阅编号（由全国组委会评阅前进行编号）：嫦娥三号软着陆轨道设计与控制策略摘要针对嫦娥三号软着陆轨道设计与控制策略问题，本文通过建立燃料最优模型、重力转换模型、统计学模型、垂直软着陆模型，对嫦娥三号软着陆的运行轨迹进行讨论，对其进行了有效的设计与控制。针对问题一和问题二，本文根据角动量守恒定律和机械能守恒

4、定律，首先计算出嫦娥三号在近月点的速度大小为v1=1693.036m/s，远月点的速度大小为V2=1614.648m/s，近月点和远月点的速度相反，并且都垂直于近月点、远月点和球心的连线。然后，我们建立燃料最优模型，并进行一定的假设，得到嫦娥三号在主减速阶段的最优轨迹。根据得到的最优轨迹方程，利用已知条件，得出软着陆过程中的水平位移为458.875km，进一步得出近月点在29.02N,19.51W，高度15Km处和远月点在29.02S,159.49E,高度100Km处。在快速调整阶段，本文建立了重力转换模型，推导出偏航角控制量和俯仰角控制量的表达式，实现发动机的推力方向与嫦娥三号速度矢量的反方

5、向保持一致。在粗避障区和静避障区阶段，利用matlabR2008a读出图形数据，建立了统计学模型，利用层次分析的思想，先确定一定的相对平坦的区域作为着陆范围，再在所选的区域中找出精确的着陆点。对于缓速下降和自由落体阶段，本文建立垂直软着陆模型，采用分段控制方式，把若干推力级别等效为一段连续的推力过程。嫦娥三号依次经过悬停、匀加速、匀减速和关机降落几个过程。 针对问题三，由于主减速阶段终端时的水平速度无法控制，进而对主减速阶段的运行轨迹产生影响。本文利用matlab画出剩余时间、偏航角和水平位移关于不同水平速度的图像。得到水平位移在小的范围内变化所造成的纬度偏角的误差很小。在精避障阶段，确定着陆

6、范围时，所选处理矩阵的大小直接影响探测器能否准确、快速的着陆。所选区域越大，能找到的最平坦的区域误差就会越小，但数据处理的时间越长。关键字：嫦娥三号 燃料最优 重力转弯 统计学 垂直软着陆 分段控制一、问题重述嫦娥三号于2013年12月2日1时30分成功发射，12月6日抵达月球轨道。嫦娥三号在着陆准备轨道上的运行质量为2.4t，其安装在下部的主减速发动机能够产生1500N到7500N的可调节推力，其比冲（即单位质量的推进剂产生的推力）为2940m/s，可以满足调整速度的控制要求。在四周安装有姿态调整发动机，在给定主减速发动机的推力方向后，能够自动通过多个发动机的脉冲组合实现各种姿态的调整控制。

7、嫦娥三号的预定着陆点为19.51W，44.12N，海拔为-2641m（见附件1）。嫦娥三号在高速飞行的情况下，要保证准确地在月球预定区域内实现软着陆，关键问题是着陆轨道与控制策略的设计。其着陆轨道设计的基本要求：着陆准备轨道为近月点15km，远月点100km的椭圆形轨道；着陆轨道为从近月点至着陆点，其软着陆过程共分为6个阶段（见附件2），要求满足每个阶段在关键点所处的状态；尽量减少软着陆过程的燃料消耗。根据上述的基本要求，请你们建立数学模型解决下面的问题：（1） 确定着陆准备轨道近月点和远月点的位置，以及嫦娥三号相应速度的大小与方向。（2）确定嫦娥三号的着陆轨道和在6个阶段的最优控制策略。（3

8、）对于你们设计的着陆轨道和控制策略做相应的误差分析和敏感性分析二、模型假设(1）地球和月球对嫦娥三号的吸引力可以忽略。（2）在着陆过程中月球的自转可以忽略。(3)嫦娥三号的着陆点服从正态分布。(4)嫦娥三号在着陆中不会出现故障。(5)月球是一个理想球体，不考虑月球非球形重力场和月球偏率的间接效应。(6)在着陆过程，嫦娥三号的质量不变。三、符号说明符号符号的意义M月球质量m嫦娥三号质量G月球表面万有引力常数V1近月点速度V2远月点速度vxo近月点时轴方向速度vyo近月点时轴方向速度vzo近月点时轴方向速度vxf终端时刻轴方向速度vyf终端时刻轴方向速度vzf终端时刻轴方向速度xo近月点轴坐标yo

9、近月点轴坐标zo近月点轴坐标xf终端时刻轴坐标yf终端时刻轴坐标zy终端时刻轴坐标trem剩余时间P主发动机推力af主发动机推力产生的加速度*偏航角*俯仰角gm月心重力加速度纬度偏角a0等效推力加速度偏移量K等效推力加速度系数四、问题分析问题一要求我们确定着陆准备轨道近月点和远月点的位置，以及嫦娥三号相应速度的大小与方向。我们认为根据角动量公式和能量守恒即可算出近月点和远月点的速度。我们认为根据主减速阶段的水平位移以及近月点高度即可求出近月点和远月点的位置。问题二要求我们确定嫦娥三号的着陆轨道和在6个阶段的最优控制策略。首先我们建立燃料次优显示制导模型，由于主减速阶段燃料消耗约占全部燃料消耗的

10、75%，该模型采用极大值原理，使主减速阶段的燃料消耗最小，并进行一定的假设，得到嫦娥三号在主减速阶段的最优轨迹。再次我们建立滑模变结构控制模型，推导出偏航角控制量和俯仰角控制量的表达式，实现发动机的推力方向与嫦娥三号速度矢量的反方向保持一致。在粗避障区和静避障区阶段，利用matlabR2008a读出图形数据，建立了统计学模型，利用层次分析的思想，先确定一定的相对平坦的区域作为着陆范围，再在所选的区域中找出精确的着陆点。对于缓速下降和自由落体阶段，本文建立垂直软着陆模型，采用分段控制方式，把若干推力级别等效为一段连续的推力过程。嫦娥三号依次经过悬停、匀加速、匀减速和关机降落几个过程。问题三要求我

11、们对设计的着陆轨道和控制策略做相应的误差分析和敏感性分析。由于主减速阶段终端时的水平速度无法控制，进而对主减速阶段的运行轨迹产生影响。本文利用matlab2008a画出剩余时间、偏航角和水平位移关于不同水平速度的图像。得到水平位移在小的范围内变化所造成的纬度偏角的误差很小。在精避障阶段，确定着陆范围时，所选处理矩阵的大小直接影响探测器能否准确、快速的着陆。所选区域越大，能找到的最平坦的区域误差就会越小，但数据处理的时间五、模型的建立与求解5.1近月点和远月点速度的求解 我们知道嫦娥三号着陆准备阶段时在近月点15km，远月点200km的椭圆轨道上运动，此时嫦娥三号完全依靠月球对它的引力做椭圆运动

12、。我们可知在运动中嫦娥三号的能量守恒，且满足椭圆运动的角动量定律。即可列出下列公式：v1r1=v2r2 （5-1）-GMm/r1+1/2mv12=-GMm/r2+1/2mv22 （5-2） 带入数据即可求得:v1=1693.036m/s v2=1614.648m/s (求解见附件程序Prog1)而近月点和远月点的速度相反，并且都垂直于近月点、远月点和球心的连线。5.2模型一的建立与求解5.2.1制动软着陆坐标系O0-X0Y0Z0，简记为O 该坐标系的原点是嫦娥三号质心在软着陆制动时刻的月心矢与标准月球椭球体表面的交点O0，O0Y0轴沿O0点与制动时刻嫦娥三号质心O1连线方向，O0X0在软着陆制

13、动时刻月球嫦娥三号运行的轨道平面内，且与O0Y0相互垂直，并指向嫦娥三号的运动方向。O0Z0轴与O0X0、O0Y0轴构成右手坐标系。利用该坐标系建立嫦娥三号软着陆过程中相对于月面的运动方程。5.2.2嫦娥三号本体坐标系O1-X1Y1Z1。该坐标系原点为嫦娥三号的质心O1，O1X1为嫦娥三号纵轴，指向头部。O1Y1在嫦娥三号主对称面内，该平面在软着陆制动时刻与制动软着陆坐标系X0O0YO平面重合，并且垂直O1X1，O1Z1轴垂直于主对称面，沿运动方向看去，O1Z1指向右方。该坐标系即我们熟悉的本体坐标系。5.2.3模型一的建立 由于嫦娥三号在主减速阶段的燃料消耗约占全部燃料消耗的75%1，因此要

14、选择一个合理的制导控制方案，使嫦娥三号在主减速阶段的燃料消耗最小。但是由于嫦娥三号是一个强非线性系统，直接求解最优控制问题比较复杂，且难以获得解析形式的表达式，通常需要利用给定初值进行迭代的方法进行计算，不利于嫦娥三号实现自主控制。因此有必要在合理假设的基础上进行简化处理：（1） 只考虑嫦娥三号的质心动力学方程和质心运动学方程，不考虑其绕质心转动动力学等方程。（2） 忽略月球自转，不考虑离心惯性力和哥氏惯性力。(3) 假设主发动机为理想安装，认为主发动机推力倾差角和推力偏差角均为0。并忽略姿态发动机推力的影响。因此，我们得到的嫦娥三号简化模型为： dvx dt=Pmcoscos dvydt=P

15、msincos-fMr2 dvzdt=-Pmsin dxdt=vx (5-3)dydt=vy dzzt=vz 5.2.4模型一的求解 我们利用极大值原理来求解最小燃料消耗。在主减速阶段开始时，主发动机开始点火。此时，制导系统的主要任务是抵消水平的初始速度，因此可以认为航偏角为180度，俯仰角为小量。可设cos=1, sin=k1+k2t带入式(5-3)得到： dvxdt=Pmcos0 dvydt=Pmk1+k2tcos0-fMr2 dvzdt=-Pmsin0 （5-4） dxdt=vt dydt=vy dzdt=vt求解得到嫦娥三号在主减速阶段的最优轨迹： vxf=aFcos0trem+vx0

16、 v yf=aFk1cos0trem+aFk2cos0(trem2)2+vy0 vzf=-aFsin0trem+vz0 ( 5-5） xf=aFcos0(trem2)2+vx0trem+x0 yf=(aFk1cos0-fMr2)(trem2)2+aFk2cos0(trem3)6+vy0trem+y0 zf=-aFsin0(trem2)2+vz0trem+z0其中：aF=Pm表示主发动机推力作用在嫦娥三号上产生的加速度。 trem表示剩余时间。有（5-5）的（1）、（3）可以推导出剩余时间的表达式：trem=(vxf-vx)2+(vzf-vz)2aF （5-6） 其中： cos0=vxf-vx(

17、vxf-vx)2+(vzf-vz)2 （5-7）r=x2+(y+R0)2+z2 （5-8）由上面近月点和远月点的速度求解可知vx0=1693m/s，vy0=0，vz0=0.由于主减速阶段结束时嫦娥三号已基本位于着陆点上方，可知此时嫦娥三号的水平速度很小，不妨假设vxf=5m/s，vyf=0，由此时嫦娥三号的速度为56m/s可得vzf=56.78m/s。将以上速度带入（5-6）、（5-7）可得： 剩余时间 trem=542s cos0=-0.9994再将trem、 cos0带入式（5-5）中的（4）可得： xf=458.875kmxf即为嫦娥三号在减速阶段的前进距离，有嫦娥三号的前进距离和起点终

18、点的高度，根据图1，在图中三角形中即可求得：图1cos=0.9655,然后求得主减速阶段的纬度偏角：=15.10°由文献上查得嫦娥三号是从纬度低的位置减速到纬度高的位置，由于在减速阶段时间很短，忽略月球自转，嫦娥三号的经度基本不变。则嫦娥三号的纬度=44.12N-15.10N=29.02N。所以嫦娥三号近月点的位置为29.02N，19.51W。而远月点的位置与近月点是对称的，所以嫦娥三号远月点的位置为29.02S，159.49E5.3模型二的建立与求解嫦娥三号在快速调整阶段的主要任务是调整姿态，将水平方向的速度减为0，即使主减速发动机的推力竖直向下。由于该阶段主要是减少水平速度，可以

19、认为该阶段的轨迹与嫦娥三号在主减速阶段的轨迹相同。我们建立重力转弯模型2，将发动机的推力方向与嫦娥三号速度矢量的反方向保持一致，实现此阶段的制动减速。因此嫦娥三号在快速调整阶段的轨迹即为式（5-5）。有（5-5）的（1）、（3）可以推导出偏航角控制量的表达式： *=0=180°-arctan(vzf-vzvxf-vx) （5-9）下面推导俯仰角控制量，根据文献，结合（5-5），可以得到垂直方向的加速度，其表达是为：dvydt=6yf-y-2vyf+2vytremtrem2 （5-10）将式（5-10）带入式（5-3），可以得到俯仰角控制量为：*=arcsin6yf-y-2vyf+2v

20、ytrem+fMr2(trem2)aFcos0(trem2) （5-11）式（5-9），（5-11）即为嫦娥三号在姿态调整过程中偏航角和俯仰角控制量的表达式，在该过程中，根据嫦娥三号运动的速度情况即可求得偏航角和俯仰角的度数，然后嫦娥三号通过控制小型姿态调整发动机来调整嫦娥三号的姿态，直到嫦娥三号的水平速度为0，速度方向竖直向下，进入粗避障阶段。5.4模型三的建立与求解5.4.1模型三的建立根据附件三和附件四所提供的图形，本文使用matlab2008a 版本的imread命令读出图形的数据，并对这些数据进行了分析和处理。本文提出分段求点的最大值的思想，即把附件三和附件四读出的高度值，按照合适的

21、精确度进行分段，然后求读出的数据在各个所分段的数目，根据这些数据再做相应处理，调整嫦娥三号的运动。因为两个阶段的方法是一样的，这里本文以第五阶段的判断过程做一次详细介绍，第四阶段不再赘述(第四阶段文件、阶段见附件)。精避障段的要求是根据附件4找出一个最平坦的着陆点，然后选择一个最好的运动路线，从现在的离地面100米的距离到达这个着陆点。很明显，从初始位置到达着陆点，最好的运动方式是直线运动，又因为开始和最终的速度都是零，据此，建立一个先加速后减速的直线运动的运动学模型，至于求解着陆点，本文提出分段求点的最大值的思想，即附件四读出的高度值，按照合适的精确度进行分段，然后求读出的数据在各个所分段的

22、数目，统计这些数据，先确定一个大致的着陆区域，其次，在这个区域中选择一个适合探测器着陆的精确区域大小，具体实现方法可以利用求解这个矩阵区域的方差，方差最小的就是我们所要求的精确区域。5.4.2模型三的求解附件五的数据信息是一个1000×1000的矩阵，本文先求出这个矩阵的最大值和最小值，利用matlab2008a求矩阵最小值为和最大值的函数得出：最小值为4，最大值为224(见附件程序Prog2)为了保证矩阵中所有数据都能够被统计，我们选择了分段起始点为4、终止点为224，并选择以4为步长进行分段，并进行统计（见附件程序Prog3）。根据常识，在分段中数目较少的不能作为着陆点，统计了前

23、5位数目，见表1表 1 位次 1 2 3 4 5 数段 9296 8892 96100 8488 100104 数目 122597 119315 109587 96633 86564 根据表中的数据，分别画出这5个数据段的三维图形（见图2，图3，这里限于篇幅原因只列出两幅图，见附件程序Prog4），通过对比这五个三维图形，可以选择一个平坦的区域，经过比对，选择矩阵行的范围为481600，列范围为：801920，构造一个新矩阵。 经过资料查找，嫦娥三号探测器的大致大小为2m×2m3，为了保证探测器有足够的空间进行着陆，因为照片的分辨率为0.1米，选择从choice2矩阵中选择30

24、15;30的矩阵，分别计算构造的30×30的矩阵的方差，方差最小所对应的那个矩阵，就是要找的着陆点，通过计算得出方差最小数值是 3.6973，根据这个方差值，找到其所对应的矩阵，记为着陆点（见附件程序Prog5），根据着陆点的数据矩阵，画出其图形(见图4)。 探测器在初始位置在1000×1000矩阵的中心位置的正上方100m处，目标位置是着陆点矩阵的中心位置，两点的连线就是嫦娥三号最优的运动轨迹（见图5）。图2图3 图4图555模型四的建立与求解5.5.1垂直软着陆模型的建立 嫦娥三号在100米悬停后，距离月面很近，且下降速度与当地水平面几乎垂直。该段我们采用简单的垂直着陆

25、方式。其中，最简单的是采用分段控制方式，如图6所示，嫦娥三号依次经过悬停、匀加速、匀减速和关机降落几个过程3。 图655.2垂直软着陆模型的求解图二给出的制导率是分段式的，即将推力按照加减速的需要分成若干级别；我们考虑把若干推力级别等效为一段连续的推力过程。这里假设等效推力加速度是线性变化的，其加速和减速过程的等效推力示意如图7所示。 图7根据图二和图三，假设高度h、合加速度a和速度v皆取向上为正，则不妨假设线性推力加速度aF1有如下形式： aF1=Kt+(a0+gm) (5-12)其中，K表示等效推力加速度的斜率，K>0表示推力逐渐增大，K<0表示推力逐渐减小; a0表示等效推力

26、加速度相对于月心重力加速度gm的初始偏移量，a0>0表示初始推力加速度大于gm，a0<0表示初始推力加速度小于gm。对式(5-12)进行两次积分，同时结合关机条件vtfl=voff,htfl=hoff,即可求得式(5-12)的系数K=12hol-hofftfl3+6(v0+voff)/tfl2 (5-13) a0=6hoff-h0ltfl2-2(2v0+voof)/tfl (5-14)式(5-13) (5-14)中下标0和f分别表示着陆段初始和关机时刻状态，下标l表示着陆段，下标off表示关机时刻。六、误差分析和敏感性分析6.1主减速阶段的误差和敏感性分析 我们在主减速阶段建立了燃

27、料最优模型，在主减速结束时嫦娥三号的速度为57m/s，此时嫦娥三号已经基本位于着陆点上方，水平速度很小，我们在模型中假设此时的水平速度为5m/s。然而在实际过程中，嫦娥三号在此时的水平速度很难控制。由于嫦娥三号的合速度为57m/s，因此，水平速度的范围为（0,57）。由式（5-6）(5-7)和式（5-5）中的（4）可以知道，水平速度的改变会影响剩余时间、偏航角和水平位移。下面图8、图9、图10分别给出了剩余时间、偏航角和水平位移关于不同水平速度的图像(画图程序见附件程序Prog14)。如下：剩余时间与水平速度的关系图像图8偏航角与水平速度的关系图像图9水平位移与水平速度的关系图像图10虽然图像

28、上水平速度的范围为（0-57），单实际上此时嫦娥三号速度方向基本向下，水平速度很小。有上面三张图可以看出，水平速度在很小的范围内变化时，在此阶段的剩余时间、偏航角、水平位移偏差很小。尤其是水平位移在（458km-459km）之间。由于在图一中，纬度偏角的两边很大，水平位移在这么小的范围内变化所造成的纬度偏角的误差很小。因此主减速段结束时的水平速度的不确定性对于纬度偏角、偏航角所产生的误差很小。6.2第五阶段的误差分析和敏感度分析精确避障区在进一步确定着陆范围时，所选范围矩阵的大小直接影响系统能否准确、快速的着陆。所选区域越大，数据量越大，找到最平坦的着陆点误差会相对较小，但因为数据处理消耗的时

29、间就越长，探测器做出反应的敏感度就会下降，反之，所选区域越小，数据量越小，找到最平坦的着陆点的误差会相对较大，但对应的数据处理消耗的时间会变小，探测器做出反应的敏感度会上升。据此，分析数据量的大小和数据处理时间以及所选区域的平坦度的关系。其中数据处理时间是利用一台联想双核、主频2.5GHz的i5 笔记本电脑运行相应的matlab程序，得到结果所消耗的时间所体现，运行时间越小，敏感度相对越大；所选区域的平坦度是由此区域多对应的矩阵的方差体现，方差越小，平坦度越大。本文选择8组实验数据，进行分析，分析结果见表（具体实现程序见附件程序Prog6Prog13）。矩阵大小150×150140×140130