发布时间 : 星期二 文章A14014002 - 陈华国 - 王德望 - 崔锦 - 图文更新完毕开始阅读
而改进后的模型是在原模型微分方程的基础之上，对参数??t?不进行人为的设定，这样使得模型更具使用价值．
原模型的微分方程为：
?dr?dt?v?2?dv?Fsincos?(t)?GM?v?dtmr2r??d??? （1） ?dt??d?Fcos?(t)2v???dt?mrr??dm??F?C?dt假设初始时刻t0?0，终端时刻tf自由．软着陆的初始条件由探测器在椭圆轨道近月点处的状态确定，即
r(0)?r0,v(0)?v0,?(0)?0,?(0)??0,m(0)?m0
为到达月面实现软着陆，显然有如下终端条件：
r(tf)?rf,v(tf)?vf,?(tf)?0
耗燃最优控制设计
将运动方程表示为状态方程的形式
?x?f(x,u)
其中系统状态变量为x??r,v,?,?,m?，控制变量u??F,??．应用Pontryagin极大值
T原理，取终端性能指标为
J??x(tf)?m(0)?m(tf)?m0?m(tf)
构造哈密顿函数为
??H(x,?,u)??Tf(x,u) （2）
其中????r,?v,??,??,?m?，满足
T19
??2?v???Fcos?2??v?2???r??3??v??22rmrr??????v???r?2??r???? ???0 (3)
???v???2?r????2?v?????r??FF???sin???cos?v?22??mmmr??由上可知，终端约束为：
?r(tf)?rf???G??v(tf)?vf?=0
??(t)?f??则横截条件为：
????G???? （4） ?(tf)????x?x?t?tf其中?为拉格朗日乘子．
由（3），（4）方程可知，沿着最优轨迹，有???0.将（1）代入（2），并注意到???0,得到
2v??v?m?m?(sin??cos??)F （5）
r2rmmrC根据Pontryagin极大值原理，最优推力为 H??rv??v??vr?2?????v?m?m?F,sin(?)?cos(?)??0??maxmmrC*F?? （6）
?m?0,?vsin(?)???cos(?)??0?mrC?m此外，由于?不受约束，根据变分法极值条件
?H?0，可得到： ???*??tan?1(综上，最优性能指标控制变量：
?vr) ??u*?F*,?*
??T20
将最优控制律代入状态方程（1）和共轭方程（3），利用初始条件和终端约束对状态方程和共轭方程进行积分，即可以得到软着陆最优轨道．此时，求最优轨道就转化成数学上对两点边值问题的求解．
参考文献
[1] 姜启源等.数学建模(第四版)[M].北京：高等教育出版社.2012.
[2] 胡良剑.孙晓君.MATLAB数学实验[M].北京：清华大学出版社.2006.
[3] 林进福.人造地球卫星的椭圆轨道方程[M].上海工程技术大学学报.第12卷第4期：
53-55，1998.
[4] 单永正.段广仁.应用非线性规划求解月球探测器软着陆最优控制问题
[M].Proceedings of the 26th Chinese Control Conference[C].2007.
21
附 录
附录1：问题二求中求解微分方程组即绘制各个参量关于时间t的关系图的gjfc.m文件的MATLAB
程序：
function y=gjfc(t,x) c=0; %c=pi/5; %c=pi/6; %c=pi/8; %c=pi/4; F=7500; %F=5000; %F=2500; %F=1500; m=2400; C=2940;
a=(F/m)*sin(c); d=(F/m)*cos(c); %u=1.633; u=4.9*10^12;
r=(15+1737.013)*1000;
y=[x(2), a-u/r^2-x(2)^2/x(1), x(2)/x(1), -x(2)*x(4)/r+d/r, 主程序为：
clear all;clc;close all; %ts=0:120; ts=0:190; %ts=220; %ts=0:800;
x0=[1752013,0,0, 9.5319e-004,2400]; [t,x]=ode45('gjfc',ts,x0); [t,x] figure(1)
plot(t,x(:,1),'o-')
xlabel('第一阶段时间t')
ylabel('第一阶段探测器离月球中心的距离r') %figure(2)
%plot(t,-x(:,2),'h') figure(3)
plot(t,x(:,3),'r+')
xlabel('第一阶段时间t')
ylabel('第一阶段探测器方位角th') %figure(4)
%plot(t,x(:,4),'go')
%xlabel('第一阶段时间t')
22
-F/C]';