打靶法含Matlab程序Word下载.docx
《打靶法含Matlab程序Word下载.docx》由会员分享，可在线阅读，更多相关《打靶法含Matlab程序Word下载.docx（7页珍藏版）》请在冰豆网上搜索。
这时候你随便给个常数减速力
到方程里飞船八成都能降落，但不是最优解。
想想整个过程，开始降落之前飞船总机械能就那么多，你需要对飞船做负功让机械能减到0。
题目里写发动机喷出翔的相对速度是一定的，直觉告诉
我飞船速度快的时候多喷一些速度慢的时候少喷一些，可以提高
做负功的效率。
但是多喷也不能超过上限7500N所以这就是一个带约束优化问题，matlab里边有专用的优化函数，用fmincon就好。
找出最优解以后把过程画出来，看看F可不可以是那5个状态量的线性组合，如果是的话就非常happy,不是的话再说。
三四阶段你可以扯点图像识别，什么二维复利叶分解找平坦区域，怎么一边下降一边根据自身状态调整路径之类的。
五六阶段还真不知道说什么。
一二阶段肯定是重点啦
（3）误差分析其实还挺难的。
可能的误差来源是地球的引力，月亮绕地球向心加速度，太阳的引力（可能会很小），对自身速度、角度的测量误差（比如你测出自身当前速度100m/s但实际上是105m/s），控制的时候F大小以及角度的误差（比如你想朝正前方向喷2000N但实际上偏了2度而且F=2010N之类）。
上一问已经求出了最优控制策略和飞船路线，把这些扰动加进去以后算出新的路线减掉理想路线求偏差，然后随便用个卡尔曼滤波器把误差给校正
AllforJoy
2014/9/1311:
14:
38
老师的思路，求大神解答给我一份呀
法）。
有限差分法）。
三、实验要求:
运用Matlab/C/C++/Java/Maple/Mathematica等其
中一种语言完成程序设计。
四、实验原理:
1打靶法:
对于线性边值问题
y"
+p（x）y'
+q（x）y=f（X）x可a,b]iy（a）=a,y（b）=P
假设L是一个微分算子使：
Ly=y”+p（x）/+q（x）y
则可得到两个微分方程:
Lyi=f（x）,yi（a）皿,y；
（a^0
-y^+pgy；
+q（x）y1=f（X）,yO"
y^a）=0
（2）
Ly2=°
y2（a）=0,y2（a）二1
二y2+p（x）y2+q（x）y2=0,y2（a）=0,y2（a）=1（3）
方程
（2）,（3）是两个二阶初值问题.假设yi是问题
（2）
的解，y是问题（3）的解，且y2（b）H0，则线性边值问题
（1）
y2（b）
的解为：
y（x）=yi（x）+卩~yi（b）y2（x）。
2.有限差分法:
基本思想是把连续的定解区域用有限个离散点构成
的网格来代替，这些离散点称作网格的节点；
把连续定解区域上的连续变量的函数用在网格上定义的离散变量函数来近似；
把原方程和定解条件中的微商用差商来近似,
分用积分和来近似，于是原微分方程和定解条件就近似地
可以得到原问题在离散点上的近似解。
然后再利用插值方法便可以从离散解得到定解问题在整个区域上的近似解。
五、实验内容:
%线性打靶法
function
[k,X,Y,wucha,P]=xxdb（dydx1,dydx2,a,b,aIpha,beta,h）n=fix（（b-a）/h）;
X=zeros（n+1,1）;
CT仁[alpha,0];
Y=zeros（n+1,length（CT1））;
丫仁zeros（n+1,length（CT1））;
Y2=zeros（n+1,length（CTI））;
X=a:
h:
b;
Y1（1,:
）=CT1;
CT2=[0,1];
Y2（1,:
）=CT2;
fork=1:
n
k1=feval（dydx1,X（k）,Y1（k,:
））x2=X（k）+h/2;
y2=Y1（k,:
）'
+k1*h/2;
k2=feval（dydx1,x2,y2）;
k3=feval（dydx1,x2,Y1（k,:
+k2*h/2）;
k4=feval（dydx1,X（k）+h,Y1（k,:
+k3*h）;
Y1（k+1,:
）=Y1（k,:
）+h*（k1'
+2*k2'
+2*k3'
+k4'
）/6,k=k+1;
end
u=Y1（:
1）
nk1=feval（dydx2,X（k）,Y2（k,:
y2=Y2（k,:
k2=feval（dydx2,x2,y2）;
k3=feval（dydx2,x2,Y2（k,:
k4=feval（dydx2,X（k）+h,Y2（k,:
Y2（k+1,:
）=Y2（k,:
v=Y2（:
1）Y=u+（beta-u（n+1））*v/v（n+1）fork=2:
n+1
wucha（k）=norm（Y（k）-Y（k-1））;
k=k+1;
end
X=X（1:
n+1）;
Y=Y（1:
n+1,:
）;
k=1:
n+1;
wucha=wucha（1:
k,:
P=[k'
X'
Y,wucha'
];
plot（X,Y（:
1）,'
ro'
X,Y1（:
g*'
X,Y2（:
1）,
xlabel（'
轴\itx'
）;
ylabel（'
轴\让y'
）
'
mp'
legend（'
是边值问题的数值解y（x）的曲线'
'
是初值问题1的数值解u（x）的曲线'
是初值问题2的数值解v（x）的曲线'
title（'
用线性打靶法求线性边值问题的数值解的图形
）
%有限差分法
function
[k,A,B1,X,Y,y,wucha,p]=yxcf（q1,q2,q3,a,b,alpha,beta,h）
n=fix（（b-a）/h）;
X=zeros（n+1,1）;
Y=zeros（n+1,1）;
A1=zeros（n,n）;
A2=zeros（n,n）;
A3=zeros（n,n）;
zeros（n,1）;
n
k1（k）=feval（q1,X（k））;
A1（k+1,k）=1+h*k1（k）/2;
k2（k）=feval（q2,X（k））;
A=zeros（n,n）;
B=
A2（k,k）=-2-（hA2）*k2（k）;
A3（k,k+1）=1-h*k1（k）/2;
k3（k）=feval（q3,X（k））;
fork=2:
B（k,1）=（hA2）*k3（k）;
B（1,1）=（hA2）*k3
（1）-（1+h*k1
（1）/2）*alpha;
B（n-1,1）=（hA2）*k3（n-1）-（1+h*k1（n-1）/2）*beta;
A=A1（1:
n-1,1:
n-1）+A2（1:
n-1）+A3（1:
n-1）;
B1=B（1:
n-1,1）;
Y=A\B1;
Y1=Y'
;
y=[alpha;
Y;
beta];
wucha（k）=norm（y（k）-y（k-1））;
y=y（1:
n+1,1）;
k=1:
plot（X,y（:
xlabel（'
）,legend（'
是边值问题的数值解y（x）的曲线'
用有限差分法求线性边值问题的数值解的图形'
）,p=[k'
y,wucha'
打靶法
3.Matlab
源代码：
创建M文件:
ys=dbf（f,a,b,alfa,beta,h,eps）
ff=@（x,y）[y
（2）,f（y
（1）,y
（2）,x）];
xvalue=a:
%）取值范围
n=length（xvalue）
s0=a-0.01;
%选取适当的s的初值
x0=[alfa,s0];
迭代初值
flag=0;
%用于判断精度
y0=rk4（ff,a,x0,h,a,b）;
ifabs（y0（1,n）-beta）<
=eps
flag=1;
y1=y0;
else
s1=s0+1;
x0=[alfa,s1];
y1=rk4（ff,a,x0,h,a,b）;
ifabs（y1（1,n）-beta）<
ifflag~=1
whileabs（y1（1,n）-beta）>
eps
s2=s1-（y1（1,n）-beta）*（s1-s0）/（y1（1,n）-y0（1,n））;
x0=[alfa,s2];
y2=rk4（ff,a,x0,h,a,b）;
s0=s1;
s1=s2;
y0=y1;
y1=y2;
yvalue=y1（1,:
ys=[xvalue'
yvalue'
x=rk4（f,t0,x0,h,a,b）
%rung-kuta
法求每个点的近似值（参考大作业一）
t=a:
%迭代区间m=length（t）;
%区间长度
t
（1）=t0;
x（:
1）=x0;
%迭代初值
fori=1:
m-1
L1=f（t（i）,x（:
i））;
L2=f（t（i）+h/2,x（:
i）'
+（h/2）*L1）;
L3=f（t（i）+h/2,x（:
+（h/2）*L2）;
L4=f（t（i）+h,x（:
+h*L3）;
i+1）=x（:
+（h/6）*（L1+2*L2+2*L3+L4）;
4.举例求二阶非线性方程的边值问题：
在matlab控制台中输入：
f=@（x,y,z）（x八2+z*x八2）;
x0l=0;
x0u=2*exp（-1）;
alfa=0;
beta=2;
h=0.01
dbf（f,x0l,x0u,y0l,y0u,h,1e-6）;
>
y=ans（:
2）;
x=ans（:
1）;
plot（x,y,'
-r'
结果：
再输入：
m=0:
0.01:
2;
n=m.*exp（-1/2*m）;
plot（
升级会员 