铅垂平面飞行弹道仿真及分析报告

.0000 .6554 1.5007 2.3227 3.1197 3.9092;
.0000 .6617 1.5134 2.3409 3.1426 3.9401;
.0000 .6698 1.5304 2.3661 3.1775 3.9835;
.0000 .6792 1.5501 2.3950 3.2162 4.0323;
3）阻力系数表
马赫数
攻角()
0
2
4
6
8
10
0.1
.4177
.4404
.5219
.6603
.8534
1.1023
0.2
.3858
.4086
.4903
.6290
.8226
1.0723
0.3
.3779
.4007
.4827
.6218
.8160
1.0666
0.4
.3785
.4015
.4838
.6234
.8184
0.6
-0.5089
-0.5227
-0.5376
-0.5520
-0.5658
-0.5791
0.7
-0.5366
-0.5502
-0.5649
-0.5790
-0.5927
-0.6058
0.8
-0.5738
-0.5871
-0.6014
-0.6153
-0.6287
-0.6415
0.9
-0.6272
-0.6407
特征面积S(m2)
特征长度L(m)
毛翼展(m)
音速SONIC(m/s)
大气密度 (kg/m3)
0.0227
1.8
0.5
343.13
1.225
使用的插值算法：
气动数据插值——等距双变元抛物线插值；
推力、重心、转动惯量等——不等距一元线性插值。
5 空气动力和空气动力矩表达式
6 编制计算程序
计算机算法采用了matlab实现,源程序见附件.
2 建立数学模型：
3计算方法研究
确定数值积分方法和积分步长
使用算法：四阶龙格—库塔法。
使用步长：h=0.005
4 原始数据：
1）初值
x=0(m)y=20.0(m)=18=18v=20(m/s)z=0(rad/s)m=52.38(kg)
2）攻角与马赫数范围(仅用于插值计算)
攻角=0~10马赫数=0.1~ 0.9
-0.6553
-0.6694
-0.6830
-0.6960
当Xg=.8896时
马赫数
攻角()
0
2
4
6
8
10
0.1
-0.6179
-0.6384
-0.6600
-0.6805
-0.6999
-0.7182
0.2
-0.6207
-0.6410
-0.6626
-0.6830
-0.7024
-0.7207
0.3
-0.6253
-0.7459
-0.7636
0.7
-0.6997
-0.7190
-0.7395
-0.7589
-0.7774
-0.7948
0.8
-0.7435
-0.7624
-0.7824
-0.8014
-0.8194
-0.8365
0.9
-0.8069
-0.8266
-0.8474
-0.8672
-0.8859
-0.9035
11）其它参数
.3787 .4018 .4846 .6249 .8209 1.0738;
.3829 .4062 .4897 .6310 .8248 1.0835;
.3855 .4091 .4934 .6363 .8358 1.0938;
.4082 .4321 .5175 .6621 .8641 1.1254;
.4947 .5192 .6073 .7571 .9672 1.2392];
-0.0053
-0.0229
-0.0391
-0.0538
-0.0693
当导弹重心变化时的修正公式：
10）阻尼力矩导数
当Xg=.9381时
马赫数
攻角()
0
2
4
6
8
10
0.1
-0.4686
-0.4829
-0.4982
-0.5130
-0.5272
-0.5409
0.2
-0.4707
-0.4850
-0.5003
程序中一些主要的变量名有:
L—特征长度;
S—特征面积;
SONIC—声速C;
RHO—大气密度 ;
h—积分步长;
程序中其他变量都是存放中间结果的变量.
7 计算结果
运行程序得到弹道曲线, 速度曲线, 攻角曲线如下:
通过运用龙格库塔法和一些插值方法求解导弹运动方程组,获得了导弹各运动参数的变化规律.通过该算法获得的变化规律比较接近真实情况.
-0.0761
-0.0969
0.6
0.0000
-0.0093
-0.0314
-0.0521
-0.0708
-0.0903
0.7
0.0000
-0.0080
-0.0286
-0.0477
-0.0650
-0.0829
0.8
0.0000
-0.0065
-0.0252
-0.0425
-0.0578
-0.0739
0.9
0.0000
0.8
.0000
.6792
1.5501
2.3950
3.2162
4.0323
0.9
.0000
.6933
1.5935
2.4706
3.3273
4.1790
5）推力数据
t(s)
.000
.15
.49
2.11
2.27
3.53
8.78
25.45
42.80
43.68
44.08
P(kgf)
331.2
614.3
505.4
.5192
.6073
.7571
.9672
1.2392
4）升力系数表
马赫数
攻角()
0
2
4
6
8
10
0.1
.0000
.6430
1.4758
2.2870
3.0713
3.8463
0.2
.0000
.6454
1.4807
2.2942
3.0814
3.8598
0.3
.0000
.6480
1.4858
2.3014
3.0915
附件:源程序
function main()
clear all
clear global
global y;
global ii;
ii=0;
h=0.005;
initl();
while y(7)>=0
ii=ii+1;
result(ii);
rk_4(8,h);
end
savedata(ii);
drawing();
%原始数据初始化
function initl()
global acx acy ajz amzaf amzwz axg ap amc agc andm andaf b L S SONIC RHO;
global y;
y=[0 20. 18 0 18 0 20 52.38];
%马赫数 ma
%攻角 alpha
t(s)
.0
2.0
2.4
ຫໍສະໝຸດ Baidu6.4
10.4
14.4
18.4
22.4
26.4
30.4
34.0
38.4
42.4
44.0
Jz(kgms)
8.35
7.88
7.86
7.81
7.78
7.75
7.73
7.71
7.70
7.70
7.69
7.69
7.69
7.69
8）导弹重心（起自头部）
t(s)
.0
2.0
2.4
10.0
18.0
根据描述飞行器在铅垂平面内运动的数学模型编制某导弹的铅垂平面无控飞行弹道仿真程序利用计算机解算初始段无控飞行弹道对初始段弹道参数的变化规律进行分析
数值实验作业
1 实验背景
实验名称：铅垂平面飞行弹道仿真及分析
实验内容与要求：根据描述飞行器在铅垂平面内运动的数学模型，编制某导弹的铅垂平面无控飞行弹道仿真程序，利用计算机解算初始段无控飞行弹道，对初始段弹道参数的变化规律进行分析。
26.0
32.0
38.0
42.0
44.0
XG(m)
.9381
.9095
.9091
.9026
.8969
.8928
.8907
.8896
.8895
.8896
9）静稳定力矩系数
马赫数
攻角()
0
2
4
6
8
10
0.1
0.0000
-0.0104
-0.0341
-0.0564
-0.0771
-0.0985
0.2
0.0000
该程序有八个函数组成,各函数之间的调用关系如下图所示.
1)子函数initl的功能是输入求解导弹运动方程组所需的原始数据.
2)子函数rk_4是四阶龙格-库塔法积分算法子函数,其中调用了子函数dery.
3)子函数dery功能是计算微分方程组的右端函数,其中调用了子函数interp.
4)插值子函数interp功能是计算所有需要插值的参数,其中调用了子函数interp11和interp33.
5)子函数interp11是不等距单变元线性插值函数,主要用于转动惯量Jz, 质心位置xg, 推力p等单变量的插值.
6)子函数interp33是等距双变元抛物线-线性插值函数,主要用于气动力系数cx, cy, 气动力矩系数 , ,等双变元参数的插值,该函数调用了子函数interp31.
7)子函数interp31是等距单变元抛物线插值函数,被interp33调用,完成子函数interp33的等距双变元抛物线-线性插值功能.
agc:质心位置变化的始末值;
andm:气动数据插值所需的Ma数的最小,最大值;
andaf: 气动数据插值所需的攻角 的最小,最大值;
y:存放积分结果的数组,该数组在程序开始时存放积分初值;
dy:存放右端函数数值的数组;
b:存放三个时间值的数组,其中b(2)存放起飞发动机工作结束时间;b(3)存放续航发动机工作结束时间;
-0.0104
-0.0341
-0.0564
-0.0770
-0.0983
0.3
0.0000
-0.0104
-0.0341
-0.0564
-0.0769
-0.0982
0.4
0.0000
-0.0105
-0.0342
-0.0564
-0.0768
-0.0979
0.5
0.0000
-0.0104
-0.0339
-0.0560
3.8731
0.4
.0000
.6512
1.4923
2.3107
3.1039
3.8891
0.5
.0000
.6554
1.5007
2.3227
3.1197
3.9092
0.6
.0000
.6617
1.5134
2.3409
3.1436
3.9401
0.7
.0000
.6698
1.5304
2.3661
3.1775
3.9835
%升力系数
acy=[.0000 .6430 1.4758 2.2870 3.0713 3.8463;
.0000 .6454 1.4807 2.2942 3.0915 3.8731;
.0000 .6480 1.4858 2.3014 3.0915 3.8731;
.0000 .6512 1.4923 2.3107 3.1039 3.8891;
关于程序中出现的数组和变量名,作如下说明:
acx:阻力系数数组;
acy: 升力系数数组;
amzaf:某一质心位置下的静稳定性导数数组;
amzwz:阻尼力矩系数导数数组;
axg:质心位置Xg随时间的变化规律;
ajz:转动惯量Jz随时间的变化规律;
ap:起飞,续航发动机的推力值;
amc: 起飞,续航发动机的燃料质量的秒流量值;
1.0700
0.5
.3787
.4018
.4846
.6249
.8209
1.0738
0.6
.3829
.4062
.4897
.6310
.8284
1.0835
0.7
.3855
.4091
.4934
.6363
.8358
1.0938
0.8
.4082
.4321
.5175
.6621
.8641
1.1254
0.9
.4947
%三个时间 b(1)为导弹离轨时间 b(2)为起飞发动机工作结束时间 b(3)为续航发动机工作结束时间
b=[0 2.1126 44.0832];
%系数表维数
n1=9;
n2=6;
%andm最小,最大值
andm(1)=0.1;
andm(2)=0.9;
%andaf最小,最大值
andaf(1)=0;
andaf(2)=10;
607.8
48.65
43.97
42.01
41.00
40.80
40.79
2.22
第一级工作结束时间：2.1126s，第二级工作结束时间：44.0832s
6）发动机质量秒流量
t(s)
0.
2.1
2.105
44.1
44.105
100
秒流量(kg/s)
2.362
2.362
0.21059
0.21059
0.
0.
7）转动惯量
-0.6455
-0.6670
-0.6874
-0.7067
-0.7249
0.4
-0.6319
-0.6521
-0.6734
-0.6937
-0.7129
-0.7310
0.5
-0.6424
-0.6624
-0.6835
-0.7036
-0.7226
-0.7406
0.6
-0.6669
-0.6866
-0.7074
-0.7272
.0000 .6933 1.5935 2.4706 3.3273 4.1790];
%推力
ap=[.000 .15 .49 2.11 2.27 3.53 8.78 25.45 42.80 43.68 44.08;
-0.5150
-0.5292
-0.5429
0.3
-0.4744
-0.4886
-0.5039
-0.5186
-0.5327
-0.5464
0.4
-0.4797
-0.4939
-0.5090
-0.5237
-0.5378
-0.5514
0.5
-0.4882
-0.5022
-0.5173
-0.5318
-0.5458
-0.5593
%阻力系数
acx=[.4177 .4404 .5219 .6603 .8534 1.1023;
.3858 .4086 .4903 .6290 .8226 1.0723;
.3779 .4007 .4827 .6218 .8160 1.0666;
.3785 .4015 .4838 .6234 .8184 1.07;