写了一个平面弹道仿真软件[MatlabApp及源码]
忆昔长别 2019-11-23原创 航空技术空间技术
中文摘要
之前写了一篇计算远程火箭弹道参数的文章,都是公式看起来不怎么直观,所以我编写了这个软件,还给计算结果加上了动画效果,希望能派上用处。本文内容包含:二体运动下的平面弹道仿真软件,使用说明,APP下载附件,Matlab源码。
关键词
MATLAB弹道仿真轨道

软件功能

卫星平面轨道入轨仿真即根据入轨点的参数,对其后的运动轨迹进行计算。可以模拟出关机点参数与轨道形状的关系。第一、第二宇宙速度等。

火箭自由飞行段弹道仿真根据主动段结束时的关机点参数,推算自由飞行段轨迹。可以观察出射程、飞行时间、再入速度与关机点参数的关系。

近地面的抛体仿真当关机点参数很小时,模型表现为水平地面的斜抛运动。

软件进行数值计算之后,会以动画效果演示飞行轨迹。

软件带有绘图保持功能,可以比较不同的初始参数的仿真结果。

软件运行效果

软件开始界面。

软件开始.PNG

左下角的开始仿真按钮被水印盖住了。

软件帮助.PNG

软件帮助及说明。

弹道导弹自由飞行段的仿真。

软件仿真-弹道.PNG

通过“绘图保持”按钮,可以比较不同关机点参数下的弹道。每次改完参数,然后点击“绘图保持”,再点击“开始仿真”。

软件仿真-多个弹道.PNG

在图形的下方有个文本框,列出本次弹道参数,包括飞行时间s,射程km,弹头再入大气层的速度mk/s。

当弹道与地面没有交点时,得到卫星轨道。

软件仿真-轨道.PNG

这里要注意设置仿真时间,不能太长,不然动画会一直持续很久。

软件仿真-多个轨道.PNG

可以比较不同的速度、弹道倾角对应的轨道情况。速度达到一定程度,弹头逃逸,离开地球。

近地面抛体运动。

软件仿真-抛体.PNG

需要点击“曲线缩放”按钮,放大画面。这样才能看到这么小的抛体轨迹。当然也可以仿真,比较多个参数的效果。

软件仿真-多个抛体.PNG


动画演示。

轨道仿真.gif

弹道仿真.gif

抛体仿真.gif

仿真软件附件:

attachment icon ballistic missile.mlappinstall 42.60KB MLAPPINSTALL 107次下载

运行次软件,需要提前下载Matlab可以安装APP的版本。我用的是Matlab2018b版本。下载之后双击自动安装,然后就可以在我的APP里找到这个应用。下图第二行的第一个。

软件安装.png


软件使用的说明


动力学模型:


目前的模型还比较简单,但是基本的估算是没问题的,尤其是弹头自由飞行段和卫星较高的轨道。以后再考虑加上地球引力摄动模型和大气阻力模型,这样的话计算结果会更精确些,尤其是低轨道和弹道的再入段部分。但是相应的,需要输入的参数也会变多。


数学模型1.PNG

这是一个中心引力场下的二体运动模型,这个是有解析解的,解析解的部分可以参考我以前写的文章。当然也可以通过数值计算求解。

因为平面看图比较直观,而且圆锥曲线的轨道是个平面。于是投影到平面下进行计算仿真。


数值计算方法:


采用变步长的四阶龙格库塔法进行求解。具体的表达式如下

求解方法.PNG

在Matlab里实现很方便,只需要就一句代码 :

 [T_out,Y_out]=ode45(<a href="/u/15908"><a href=/u/15908 target="_blank" data-tag="nkcsource" data-type="at">@daod</a></a>an,[t0 T],Y0,options)



Matlab独立仿真代码源码:


贴到Matlab里就能运行,无需安装。想改的话自己改吧~

%平面弹道导弹被动段弹道仿真(中心引力体运动)
clc;clear;close all;

T=2e4;%最长仿真时间s
higth_k = 10e3;%关机点高度m
speed_k = 7.9e3;%关机点速度大小m/s
theta_k = 50*pi/180;%关机点弹道倾角rad
latitude_k = -50*pi/180;%关机点纬度deg

t0=0;%仿真开始时间s
RE=6371e3;%地球半径m
rk=RE+higth_k;
Y0=[speed_k*sin(theta_k+latitude_k),-speed_k*cos(theta_k+latitude_k),rk*cos(latitude_k),rk*sin(latitude_k)]';
options=odeset('Reltol',1e-4,'AbsTol',1e-4,'Events',@eventfun);
[T_out,Y_out]=ode45(<a href="/u/15908"><a href=/u/15908 target="_blank" data-tag="nkcsource" data-type="at">@daod</a></a>an,[t0 T],Y0,options);

