嫦娥三软着陆轨道设计与控制策略.pdf

嫦娥三软着陆轨道设计与控制策略 集团文件发布号：（9816-UATWW-MWUB-WUNN-INNUL-DTY- 嫦娥三号软着陆轨道设计与控制策略 摘要 本文以嫦娥三号登月为背景，研究的是嫦娥三号软着陆轨道设计与最优控制策略问题根据动力学相关原理，建立了嫦娥三号软着陆轨迹模型，得到软着陆过程中各阶段的最优控制策略 针对问题一，通过已知条件求解主减速阶段运动过程，通过水平位移量反推近月点位置建立模型一确定近月点和远月点的位置，以及嫦娥三号速度大小与方向首先以月球中心为坐标原点建立空间坐标系，根据计算的作用力可知地球影响较小，故忽略不计然后将嫦娥三号软着陆看作抛物线的运动过程，计算在最大推力下的减速运动，求得月面偏移距离为 462.4km，由此计算出偏移角度为15.25°从而得出近月点和远月点的经纬度分别为（34.76°W，44.12°N）和（34.76°E，44.12°S）最后在软着陆的椭圆轨道上，由动力势能和重力势能的变化，计算出嫦娥三号在远月点和近月点的速度分别为 1700/和 1615/，沿轨道切线方向 针对问题二，我们根据牛顿第二定律，以每个阶段初始点以及终止点的状态作为约束，以燃料消耗最少作为优化目标，可以建立全局最优模型。

而通过将轨迹离散化，进行逐步迭代从而求得每个阶段的水平位移，并分别得到软着陆过程中六个过程中的着陆轨迹方程以及其对应的最优控制策略而在粗避障以及精避障阶段，我们将所给的数字高程图均分为 9 块，综合考量每一块的相对高程差和平坦度指标来选取最佳着陆点在粗避障阶段，根据燃料消耗最少的目标，选择把先将主减速发动机关闭，在进行一段时间匀加速直线运动后再打开发动机，进行减速直线运动作为最优的控制策略  针对问题三，首先我们改变近月点处到月表的距离和减速发动机的推力这两个因素，对嫦娥三号处的水平位移、燃料消耗等等因素进行灵敏度以及误差的分析，可以观察到近月点离月表的距离与水平位移和燃料消耗均呈线性正相关，同时注意到减速发动机的推力与水平位移呈线性负相关，与该燃料消耗却又呈线性正相关，这也与常识相符合由于嫦娥三号在主减速段水平位移最大，因此我们选取该段从对近月点离月表的距离和减速发动机提供的推力变化两个变量来对模型进行阶段的误差分析，通过计算每个阶段时间的相对误差对最优化后的模型进行误差分析 最后，本文对所建立的模型进行评价，指出优缺点并提出改进的方向 关键词：抛物线；最优控制；线性正相关；相似度 1.04% 一、 问题重述 嫦娥三号于 2013 年 12 月 2 日 1 时 30 分成功发射，12 月 6 日抵达月球轨道。

嫦娥三号在着陆准备轨道上的运行质量为 2.4t，其安装在下部的主减速发动机能够产生 1500N 到 7500N 的可调节推力，其比冲（即单位质量的推进剂产生的推力）为 2940m/s，可以满足调整速度的控制要求在四周安装有姿态调整发动机，在给定主减速发动机的推力方向后，能够自动通过多个发动机的脉冲组合实现各种姿态的调整控制嫦娥三号的预定着陆点为 19.51W，44.12N，海拔为-2641m 嫦娥三号在高速飞行的情况下，要保证准确地在月球预定区域内实现软着陆，关键问题是着陆轨道与控制策略的设计其着陆轨道设计的基本要求：着陆准备轨道为近月点 15km，远月点 100km 的椭圆形轨道；着陆轨道为从近月点至着陆点，其软着陆过程共分为 6 个阶段要求满足每个阶段在关键点所处的状态；尽量减少软着陆过程的燃料消耗根据上述的基本要求，请你们建立数学模型解决下面的问题：  1. 确定着陆准备轨道近月点和远月点的位置，以及嫦娥三号相应速度的大小与方向 2. 确定嫦娥三号的着陆轨道和在 6 个阶段的最优控制策略 3. 对于你们设计的着陆轨道和控制策略做相应的误差分析和敏感性分析 二、 问题分析 2.1 问题一的分析 针对问题一，首先我们需要根据预定的着陆点的经纬度确定轨道，然后通过抛物线的运动计算出在月球着陆时的水平路程，然后计算出偏移角度，据此确定近月点的经纬度，而嫦娥三号的着陆轨道为过月球中心点的椭圆轨道，所以远月点的经纬度和近月点对称，则可以由近月点计算出远月点的经纬度。

