嫦娥三号软着陆轨道设计与控制策略:建模、优化与仿真

📅 发布时间:2026/10/12 1:46:59
嫦娥三号软着陆轨道设计与控制策略:建模、优化与仿真
简介《2014年全国大学生数学建模竞赛A题嫦娥三号软着陆轨道设计与控制策略》是一份完整的赛题论文面向数学建模参赛者、航天轨道设计爱好者及毕业设计学生。论文以嫦娥三号月球软着陆为背景完整呈现六阶段软着陆过程的建模思路采用逆向推理与微元分析求解近月点与远月点位置、速度及方向利用着陆器质心动力学方程与模拟退火算法优化燃料消耗并结合图像分区与优选法处理粗避障、精避障阶段最终给出安全区域和软着陆轨道。资源包共1个doc文件约771KB内容涵盖摘要、问题重述、模型假设、受力分析与公式推导、结果图表以及误差和敏感性分析结构清晰便于直接阅读和二次修改。目前已有88人学习浏览对于需要完整赛题论文模板或系统理解嫦娥三号轨道控制模型的读者具有较高的参考价值。1. 2014国赛A题一道把轨道设计与控制策略捆在一起的软着陆优化题2014年全国大学数学建模竞赛A题嫦娥三号软着陆轨道设计与控制策略是我见过的竞赛题里少数把轨道设计和控制策略两件事同时推到台前的题目。拿到题表那一晚好多人先去找嫦娥三号的“标准答案”想从公开资料里翻出真实落点坐标直接交差可出题人早已把近月点、远月点的高度和速度、动力下降各阶段的高度范围摊在表里了。真正要交的是一套能从这些边界条件出发、把燃料最省或落点精度最高的优化模型做出来再把推力、姿态控制策略落到每个阶段的方案。无论你在备战国赛、辅导建模队伍还是研究着陆制导算法这条“建模、求解、验证、成文”的链条都能直接拆走复用。2. 先把软着陆问题拆成能算的数学建模坐标系、阶段划分与动力学方程软着陆轨道设计这道题最麻烦的不是公式复杂而是物理图像不统一。环月段是几百公里高的椭圆轨道近月点之后进入动力下降高度从十几公里一路降到零。同一个问题里同时出现轨道力学和短周期控制两套时间尺度第一步必须把坐标系和阶段边界敲定否则后面所有计算都是各说各话。2.1 坐标系怎么选月固系看飞行、月心惯性系看轨道我习惯把问题拆成两个坐标系来用。环月段从远月点到近月点用月心惯性系因为这一段主要受月球中心引力场支配直接写成开普勒轨道就能描述动力下降段从近月点到月面改到月球固连坐标系坐标轴跟着月球自转走着陆点的经纬度固定重力和地形方向都好表达控制律写起来也直观。两套坐标系的转换不复杂转动角速度就是月球自转角速度量级在 2.66e-6 rad/s。动力下降段全程只有十二分钟左右自转带来的切向速度变化不大但如果你把近月点的速度方向反推错后面全盘皆错。我见过不少队伍把近月点速度当成“水平向东”直接写进初值结果轨道反演出来跟题表给出的远月点速度对不上问题就出在坐标转换时漏了牵连速度。重力模型方面软着陆高度从 15 km 到 0相对月面半径 1738 km 来说很小把月球当球体、用中心引力场就够了。要不要加 J2 摄动项从竞赛获奖角度看没必要加了不仅多一组公式还容易论文里冒出“我们考虑了 J2 但影响小于 1%”这种自我拆台的话。我一般会在正文里写明忽略月球的扁率摄动同时保留科氏项再用一句话给数值评估科氏加速度约 10^-4 m/s^2 量级对比主推力产生的加速度小三个数量级以上。后面可以再补充一句。这样既有理论深度又不会让模型重量失控。2.2 六个阶段怎么切近月点、主减速、快速调整到缓速下降2014年A题的题表把动力下降过程分成了主减速段、快速调整段、粗避障段、精避障段、缓速下降段等阶段每个阶段高度范围和控制目标都不一样。这个划分直接决定建模粒度我按常见做法把主减速段的起点定在近月点高度约 15 km、速度约 1.7 km/s这一段是燃料消耗的大头任务是消掉大部分水平速度把高度和速度压到 3 km、60 m/s 这个量级快速调整段把飞行器姿态从接近水平的大倾斜拉回竖直同时继续减速。粗避障段用光学成像在几百米高度上选没有大坑的着陆区精避障段在几十米高度上再挑一次落点缓速下降段从百米高度附近垂直下降到落地前把速度归零。注意每个阶段的切换条件应该用高度阈值而不是时间因为轨道偏差会让同一时刻的高度完全不同我在论文里把阶段切换写成 if h h_switch 这种形式做仿真时直接从上一步的终点状态接下一步的初值程序结构会清爽很多。阶段划分表放在论文里还有个好处就是评委一眼能看出你完整覆盖了题表给的全部阶段。2.3 动力学方程与变质量模型软着陆终端的三个约束状态向量取月固系下的着陆器位置 r、速度 v 和质量 m动力学方程写如下形式r_dot vv_dot F / m * u_vec - mu * r / |r|^3 - 2 Ω × v - Ω × (Ω × r)m_dot -|F| / (Isp * g0)其中 F 是推力大小u_vec 是推力方向单位矢量Isp 是比冲g0 是地球海平面引力加速度。题目表格给出发动机推力和比冲后用第三个式子可以换算质量消耗燃料最优自然等价于终端质量尽量大。这个变质量项很多新手容易漏把质量当常数算算出来的最优轨迹阶段时间误差能到几十秒论文验算时特别难看。终端软着陆条件有三个高度等于月球半径更细一点是等于月球半径加 4 m之后关机靠着陆腿缓冲速度 v 等于零姿态角接近竖直。前两个直接写进约束第三个在二维模型里表现为推力方向与当地垂直方向夹角小于某个阈值。三条一个都不能少缺了水平速度约束优化器一定会“聪明”地给你留一点点水平速度因为这样能省燃料落地却会侧翻。目标函数我建议直接选燃料最优写成 max m(t_f) 或等价地 min ∫F dt。为什么不是时间最短因为工程约束里落点精度和燃料余量优先于着陆时间而且燃料最优在数学上更好写评审时也更容易讲清楚重力损失、推力损失这些物理概念。3. 用数值方法把最优点着轨道算出来直接离散化、初值猜测与跟踪控制模型搭好之后下一关是怎么把最优控制问题解出来。很多队伍在这里卡住因为最优控制听着高深实际是“选什么算法、初值怎么猜、算出来怎么验证”三个工程问题。3.1 为什么国赛里推荐直接法而不是变分法求解最优控制有两条路线。间接法走变分推导最优性条件得到两点边值问题再用打靶法去解共轭变量初值直接法把状态和控制离散成有限个变量把最优控制转成非线性规划问题交给求解器。我自己的经验是国赛三天时间不要碰间接法公式推导一上午一个共轭变量初值猜不好解就直接发散成一条完全不落地的抛物线。直接法则直观得多控制量参数化成若干段常值推力状态用数值积分推进终点约束写成等式约束丢给 SLSQP 或 IPOPT 去迭代。它的好处是每个变量都有明确物理意义出了问题能顺藤摸瓜查是推力越界、质量飞了还是约束里高度量纲写错。下面给一个能抄的最小算例理解了这个结构往二维、三维扩展就是加状态和控制维度的事。3.2 一个能跑起来的最小案例一维垂直软着陆的离散推力优化import numpy as np from scipy.optimize import minimize # 物理参数全用 SI 制避免量纲混乱 mu 4.904e12 # 月球引力常数m^3/s^2 R_m 1738e3 # 月球平均半径m g0 9.80665 # 地球引力加速度m/s^2 Isp 2940 # 比冲s数值以题表为准 P_max 7500.0 # 最大推力N m_fuel 2400.0 # 着陆器初始质量kg h0 15e3 # 初始高度m v0 -50.0 # 初始垂直速度m/s T_total 80.0 # 总飞行时间s N 20 # 控制离散段数 # 优化变量做归一化推力、时间都映射到 [0,1] # 这就是数据规范化处理避免推力 7500 N 和时间 80 s 数量级差太多让梯度失衡 def simulate(x): P x[:N] * P_max # 每段推力大小 dt x[N] * T_total / N # 每段时长 h, v, m h0, v0, m_fuel for k in range(N): a_thrust P[k] / m # 推力加速度 a_grav mu / (R_m h)**2 # 引力加速度 acc a_thrust - a_grav # 半隐式欧拉先更新速度再更新位置稳定性比显式欧拉好 v v acc * dt h h v * dt m m - P[k] / (Isp * g0) * dt if m 100 or h 0: break # 燃料耗尽或提前撞月跳出 return h, v, m def objective(x): h, v, m simulate(x) return -m # 燃料最优 终端质量最大 def constraints(x): h, v, m simulate(x) return [h, v] # 终端高度与速度同时归零 # 初值全程 70% 推力时间取 80s x0 np.hstack([0.7 * np.ones(N), 0.8]) bounds [(0, 1)] * N [(0.5, 1.0)] cons {type: eq, fun: constraints} res minimize(objective, x0, methodSLSQP, boundsbounds, constraintscons, options{ftol: 1e-8, maxiter: 300}) # 看一下最优推力和终端状态 P_opt res.x[:N] * P_max h_opt, v_opt, m_opt simulate(res.x) print(f终端高度 {h_opt:.2f} m终端速度 {v_opt:.2f} m/s剩余质量 {m_opt:.1f} kg)代码里的关键点有三处。第一推力变量用 x[:N] 再乘 P_max时间变量用 x[N] 再乘 T_total把优化变量全部归一化到 0 到 1这就是典型的数据规范化处理。SLSQP 靠有限差分算梯度如果变量物理量纲差异过大比如推力 7500 和时间 80数值梯度会被大量级的变量淹没收敛极慢归一化之后明显稳了。第二半隐式欧拉先更新速度再更新位置比普通欧拉在重力场里不容易出现能量漂移但对这道题来说精度还是不够真正算题时建议换成 scipy.integrate.solve_ivp 的 RK45然后把目标函数和约束改成对仿真终态的取值。第三约束函数返回的是 [h, v] 两个值SLSQP 要求所有等式约束同时为 0这正好对应终端高度和速度归零。需要注意这里的一维版本只是为了让你理解“控制离散化加约束优化”这个套路。真正做 A 题时初始速度要换成题目给出的近月点速度在径向的分量总飞行时间也要放开作为优化变量之一。扩展到二维时状态变成 [x, z, vx, vz, m]控制量变成推力大小和推力方向角两个变量约束再加水平速度 vx 等于 0以及落点 x 落在预选着陆区范围内。框架完全一样只是状态和约束多几行。3.3 把最优轨道变成控制基准跟踪误差与反馈修正离线优化出来的是一组开环最优轨迹工程上叫标称轨道。真实飞行里发动机推力有偏差、导航有噪声照着标称推力曲线开环执行落点误差可能到公里级。国赛论文里这一步常用做法是设计闭环跟踪控制器把标称轨道的状态存成时间表实际状态与标称状态的偏差喂给一个 PD 反馈对推力指令做修正。写成控制律就是 P_cmd P_nom kp * (h_nom - h) kd * (v_nom - v)kp、kd 通过仿真试凑。更漂亮的做法是用 LQR 求最优增益但国赛论文中 PD 加上一组控制策略评估的扰动实验已经能拿到很好的评价。这一段把轨道设计和控制策略连成了线轨道是开环最优解控制是闭环稳定性两者缺一不可。很多队伍只算轨道不谈控制或者只谈 PID 不交代跟踪的是哪条轨道评审一眼就看出两张皮。4. 控制策略怎么落地六个阶段的推力剖面与姿态时序轨道算出来只是上半场下半场要把六个阶段分别给出能执行的推力大小和姿态时序。这是论文里最有工程味道的部分也是新手最容易写得像天书的部分。4.1 主减速段固定推力加姿态角剖面主减速段从近月点开始发动机通常工作在最大推力状态控制量主要是推力方向。把推力方向角参数化成关于时间的三次多项式ψ(t) a0 a1t a2t^2 a3*t^3系数当成优化变量终端满足高度、速度和水平速度条件。为什么用多项式而不是逐点离散因为变量数量少迭代更快论文里画出来的姿态角曲线也平滑缺点是姿态角变化率最大值得检验不然曲线再平滑姿态机也追不上。这一段的物理直觉是初始时水平速度占主导推力方向要大幅朝后倾斜把水平分量消掉随着速度降下来推力方向逐渐转向垂直向上抵抗重力。落地的高度-速度走廊曲线画出来应该是一条从右下向左上平滑逼近零点的弧线。我一般会在论文里突出一个结论主减速段推力大小恒定方向角单调变化这既符合发动机工作特性又方便工程实现。注意这段话里提到的推力方向角多项式是一种常见做法不代表真题标准答案。4.2 快速调整段和避障段推力大小与方向联合调节到了 3 km 以下制动不再是唯一目标位置修正开始上场。快速调整段要把姿态从倾斜拉回接近竖直同时把垂直速度限制在安全范围粗避障段和精避障段则要回答“落到哪里”的问题。粗避障段的常见简化是用光学成像数据生成网格高程计算每个候选点的坡度、粗糙度、高程标准差避开坡度超过阈值或者地形起伏太大的区域在剩余区域里选一个目标。这时候控制量从“消速度”变成“落位置”期望加速度拆成垂向减速项和水平位移项。精避障段的逻辑也一样只是高度更低、约束更紧。我见过不少队伍在这里钻进死胡同试图仿真相机成像、阴影检测、坑沿识别写了几百行代码还是没结果。其实对国赛论文来说用“坡度加粗糙度筛选落点”这个简化完全够用把避障当作落点约束放进优化问题里最后在论文里说明“简化了对地形数据的处理但保留了避障对轨道的影响”。这一点放到第5章避坑部分还会展开。4.3 控制策略参数表写进论文直接抄的时序设计这一节把六个阶段的控制参数整理成表放进论文可以直接用。表里的高度和速度范围是常规量级用来示意相对关系具体数值以你们手上的题表为准不要照抄。阶段高度范围量级主要控制量推力策略终端条件主减速段15 km → 3 km推力方向角最大推力常值高度约 3 km速度约 60 m/s快速调整段3 km → 2 km推力方向、推力大小推力随偏差调节姿态接近竖直速度继续下降粗避障段2 km → 100 m推力大小、水平位移按落点筛选结果机动锁定安全着陆区精避障段100 m → 30 m推力大小、水平位移垂直控制器跟踪水平速度趋近于零缓速下降段30 m → 4 m推力大小PD 跟踪速度剖面速度趋近于零高度 4 m表后的一个重要参数是推力调节范围。变推力发动机不可能瞬间从 7500 N 跳到 1500 N每一次推力指令变化都要检查变化率限制。这个限制在题表里不一定直接给出但对评审来说你的控制策略如果不能在实际发动机上执行整套方案就要打折扣。我在论文里习惯加一句话控制指令按照变推力发动机可实现的斜率做平滑处理并以仿真中的推力变化率作为约束验证。这里也可以嵌入一句关于“控制策略评估”的话——对推力变化率越界做评估能直接看出策略的可执行性。5. 复现这道题时最容易翻车的五个坑与排查思路以下是我反复在赛题复现里见到的五个坑按“现象、原因、解决”写每一条都是血泪经验。5.1 单位不统一数据规范化处理没做现象题目给的是 km 和 km/s引力常数用的是 m 制代进动力学方程发现加速度差得离谱轨道直接飞到几万公里外。原因竞赛题表里的单位很良心常用工程单位但你自己一混合就出乱子。真正的问题不在换算而在换算完之后位置 15000 m 和推力 7500 N 的量级差异会让优化求解器无所适从。解决先建立一张单位换算表贴在论文附件里所有计算统一进 SI 制进入优化器之前再做归一化这就是做数据规范化处理。别嫌这一步麻烦它能省掉你后面一半的 Debug 时间。5.2 近月点速度方向搞反轨道和题表对不上现象用题表给的近月点速度算椭圆轨道总发现远月点的速度大小对不上或者近月点高度对不上。原因近月点位于椭圆轨道的近月点附近速度方向几乎是当地水平方向但题表给的是速度大小方向要自己判断。很多队伍想当然地“垂直向下”算出来当然是错的。解决先用能量守恒和角动量守恒验证初值。在轨道近月点速度方向与位置矢量近似垂直可以根据这个几何关系确定速度方向再用开普勒方程正算一遍远月点速度和题表对照误差在 1% 以内才能继续往下做。这一步验证做完后面就不会白忙。5.3 积分器步长太糙结果发散现象用显式欧拉法步长取 0.1 s动力下降段总长 600 s算到后面位置误差漂到几百米优化器怎么调都不收敛。原因显式欧拉是一阶方法误差随时间线性累积重力场里还会出现能量漂移。时间步长取 0.1 s 看着很细但对十几分钟的飞行来说完全不够。解决换用四阶 Runge-Kutta 或自适应步长求解器MATLAB 用 ode45Python 用 scipy.integrate.solve_ivp相对容差设 1e-6绝对容差按位置和速度的量级分别设置。积分器换掉之后很多“玄学”的发散现象会突然消失。5.4 最优轨道和控制策略写成了两张皮现象论文里先画一条漂亮的软着陆轨道再贴一个 PID 控制框图但问到“控制器跟踪的是哪条轨道”就答不上来。原因轨道设计阶段只做了开环优化控制阶段没有把标称轨道带进去做闭环仿真。两条线各算各的论文自然没有说服力。解决让控制器明确跟踪第 3 章的标称轨道。闭环仿真里加入推力偏差和初始速度扰动输出跟踪误差曲线并画在目标轨道旁边作为对照。这就完成了一次控制策略评估控制器能不能把状态拉回标称轨道落点偏差和燃料余量是多少一张图全说清楚。5.5 避障段不做简化直接卡死现象意识到避障要处理地形数据开始写图像处理代码忙了一晚上仿真的机时全花在假地形上最后的论文里只字未提避障。原因建模时想当然地认为避障必须仿真光学相机忘了竞赛的目的是评估你的建模能力和控制策略不是让你做一套完整的计算机视觉系统。解决把避障转化为落点约束。用“坡度、粗糙度、高程标准差”做筛子选出安全落点区域再把你所在区域的边界作为不等式约束放进优化。论文里写“简化地形感知为区域筛选重点研究避障机动对燃料消耗的影响”既诚实又聚焦。控制策略部分只保留水平位移修正这样论文逻辑反而更完整。6. 从算例到一篇国赛论文图表组织、鲁棒性验证与提分技巧模型和仿真做完了最后 48 小时的重心要转移到“怎么把结果讲成一个完整故事”。国赛评委看论文的速度比你想象中快图表就是你的第一印象。6.1 论文图表怎么组织一图胜千言至少要画五类图软着陆轨道的二维剖面图、高度-速度曲线、推力随时间变化曲线、质量消耗曲线、控制策略的时序图。轨道剖面图用高度作纵轴、水平距离作横轴把各阶段切换点标在曲线上评委能直接看到六个阶段的衔接高度-速度曲线用来验证终端约束推力曲线要能看出主减速段常推力、避障段变推力的特征。表格方面初始状态与终端状态对照表必放列出每个物理量的初值、最优值、终端残差。表格里所有数值统一保留四位有效数字单位写在表头不要在每个格子里重复写。Word 排版记住两个习惯图题放图下方居中表题放表上方居中论文正文里先引用图号再放图。这其实是很多拿奖队伍和拿省奖队伍的直观差距。6.2 验证方法残差、终端误差和扰动实验验证分三个层次。第一层把优化得到的最优控制序列代回原始动力学方程用更高精度的积分器重新推一遍比较重推轨迹和优化轨迹的状态偏差偏差小于 1e-3 量级才算可信。第二层把终端约束的误差做成一张小表高度误差、垂直速度误差、水平速度误差全部列出来让评委自己看数量级。第三层做扰动实验这就是控制策略评估的标准动作对初始高度、初始速度、推力系数各加 ±10% 的扰动重跑闭环仿真记录落点偏差和燃料消耗变化证明策略有鲁棒性。三张表放上去论文的技术含量立刻不一样。6.3 一个让论文提档的进阶技巧参数扫描与灵敏度分析最后一个提分技巧是做灵敏度扫描。对近月点高度 h0、最大推力 P_max、着陆器初始质量 m0 三个参数分别做 0.9 倍到 1.1 倍的变化每种情况重算一次最优燃料消耗画成柱状图。你大概率会发现燃料消耗对初始速度最敏感对推力上限不太敏感这个结论本身就能写出一段“控制策略评估”性质的讨论。我自己的习惯是交卷前留一小时做一张全文假设表把每个简化假设对应的误差来源和影响量级写清楚这个动作的性价比比多调一版参数高得多。做这类赛题真正的分水岭从来不是会不会列公式而是能不能在四天之内的每道弯上及时踩住刹车。这篇笔记把最常翻车的五个地方提前标了出来希望帮到你。本文还有配套的精品资源点击获取