飞行方案大作业1.docx
《飞行方案大作业1.docx》由会员分享，可在线阅读，更多相关《飞行方案大作业1.docx（17页珍藏版）》请在冰点文库上搜索。
飞行方案大作业1
航天飞行动力学大作业
韩谨阳
2015300464
1、方案飞行2、弹道设计3、卫星摄动与机动
第一部分
飞行方案
第三部分
卫星的摄动与机动
第二部分
弹道设计
飞行方案大作业
问题描述
在已知导弹质量、转动惯量、发动机推力等参数的情况下，导弹分为三个飞行方案，即三个阶段飞行。
阶段一：
飞行距离在x<9100m，采用追踪法，其中方案高度与距离的关系、方案弹道倾角与高度的关系如下：
（1）
J=2000cos（0.0003141.1x）5000
、z=k.：
（H-H*）+L（H-H*）
阶段二：
飞行距离在24000mx9100m，采用追踪法，其中方案高度与距离的关系、方案弹道倾角与高度的关系、导弹因燃料消耗而质量改变参数如下：
H*=3050m
*「
（2）
^-k：
（H-H）+kH
ms=0.46kg/s（3）
阶段三：
飞行方案x24000m&&y・0,而最终目标位置为xm=30000m
采用比例导引法
（4）
叵=3
dtdt
vJo=k（q-q。
）
「z二心（==*）+“（」-扌）
要求:
1）计算纵向理想弹道，给出采用瞬时平衡假设mzZ“mzz：
z=0时所有纵向参数随时间的
变化曲线。
2）不考虑气动力下洗影响，计算飞行器沿理想弹道飞行时，你认为可以作为特性点的
以上点处的纵向短周期扰动运动的动力系数，并分析其在特性点处的自由扰动的稳定性,
以及计算在各个特性点处弹体传递函数W.（s）,W^y.（s）,W（s）。
二、建立模型
基于“瞬时平衡”假设，导弹在铅垂平面内运动的质心运动方程组为:
dV
mPcos:
b-Xb-mgsin二dtbb
、,d二
mVPsin:
bYb-mgcos=dt
mjz、z二0
从而可得到导弹攻角与弹道倾角之间的关系
（6）
mTz
mZ
其中假设公式
（1）的'匕=k]-二*）+心（扌-丁）中的心=-9,k寸=-0.5；
2）公式
又因为阶段二需要考虑导弹质量随时间的变化，因此阶段二的模型需要联立公式（
（5）、公式（6）、公式（7）
最后一阶段，因为利用了比例导引法
公式（4）的k=2，可得导弹到达目标的相对微分方程为
r四
dt
而导引率
*
旦★生、其中k=2；
dtdt
因为第三阶段的初始参数及终点坐标均为直角坐标系，由下图可知将
x=30000_rcosq
y=rsinq
代入到公式（4），得到直角坐标系下的微分方程组
另外补充方程法向平衡方程：
dxdr.dq
cosqrsinq-
dtdtdt
dydr..dq
sinq-rsinq-
dtdtdt
de
mVPsin:
:
Y-mgcosv
三、算法实现
编程使用MATLAB软件，并运用欧拉方程解微分方程，即ode45函数;
四、程序源代码
*************************
阶段
******************************
functiondy=jieduan1（t,y）dy=zeros（4,1）;
m=320;
g=9.8;
P=2000;
q=0.5*1.2495*（（288.15-0.0065*y（4））/288.15）.A4.2558*y
（1）.A2;k=-9;
dk=-0.5;
Hi=2000*cos（0.000314*1.1*y（3））+5000;
dHi=-2000*0.000314*1.1*sin（y（3））;
delta=k*（y（4）-Hi）+dk*（dy（3）-dHi）;
alpha=0.34*delta;
Xb=（0.2+0.005*alphaA2）*q*0.45;
Yb=（0.25*alpha+0.05*delta）*q*0.45;
dy=zeros（4,1）;
dy（l）=p*cos（alpha）/m-Xb/m-g*sin（y
（2））;
dy
（2）=P*sin（alpha）/m/y
（1）+Yb/m/y
（1）-g*cos（y
（2））/y
（1）;
dy（3）=y
（1）*cos（y
（2））;
dy（4）=y
（1）*sin（y
（2））;
end
******************************段二******************************
functiondy=jieduan2（t,y）
dy=zeros（4,1）;
m=320-0.46*t;
g=9.8;
P=2000;
q=0.5*1.2495*（（288.15-0.0065*y（4））/288.15）.A4.2558*y
（1）.A2;
k=-0.25;
Hi=3050;
delta=k*（y（4）-Hi）;
alpha=0.34*delta;
Xb=（0.2+0.005*alphaA2）*q*0.45;
Yb=（0.25*alpha+0.05*delta）*q*0.45;
dy
（1）=P*cos（alpha/180*pi）/m-Xb/m-g*sin（y
（2）/180*pi）;
dy
（2）=P*sin（alpha/180*pi）/m/y
（1）+Yb/m/y
（1）-g*cos（y
（2）/180*pi）/y
（1）;
dy（3）=y
（1）*cos（y
（2）/180*pi）;
dy（4）=y
（1）*sin（y
（2）/180*pi）;
end
*******************************[ZAv亠********************************
functiondy=jieduan3（t,y）v=y（4）;
k=10;
m=285.04-0.46*t;
q0=-atan（3050/6000）;
g=9.8;
q1=0.5*1.2495*（（288.15-0.0065*y
（2））/288.15）.人4.2558*丫（4）.人2;
k1=10;
dk1=0.05;
dy=zeros（4,1）;
r=sqrt（y
（1）A2+y
（2）人2）;
q=atan（y
（2）/（y
（1）-30000））;
elta=q-y（3）;
dr=-v*cos（elta）;
tht=q0+k*（q-q0）;
dq=v/r*sin（elta）;
dtht=k*dq;
delta=k1*（y（3）-tht）+dk1*（dy（3）-dtht）;
alpha=0.34*delta;
dy
（1）=-dr*cos（q）+r*sin（q）*dq;
dy
（2）=-dr*sin（q）-r*cos（q）*dq;
Yb=（0.25*alpha+0.05*delta）*q1*0.45;
dy（3）=（2000*sin（alpha）/m+Yb/m-g*cos（y（3）））/v;
y（4）=v;
end
main
***********************************
************************************
m
（1）=287.2204;%导弹质量
P=2000;%发动机推力
g=9.8;
k=5;
det
（1）=0.045;
a
（1）=0.6186;
sit
（1）=-0.000002024;
V
（1）=217.2867;
x
（1）=24000;
H
（1）=3071;
H1
（1）=3050;
S=0.45;
L=2.5;
k仁-0.14;
k2=-0.06;
sit1
（1）=sit
（1）;
p0=1.2495;
T0=288.15;
T
（1）=T0-0.0065*H
（1）;
p
（1）=p0*（T
（1）/T0）A4.25588;
q
（1）=1/2*p
（1）*V
（1）A2;%
Cx
（1）=0.2+0.005*a
（1）A2;
Cy
（1）=0.25*a
（1）+0.05*det
（1）*180/pi;
Y
（1）=Cy
（1）*q
（1）*S;
X
（1）=Cx
（1）*q
（1）*S;
SIT
（1）=（P*sind（a
（1））+（Y
（1）-m
（1）*g*cos（sit
（1））））/m
（1）/V
（1）;
Q
（1）=atan（-H
（1）/（30000-x
（1）））+pi;
r
（1）=6708.2039;
R
（1）=-V
（1）*cos（Q
（1））;
n
（1）=Q
（1）+pi;
SIT1
（1）=k/r
（1）*（V
（1）*sin（n
（1）））;
mza=-0.1;
mzdet=0.024;
t=0;
i=0;
dt=0.01;
ms=0.46;
whileH>0&H1>0
初始速度初始位置
初始高度
参考面积
参考长度
大气密度计算公式
升力系数
俯仰力矩系数对攻角的偏导数俯仰力矩系数对舵偏角的偏导数
质量秒消耗量
运用迭代法求解
i=i+1;
t=t+dt;
det（i+1）=k1*（sit（i）-sit1（i））+k2*（SIT（i）-SIT1（i））;
a（i+1）=-mzdet/mza*det（i）*180/pi;
Cy（i+1）=0.25*a（i）+0.05*det（i）*180/pi;
Cx（i+1）=0.2+0.005*a（iF2;
Y（i+1）=Cy（i）*q（i）*S;
X（i+1）=Cx（i）*q（i）*S;
m（i+1）=m（i）-ms*dt;
sit（i+1）=sit（i）+（P*sind（a（i））+（Y（i）-m（i）*g*cos（sit（i））））/m（i）/V（i）*dt;
V（i+1）=V（i）+（P*cosd（a（i））-（X（i）+m（i）*g*sin（sit（i））））/m（i）*dt;
x（i+1）=x（i）+V（i）*cos（sit（i））*dt;
H（i+1）=H（i）+V（i）*sin（sit（i））*dt;
Q（i+1）=atan（-H（i）/（30000-x（i）））+pi;
sit1（i+1）=k*（Q（i）-Q
（1））;
H1（i+1）=H（i）+V（i）*sin（sit1（i））;
SIT（i+1）=（sit（i+1）-sit（i））/dt;
r（i+1）=（H（i）A2+（30000-x（i））A2）A（1/2）;
R（i+1）=（r（i+1）-r（i））/dt;
n（i+1）=acos（-R（i）/V（i））+pi;
SIT1（i+1）=k/r（i）*（V（i）*sin（n（i）））;
T（i+1）=T0-0.0065*H（i+1）;
p（i+1）=p0*（T（i+1）/T0）A4.25588;
q（i+1）=1/2*p（i+1）*V（i+1）A2;
end
plot（x,H）;
holdon
[t,y]=ode45（'jieduan1',[039.0564],[250007000]）;
plot（y（:
3）,y（:
4））;
holdon
[t,y]=ode45（'jieduan2',[39.0564115],[192.768-0.00991002998.71]）;
plot（y（:
3）,y（:
4））;
其中每一段的初始值，均为上阶段的结束值
所以每一阶段计算结束后，需要再给出所有数据的结果，找到每一段距离相对应的数据，即
为初始值。
五、结果分析
制出导弹三个阶段的飞行轨迹如图
（1）
图
（1）图
（2）是第一阶段纵向参数随时间的变化曲线;
'■曲栈”小上*广二#随曲紡
JW易％1-1
Jtt
■
01
Jhi
&
i.\F.
i
1\/
*4-St*94<■
L|
”1}l|#£kJpJ*4
*十一丁」i_u■宦■
KB
Paul*•■"
MilM-yrUiffl,
r«
rimJ-th■.*_-..a
：
■*b・•Sfrx]fe：
■3
net
■MA*
f9M
1
«K
TtWi
■IM
r^j
fPt
JttC>
in（i
q
T即4喜
U1ipA-1?
n«
图
（2）
图（
3）时第二阶段纵向飞行参数随时间的变化曲线
由图
（1）导弹在第一阶段，从初始高度7000m,开始下降飞行，在距离9100m时,
开始变为登高飞行，距离达到24000m至目标30000m这一阶段为导弹的下降寻找目标阶段；
由图
（2）得，第二阶段的飞行速度先增加后减小，在第一阶段末尾阶段速度减小至192.768m/s；
弹道倾角先减小后增加，海拔高度随时间的增加而减小；
0度和3050m之
由图（3）得，第三阶段为登高飞行，所以弹道倾角和海拔高度分别在间振荡，而速度也基本在140m/s至150m/s之间徘徊；
六、特性点的动力系数、传函
分别取特性点1：
x=0时；
特性点2：
x=9100时；
特性点3:
x=24000时；
特性点4:
x=30000时
由纵向自由扰动的稳定性条件-a22a34.0即纵向自由扰动运动稳定。
根据以下公式：
Ka=
a24
（MZ）o=（mZqSL）o（叭=（mqsL）。
a22
z0
何）0
Z0
a25
z0
（PYv）
（mV）o
（YZ）0
（mV）o
得到以下值:
特性点1
7000
0.6012
250
-2113.912
507.33
特性点2
9100
0.9273
192
-3676.815
882.4357
特性点3
24000
0.92311
206
-1974.445
473.86694
特性点4
30000
1.2495
195
-2672.563
641.4152
特性点1
0.00634921
6.7108
-1.6101
0.05142
0.0052848
-0.0176
0.386
0.0115
特性点2
0.006349
11.672
-2.8014
0.07096
0.009191
-0.02622
0.2927
0.01131
特性点3
0.006349
6.26808
-1.50434
0.04968
0.004936
-0.01668
0.39941
0.011189
特性点4
0.00634921
8.48431
-2.036
0.073581
0.0084
-0.02608
0.3433
0.01372
特性点1的传递函数:
W（s）=-o.017625严86013a1）
Js（0.149s20.0086s1）
wy（s）
0.017625（s0.0063492）
s（0.149s20.0086s1）
0.4496（1一0.0012667（s0.006349）
（0.149s20.00861）
特征点2的传递函数:
w（s）=-
0.02622（0.2927s1）
s（0.0857s20.0066s1）
0.02622（1-0.002689（s0.0063492））
2
s（0.0857s0.0066s1）
0.71035（1-0.0022667（s0.006349）
（0.857s20.0061）
特性点3的传递函数:
W■/s）=-
0.01685（0.39s1）
2
s（0.159s20.0088s1）
wy（s）
0.016（10.00131（s0.0063492）
s（0.159s20.0089s1）
0.33544（10.00131（s0.006349）
（0.1595s20.00891）
特征点4的传递函数:
W（s）=-
0.02608（0.3433a1）
s（0.1178s20.00942s1）
Wy（s）
0.02608（1+0.001419s0.0063492））
s（0.1178s20.0093s1）
0.5189（1+0.00147（s0.006349）
（0.158s20.009421）
升级会员 