最后因为在着陆轨道上卫星的能量守恒，则可以通过势能和动能的转换来计算嫦娥三号的速度和方向 2.2 问题二的分析 针对问题二，确定嫦娥三号着陆轨道以及 6 个阶段的最优控制策略时，始终要满足燃料消耗最小原则在问题 1 中已经对近月点的运动情况进行求解，近月点和主减速段终值点的位置、速度及发动机推力大小均已知，在此基础上，给定准备轨道、主减速段最优控制策略快速调整段主要是对探测器姿态进行调整，采取与主减速同样的建模方法，得到该段质心动力学方程，在满足约束条件及阶段要求下给出具体最优控制策略对于粗避障段，首先对其数字高程图进行划分，对每个区域的平坦程度进行分析，取最平坦区域作为着陆大致范围同样需建立动力学模型对运动轨迹进行描述，考虑燃料消耗最少的目标，选择先匀加速后又在恒定推力下减速至 0 2.3 问题三的分析  针对问题三，我们还是通过改变近月点离月球表面的距离和减速发动机提供的两个变量，通过 MATLAB 进行编程，确定对嫦娥三号在水平位移、燃料消耗等参数方面的灵敏度如何，进而求得近月点离月表的距离与水平位移和燃料消耗之间的相对关系以及更重要的是减速发动机提供的推力与水平位移及燃料消耗之间的相对关系。

最后根据嫦娥三号在主减速段水平位移是最大的，进而我们选取该段从对近月点离月球表面的距离变化和减速发动机提供的推力变化两个角度来进行我们对于模型整体的误差的分析 三、 模型假设 1. 忽略地球等星体引力的影响； 2. 假设嫦娥三号在运行过程中不出现故障； 3. 假设飞行器变化姿态是瞬间完成的； 4. 不考虑任何摩擦力； 5. 假设月球为球体，半径以平均半径为准，并且引力场分布均匀 四、 符号说明 符号 符号说明 Δm 由燃料消耗导致的嫦娥三号减少质量 ℎ1 嫦娥三号主减速阶段初始高度 Fmax 发动机推力上限 M 嫦娥三号净重 ℎ1 主减速阶段终值点速度 m 嫦娥三号的质量 m0 嫦娥三号初始质量 α 发动机推力方向在 xoy 平面的投影与 x 轴的夹角 五、 模型的建立与求解 5.1 问题一 5.1.1 模型的建立与分析  由万有引力公式F =ℎℎℎℎ2计算 F，再由牛顿第二定律F = ma计算地球和月球在近月点和远月点处的重力加速度 表 1 地球及月球在近月点和远月点的重力加速度（单位：m/s2） 近月点 远月点 地球 0.0031 0.0024 月球 1.5966 1.4523 由上表可知，地球在近月点和远月点的重力加速度数值很小，即地球对嫦娥三号与月球影响很小，故可忽略不计。

所以本模型只考虑月球对嫦娥三号的影响 根据附件内的软着陆过程示意图，即嫦娥三号将在近月点 15 公里处以抛物线下降，将其看作匀减速运动过程利用 MATLAB 绘制嫦娥三号绕月飞行的三维动态图，更直观的反应嫦娥三号的环月飞行，如图 1 图 1 同时，由附件内嫦娥三号着陆区域和着陆点示意图可知，只要保证嫦娥三号的着陆区域在虹湾着陆区，即认为着陆成功 为保证嫦娥三号以最大概率降落到精准的着陆点和虹湾着陆区，经过分析后选择以北纬 44.12°作为软着陆的绕月轨道在这种确定纬度的绕月轨道中，月球对嫦娥三号的万有引力，可以分解为两个方向一个是绕月的向心力，一个是与绕行面相切的力，则选择最终状态为绕赤道运行更为准确根据实际分析，嫦娥三号的绕月平面应与南北极轴重合 5.1.2 模型的求解 经查阅资料，嫦娥三号主发动机能够产生从 1500N 至 7500N 的可调节推力 主减速阶段看作平抛运动： 起始速度：ℎℎ=1.7×103ℎ/ℎ  加速度的取值范围：0.625 ≤ℎℎ≤ 3.125m/ℎ2 平抛产生的距离：x =ℎℎ22ℎℎ=462.4ℎℎ 图 2 结合图 2，得到准备着陆的点与软着陆点相差 15.25°，即可算出近月点的经纬度，同时根据对称性，求得远月点的经纬度。