figure(1);
%绘图
h1=plot(Y_out(1,3)*1e-3,Y_out(1,4)*1e-3,'-.ob');
a1=gca;
grid on;
xlabel('x/km');
ylabel('z/km');
%绘制地球
q=0:0.01:2*pi;
x_earth = RE*cos(q);
z_earth = RE*sin(q);
hold on;
h2=plot([0,Y0(3)]*1e-3,[0,Y0(4)]*1e-3,'-ok');
h3=plot(x_earth*1e-3,z_earth*1e-3,'k','LineWidth',2);
axis equal;
for i=1:1:length(T_out)
    set(h1,'xdata',Y_out(1:i,3)*1e-3,'ydata',Y_out(1:i,4)*1e-3);
    pause(0.1);
end

latitude_end = atan2(Y_out(end,4),Y_out(end,3));
distance_fly = abs(latitude_end - latitude_k)*RE*1e-3;
fprintf('弹头飞行距离为:%.3f km\n',distance_fly);

%事件检测函数,数值积分终止条件
function [value,isterminal,direction] = eventfun(T_out,Y_out)
RE=6371e3;%地球半径m
value=sqrt(Y_out(3)*Y_out(3)+Y_out(4)*Y_out(4))-RE+T_out*0; %触发时间,当其值为0的时候,时间会触发 到地球表面的距离为0
isterminal=1;   %设为1时会,触发时间会停止求解器,设0时触发不影响工作
direction=-1;    %触发方向设1时是上升触发,设-1是下降触发,设0是双向触发
end

%二体运动方程
function dy=daodan(t,Y)
dy=zeros(4,1);%Y=[Vx,Vz,x,z]'m
GM=3.98603e14;%地球引力常数
R=sqrt(Y(3)*Y(3)+Y(4)*Y(4))+t*0;%到地心距离m
dy([1,2])=-GM/R^3*[Y(3),Y(4)]';%动力学方程m/s^2
dy([3,4])=Y([1,2]);%运动学方程m/s
end


[修改于 9 个月前 - 2019-11-24 01:12:22]

来自:航天航空 / 空间技术航天航空 / 航空技术
 
7
2019-11-23 5:52:56
忆昔长别(作者)
1楼

GIF图片全挂了。。。。。

[修改于 9 个月前 - 2019-11-23 18:08:52]

评论(1)折叠评论
加载评论中,请稍候...
折叠评论
2020-5-5 8:07:08
2楼

封装成exe好了

折叠评论
加载评论中,请稍候...
折叠评论
2020-05-06 11:16:55
2020-5-6 11:16:55
3楼

建议楼主将该帖子设置一定的可见性,因为可能涉嫌违反联合国2321号决议。。


评论(2)折叠评论
加载评论中,请稍候...
折叠评论

想参与大家的讨论?现在就 登录 或者 注册

%7B%22isDisplay%22%3Atrue%7D

仅供内部学术交流或培训使用,请先保存到本地。本内容不代表科创观点,未经原作者同意,请勿转载。

插入资源
全部
图片
视频
音频
附件
全部
未使用
已使用
正在上传
空空如也~
上传中..{{f.progress}}%
处理中..
上传失败,点击重试
{{f.name}}
空空如也~
(视频){{r.oname}}
{{selectedResourcesId.indexOf(r.rid) + 1}}
处理中..
处理失败
插入表情
我的表情
共享表情
Emoji
上传
注意事项
最大尺寸100px,超过会被压缩。为保证效果,建议上传前自行处理。
建议上传自己DIY的表情,严禁上传侵权内容。
点击重试等待上传{{s.progress}}%处理中...已上传
空空如也~
草稿箱
加载中...
此处只插入正文,如果要使用草稿中的其余内容,请点击继续创作。
{{fromNow(d.toc)}}
{{getDraftInfo(d)}}
标题:{{d.t}}
内容:{{d.c}}
继续创作
删除插入插入
{{forum.displayName}}
{{forum.countThreads}}
篇文章,
{{forum.countPosts}}
条回复
{{forum.description || "暂无简介"}}
ID: {{user.uid}}
学术分隐藏
{{submitted?"":"投诉"}}
请选择违规类型:
{{reason.description}}
支持的图片格式:jpg, jpeg, png