由附件内信息可知距离月球表面 15km 时，速度大小为 1700m/s，将此速度看作近月点速度，从近月点到远月点可看做是重力势能和动能转化的过程，而远月点距离地球表面为 100km，可以计算出重力势能的变化，计算出远月点的速度 ℎℎ2=ℎℎ1−ℎℎ2−ℎℎ1 ℎℎ1=ℎℎ122 ℎℎ2=ℎℎ222 ℎℎ2=ℎℎ112− (ℎℎ月ℎℎ+ℎ−ℎℎ月ℎℎ+ℎ) 由此可得到近月点、远月点的速度以及经纬度： 表 2 近月点、远月点的速度以及经纬度 近月点 远月点 速度 1700m/s 1615m/s 经纬度 (34.76°W,44.12°N) (34.76°E，44.12°S) 5.2 问题二 5.2.1 着陆准备阶段 嫦娥三号在着陆准备轨道上绕月球做椭圆运动，当且仅当其处于近月点时，恰好刚刚进入着陆轨迹基于问题 1 中求解结果可知，在该模型坐标系下，以最小燃料消耗为目标，嫦娥三号在近月点处的飞行状态如下表所示： 表 3 嫦娥三号在近月点处飞行状态数据表 推力方向 速度大小 速度方向 位置坐标 vA 的反向 1692m/s 切于 A 点且α19.51°W31.68° =9.654° N 根据表 3 中相应数据，对嫦娥三号的着陆轨道初始点位置进行选定，其着陆准备轨道必然在近月点与远月点经纬度及着陆点位置三点所构成的平面之内，且由已知条件近月点海拔 15 千米、远月点海拔 100 千米和月球形状扁率 1/963.7256，可以完全确定着陆准备的椭圆轨道，在该轨道近月点处，按照表中飞行状态数据作为最优控制策略，对嫦娥三号进行控制。

5.2.2 主减速阶段最优控制模型 模型的建立 设计主减速段的控制策略时，需根据燃料消耗最小原则进行设计，为此，定义燃料消耗的性能指标[1]： minℎℎ=∫ℎ′(ℎ)ℎℎℎ0（1） 其中，m’为单位时间燃料消耗；Δm 为由燃料消耗导致的嫦娥三号减少质量 由题目可知比冲关系： ℎ推=ℎ′ℎℎ=ℎℎℎℎℎℎ（2） 嫦娥三号的质量变化： m =ℎ0−ℎℎ（3） 其中，m 为嫦娥三号的质量；m0为嫦娥三号初始质量 由于嫦娥三号着陆方式符合重力转弯软着陆[1]的情况，即ℎ推的方向与下降速度方向相反，α为发动机推力方向在 xoy 平面的投影与 x 轴的夹角 将推力分别沿 x，y 方向进行分解，根据牛顿第二定律可得： −ℎ推ℎℎℎℎ=ℎℎℎ（4） mℎ月−ℎ推ℎℎℎℎ=ℎℎℎ（5） 将（2）、（3）带入（4）、（5）中化简可得：  ℎℎ=−ℎ推ℎℎℎℎℎ0−ℎ（6） ℎℎ=ℎ月−ℎ推ℎℎℎℎℎ0−ℎ（7） 利用（6）、（7）并结合质点运动的微分方程，最终得到嫦娥三号软着陆主减速阶段轨迹的微分方程： {ℎ2ℎℎℎ2=−ℎ推ℎℎℎℎℎ0−ℎℎ2ℎℎℎ2=ℎ月−ℎ推ℎℎℎℎℎ0−ℎ（8） 而其中约束条件分为边值约束以及过程约束。

边值约束为： a. 初值约束：ℎ|ℎ=0=ℎ1，ℎℎℎℎ|ℎ=0=ℎℎ，ℎℎℎℎ|ℎ=0=0 b. 终值约束：ℎ|ℎ=ℎ1=ℎ2，ℎℎℎℎ|ℎ=ℎ1=ℎ1 过程约束为： a. 推力约束：0 ≤ℎ推≤ℎℎℎℎ b. 高度约束：h ≥ r c. 质量约束：M ≤ m ≤ℎ0 其中，ℎ1为嫦娥三号主减速阶段初始高度；ℎ2为主减速阶段终值点的高度；ℎ1为主减速阶段终值点速度；t1为主减速阶段时间；Fmax为发动机推力上限；r 为月球半径；M 为嫦娥三号净重 原则为燃料消耗最小原则：m = ∫ℎ推ℎℎℎℎℎ10最小 模型的求解 利用 MATLAB 软件进行求解，求解结果如表 4，着陆轨道如图 3 表 4 主减速阶段最优控制策略 推力 方向 燃料消耗 时间 水平位移 末速度 7500N v 的反向 1075kg 421.2s 377095m 57.15m/s  图 3 主减速阶段运动轨迹 5.2.3 快速调整阶段最优控制模型 模型的建立 采用与模型 1 同样的建模方法，可得到嫦娥三号快速调整阶段动力学方程 {ℎ2ℎℎℎ2=−ℎ推ℎℎℎℎℎ1−ℎℎ2ℎℎℎ2=ℎ月−ℎ推ℎℎℎℎℎ1−ℎ（9） 其中 m1主减速阶段终止点嫦娥三号的质量。

边值约束条件： a) 初值约束：ℎ|ℎ=ℎ1=ℎ2，ℎ|ℎ=ℎ1=ℎ2，ℎℎℎℎ|ℎ=ℎ1=ℎℎ1，ℎℎℎℎ|ℎ=ℎ1=ℎℎ1，ℎ|ℎ=ℎ1=ℎ1 b) 终值约束：ℎℎℎℎ|ℎ=ℎ1+ℎ2=0，ℎ|ℎ=ℎ1+ℎ2=90°，ℎ|ℎ=ℎ1=ℎ3 其中，h3为嫦娥三号在快速调整段终点的高度；t2为快速调整段时间；vx1、vy1分别为嫦娥三号在主减速段终点的速度在 x 轴和 y 轴的分量 整个软着陆过程中过程约束均为模型一中的过程约束条件 燃料最少原则：m = ∫ℎ推ℎℎℎℎℎ1+ℎ2ℎ1最小 5.2.3.2 模型求解： 基于模型 1 相同的求解方法，运用仿真迭代运算，由于初始角度α1确定，通过改变推力推 F 的大小，观察燃料消耗量Δm 变化，利用 MATLAB 软件得到快速调整段推力与各变量之间的关系如下： 图 4 快速调整阶段推力与各变量关系 由题可知，在快速调整段需满足水平末速度为 0,发动机推力方向向下的条件，观察图 2 中推力与水平末速度的关系可以发现，推力至少要大于 5000N 才能够保证水平末速度为 0，而推力大致在 5000N～6000N 之间时，燃油消耗量成明显 快速上升趋势，此推力阶段，下落时间变化也与燃料消耗趋势基本一致，水平位移变化量也近似最大。

基于变化趋势的大致分析能够确定在 4500N～5500N 之间将存在一个最小推力，使得水平末速度恰好为 0，此时燃料消耗量即为最小值，现利用燃料消耗最小原则，利用 MATLAB 软件对推力进行具体数值确定： 表 5 发动机推力与各变量关系 推力 燃料消耗 末速度 水平末速度 时间 水平位移 4501 25.87 23.92 6.96 16.9 307.14 4502 25.88 23.90 6.96 16.9 307.09 … … … … … … 5082 40.79 0.84 0.02 23.6 289.13 5083 42.19 0.54 0.01 24.4 289.06 5084 44.79 0.20 0.00 25.9 288.99 5085 46.53 0.22 0.00 26.9 288.92 … … … … … … 5499 904.90 1.19 0.00 483.8 262.18 5500 925.83 1.17 0.00 494.9 262.07 由表中数据可知，当推力为 5083N 时，水平末速度首次达到 0m/s，此时燃料消耗量为 42.19kg，末状态合速度大小为 0.554/ms，水平位移 289.06m，此阶段运行时间为 24.4s，且在重力转弯软着陆情况下，主减速发动机产生的可调节推力方向始终与合速度方向相反。

根据以上数据，利用 MATLAB 画出快速调整阶段的轨迹图如下 图 5 快速调整阶段轨迹 因此，最优控制策略为 推力 燃料消耗 时间 水平位移 末速度 5083 42.19 24.4 289 0.554 表 6 快速调整阶段最优控制策略 5.2.4 粗避障阶段最优控制模型 5.2.4.1 高程图分析  由于粗避障阶段主要是对大致着陆方位进行初步确定，所以可将原始图片分割成九块小区域，以便于嫦娥三号进一步缩小着陆范围，分割区域及编号如下图所示 图 62400m 高程图分块 为衡量各区域高程差相对于整体区域的分布情况，定义相对高程差作为评价指标： 相对高程=|该区域高程均值-总体高程均值|/总体高程均值 嫦娥三号对着陆区域的筛选目的是为了避开大的陨石坑，而对于陨石坑存在与否的判定，主要是通过所拍照片的各区域高程差的大小来反映，因此，需利用MATLAB 软件来计算每一块区域关于高程差的相关统计量： 表 7 不同区域高程差统计 区域 1 2 3 4 5 6 7 8 9 均值 117.67 117.77 94.13 96.39 106.08 115.33 92.69 108.84 95.87 极差 158 218 60 76 110 173 75 126 76 标准差 30.39 45.49 5.63 10.12 20.12 37.47 8.88 24.59 11.36 相对高程差 0.12 0.12 0.10 0.08 0.01 0.10 0.12 0.04 0.09 针对以上相关统计量的结果进行分析： i. 均值能够反应各个区域平均的高程差情况，但是对于凹凸差异较大的区域而言，不能够作为月面情况的描述指标； ii. 极差能够直观反映出各区域陨石坑凹凸高度变化范围，但是不能对区域整体平滑程度做出判断。

iii. 标准差刻画的是高程波动情况，即该区域高程值的波动； iv. 相对高程刻画的是高程的平均程度  后两个指标都能很好描述该区域是否适合着陆由于这两个指标本身表征的含义不具有统一性，但对比标准统一，都是为了反映各区域高程的整体情况，因此，需把这两个指标进行归一化，采取线性加权的方式，综合成一个平坦度指标且平坦度越小，就越适合作为着陆点 进行归一化处理后得到以下结果： 表 8 归一化处理后不同区域高程差统计 区域 1 2 3 4 5 6 7 8 9 标准差 0.62 1.00 0.00 0.11 0.36 0.80 0.08 0.48 0.14 相对高程差 0.99 1.00 0.83 0.64 0.00 0.79 0.96 0.24 0.68 平坦度 1.61 2.00 0.83 0.75 0.36 1.59 1.04 0.71 0.83 由表中平坦度大小可知，区域 5 的平坦度最小，则表明区域 5 最适合最为着陆点所在范围因此，取区域 5 的中心位置作为粗避障段终值点坐标，由题中所给附件 3 可知，区域 5 中心点位置坐标为（1150,1150），与整个区域中心点重合，因此，粗避障阶段总的水平位移为 0。

5.2.4.2 模型建立 建立嫦娥三号粗避障阶段质心动力学方程： {ℎ2ℎℎℎ2=−ℎ推ℎℎℎℎℎ2−ℎℎ2ℎℎℎ2=ℎ月−ℎ推ℎℎℎℎℎ2−ℎ（10） 其中：m2为快速调整段结束后嫦娥三号的质量 边值约束条件： a) 初值约束：ℎ|ℎ=ℎ1+ℎ2=ℎ3，，，ℎℎℎℎ|ℎ=ℎ1+ℎ2=ℎℎ2， b) 终值约束：ℎℎℎℎ|ℎ=ℎ1+ℎ2+ℎ3=0，ℎ|ℎ=ℎ1+ℎ2+ℎ3=ℎ4，  其中，h4为嫦娥三号在粗避障段终点的高度；t3为粗避障段时间；vy2为嫦娥三号在快速调整段终点的速度在 y 轴的分量 整个软着陆过程中过程约束均为模型一中的过程约束条件 燃料最少原则：m = ∫ℎ推ℎℎℎℎℎ1+ℎ2+ℎ3ℎ1+ℎ2最小 5.2.4.3 模型求解 ① 在此情况下，匀加速直线运动阶段无燃料损失，仅当启动发动机推力时，才会产生燃料的消耗由于匀加速直线运动的时间存在不确定性，因此，需根据相关文献[2]中粗避障段下落时间的理论值，对均加速直线运动时间进行大致设定，进而得到匀加速直线运动的末状态 ② 采用模型 1 相同的求解方法，在假设推力为恒力的条件下，通过计算机迭代仿真求得该情况下燃料的总消耗量及下落总时间。

③ 选取不同的匀加速直线运动时间，进而得到各种情况下的最小燃料消耗量，从而近似确定最小的燃料消耗量及下落时间 利用 MATLAB 软件在选取不同匀加速直线运动时间的情况下，观察燃料消耗量、恒定推力等运动过程描述量的变化情况，得到如下表格： 表 9 运动描述随匀加速时间变化 匀加速时间 总时间 末速度 剩余质量 燃料消耗 恒定推力 25 108.3 0.014 1203.99 74.740 2637.90 26 104.5 0.057 1206.55 72.179 2703.26 30 91.7 0.034 1215.15 63.579 3029.56 35 79.3 0.001 1223.53 55.195 3663.14 40 69.7 0.036 1230.09 48.639 4814.74 45 62.1 0.173 1235.35 43.380 7458.25 在不断改变均加速直线运动时间的情况下，匀加速度直线运动时间取 45s时，燃料消耗量最少为 43.380kg，此时运动总时间最少为 62.1s，推力最大为 7458.25N根据以上变化规律进行推断，若增加匀加速直线运动时间不断，恒定推力的取值也随之增加，而运动总时间不断减少，同时消耗燃料质量也有所下降。

因此最优控制策略为： 表 10 粗避障阶段最优控制策略 时间段 运动状态 推力 方向 0~45 匀加速直线 0 45~62.1 减速直线 7458.25 v 的反向 此策略下，燃料消耗量为 43.380kg 5.2.5 精避障段最优控制模型 5.2.5.1 高程图分析 处理方式同 2400m 处高程图，得出归一化处理后各指标数值如下 表 11 归一化后精避障段各指标 区域 1 2 3 4 5 6 7 8 9 标准差 0.161 0.047 0.183 1 0.589 0.137 0.455 0.262 0 相对高程 0.052 0.250 0.075 0 1 0.385 0.497 0.210 0.415 平坦度 0.214 0.298 0.258 1 1.589 0.521 0.952 0.472 0.415 由表中平坦度大小可知，区域 1 的平坦度最小，则表明区域 1 最适合最为着陆点所在范围因此取区域 1 的中心位置作为精避障段终值点坐标，由题中所给附件 4 中高程图可知，区域 1 中心点位置坐标为（200,800），整个区域中心点坐标为（500,500） 以区域 1 的中心点和整个区域中心点之间的距离作为粗避障阶段总的水平位移，根据两中心点坐标，利用两点间距离公式得到两点在高程图上的距离为 424.3，由题可知，该高程的水平分辨率为 0.1m/像素，则可得实际距离 42.43m。

 5.2.5.2 模型建立 建立嫦娥三号粗避障阶段质心动力学方程： {ℎ2ℎℎℎ2=−ℎ推ℎℎℎℎℎ3−ℎℎ2ℎℎℎ2=ℎ月−ℎ推ℎℎℎℎℎ3−ℎ（11） 其中：m3为粗避障段结束后嫦娥三号的质量 边值约束条件： a) 初值约束：ℎ|ℎ=ℎ1+ℎ2+ℎ3=ℎ4，ℎ|ℎ=ℎ1+ℎ2+ℎ3=ℎ1+ℎ2，ℎℎℎℎ|ℎ=ℎ1+ℎ2+ℎ3=0，ℎℎℎℎ|ℎ=ℎ1+ℎ2+ℎ3=0，ℎ|ℎ=ℎ1+ℎ2+ℎ3=0 b) 终值约束：ℎℎℎℎ|ℎ=ℎ1+ℎ2+ℎ3+ℎ4=0，ℎ|ℎ=ℎ1+ℎ2+ℎ3+ℎ4=ℎ5，ℎ|ℎ=ℎ1+ℎ2+ℎ3+ℎ4=ℎ1+ℎ2+ℎ4 其中，h5为嫦娥三号在精避障段终点的高度；t4为精避障段时间 整个软着陆过程中过程约束均为模型一中的过程约束条件 燃料最少原则：m = ∫ℎ推ℎℎℎℎℎ1+ℎ2+ℎ3+ℎ4ℎ1+ℎ2+ℎ3最小 5.2.5.3 模型求解 ① 采用模型 1 相同的求解方法，在假设推力为恒力的条件下，通过计算机迭代仿真求对该情况下燃料的总消耗量进行求解，但由于精避障阶段运动的距离小，从而可以忽略由于燃料消耗而导致嫦娥三号质量的减少，即假设其在该阶段的质量保持不变，从而对模型进行进一步化简。

② 设定该阶段前一部分的运行时间，改变推力的大小，找到满足精避障段结束后水平速度为 0 的参数，进而确定前一部分运动结束后的状态选取不同的前一部分运行时间，分别求得到各自满足要求的推力大小，根据相关文献[2]中粗避障段下落时间的理论值大致为 22s，从而近似确定推力大小及其运动状态，得到精避障阶段运动轨迹图：  图 7 精避障段轨迹 根据求解结果，最优控制策略为 表 12 精避障阶段最优控制策略 下落时间 推力大小 方向 4.63 2445.37 v 的同向 4.63 2445.37 v 的反向 5.2.6 缓速下降段最优控制模型 5.2.6.1 模型的建立 由于缓速下降段α=90°，因此，可根据粗避障段的化简模型，同理可得嫦娥三号质心动力学方程： ℎ2ℎℎℎ2=ℎ月−ℎ推ℎℎℎℎℎ4−ℎ（12） 其中 m4为精避障段结束后嫦娥三号的质量 边值约束条件： a) 初值约束：ℎ|ℎ=ℎ1+ℎ2+ℎ3+ℎ4=ℎ5，，，ℎℎℎℎ|ℎ=ℎ1+ℎ2+ℎ3+ℎ4=ℎℎ4 b) 终值约束：ℎℎℎℎ|ℎ=ℎ1+ℎ2+ℎ3+ℎ4+ℎ5=0，ℎ|ℎ=ℎ1+ℎ2+ℎ3+ℎ4+ℎ5=ℎ6 其中，h5为嫦娥三号在缓速下降终点的高度；t5为缓速下降段时间。

整个软着陆过程中过程约束均为模型一中的过程约束条件 燃料最少原则：m = ∫ℎ推ℎℎℎℎℎ1+ℎ2+ℎ3+ℎ4+ℎ5ℎ1+ℎ2+ℎ3+ℎ4最小 5.2.6.2 模型求解 同粗避障阶段的化简模型求解方法相同，可得到缓速下降段过程最优控制策略： 表 13 缓速下降阶段最优控制策略 下落时间 推力 方向 3.44 5416.77 竖直向上 5.2.7 整体分析  整个软着陆过程燃料消耗情况如下 表 14 软着陆过程各阶段燃料消耗 阶段 着陆准备 主减速 快速调整 粗避障 精避障 缓速下降 消耗量 0 1075 42.19 43.38 0 0 总消耗量 1160.57 根据燃料消耗情况进行分析，主减速段消耗量最大，对于准备着陆轨道、精屏障段及缓速下降段而言，几乎不产生燃料的损耗 5.3 问题三 5.3.1 软着陆过程分析 综合第二问对于软着陆过程的分析，可有： 利用 MATLAB 软件做出着落过程轨迹图 图 8 通过着陆轨迹图直观反映嫦娥三号着陆的大致曲线 5.3.2 误差分析 1) 时间误差分析 表 15 实际值与模型值的各阶段时间相对误差 阶段 主减速阶段 快速调整阶段 粗避障段 精避障段 缓冲阶段 实际数值 487 16 125 22 19 模型数值 420.9 23.8 62.0 9.28 3.46 相对误差 0.136 0.530 0.510 0.582 0.809 根据上表的相对误差我们可以得到，主减速段的相对误差的最小值为 0.136，同时缓冲阶段的相对误差为 0.809，所以在各个下落阶段时间软着陆过程是一个变量，会受到多方面因素的影响，例如推力、阻力、下落高度的影响。

需要指出的 是，本模型在求解的过程中，考虑的外界因素较少，因此造成了后面的三个阶段下落时间远远小于实际值 2) 水平位移因推力改变产生的相对误差 表 16 主减速阶段推力改变时的水平位移的相对误差 推力 7000 7010 7020 7030 7040 7050 相对误差 0.0511 0.0501 0.0492 0.0488 0.0472 0.0462 同时用趋势图更直观展现二者的关系，即主减速阶段推力对于水平位移相对误差的影响 图 9 根据图像可以得出，主减速阶段的推力在逐渐增大的时候相对误差就会逐渐降低，呈现出负相关而在主减速段最优策略选择时，取推力为 7500N,那么此时的位移的相对误差变为 0故有，前问中所建立的最优控制模型，相对误差较小，取值具有合理性 3) 下落高度改变对于水平位移产生的相对误差 表 17 下落距离 11985 11990 11995 12000 12005 12010 12015 相对误差 0.00002 0.00002 0.00001 0.00000 0.00001 0.00003 0.0003 从表中的数据可以看出，主减速阶段的下落高度对于水平位移的相对的误差影响很小，而且随着下落距离的增加它的变化也有规律性。

所以根据趋势图有： 图 10 首先，在下落高度小于 1.2×104米时，相对误差和下落的高度基本呈现负相关趋势  其次，在下落高度大于 1.2×104米时，相对误差和下落高度基本呈现出正相关趋势 5.3.3 敏感性分析 1） 由于软着陆的过程中下落的高度和嫦娥三号的运动状态之间具有很高的相关性，故我们通过主减速阶段的例子，改变下落高度，观察其他变量的情况即有: 图 11 下落高度同其他变量之间的关系 因此，根据、燃料消耗量等等，与水平方向的夹角基本呈现出线性增大趋势，最终的质量以及末速度和下落的高度呈现出负相关，这与常识相符合，同时这也可以说明下落高度的改变对于其他的各个变量的数值都有较大的影响 2） 其次根据发动机推力时影响嫦娥三号下落高度以及飞行状态的重要因素，故通过减速阶段的例子，通过改变推力 F,观察其他变量的变化情况 图 12 通过趋势图可以比较直观的看出来推力变化与各个变量之间存在着的较强的关联性，推力的变化也就能够影响其余变量的数值了 六、 模型的评价与推广 6.1 模型的优点 1.模型一分别以着陆点的经纬度作为准备着陆轨道，选取精度不变的轨道处于稳定状态，不需要产生推力，保证了燃料消耗的最优； 2.采用动力学模型对嫦娥三号软着陆过程轨迹进行描述，具有严谨的推导过程，且具建立的模型具有很强的推广价值。

6.2 模型的缺点 1.模型建立过程中考虑的因素较少，在处理问题时可能存在一些误差； 2.迭代仿真计算方法求解过程略微繁琐，结果不够精确  6.3 模型的推广 本文主要运用到动力学相关知识以及仿真迭代求解方法，在宏观、微观经济，射击模拟，汽车驾驶训练，军事训练等都有广泛的应用，具有很强的推广意义 七、 参考文献 [1]王鹏基，张熇，曲广吉月球软着陆重力转弯轨道设计与分析[C].//中国宇航学会深空探测技术专业委员会第二届学术会议.2005； [2]张洪华，关轶峰，黄翔宇，嫦娥三号着陆器动力下降的制导导航与控制[J].中国科学:技术科学,2014,4:006； [3]张仲满月球软着陆的制导算法研究[D]，哈尔滨，哈尔滨工业大学，2013 八、 附录 附录一：MATLAB2015a 求解主减速阶段最优控制策略及轨迹 clc;clear;closeall; g=1.633;%òá|óùè m0=2.4*10^3;%àD3êêá Ve=2940;%±è3 theta=9.654*pi/180;%3ùèó·òμD F=7500;%íá| V0=1692.464;%üèμ3ùè t=0;%3êê± T=0.1;%ê±23¤ Vx0=V0*cos(-theta);%3ùè Vy0=V0*sin(-theta);%êú±3ùè  Ay0=g-F*sin(-theta)/(m0-F/Ve*t);%êú±3óùè Ax0=-F*cos(-theta)/(m0-F/Ve*t);%3óùè count=0; X_res=Vx0*t+0.5*Ax0*t^2; Y_res=Vy0*t+0.5*Ax0*t^2; Result=[]; h=12000;%÷ù×μàà %%μü′úó·aùèoí·aò while(Y_res