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

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

软件安装.png


软件使用的说明


动力学模型:


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


数学模型1.png

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

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


数值计算方法:


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

求解方法.png

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

 [T_out,Y_out]=ode45(<a href=/l?t=L3UvMTU5MDg%3D target="_blank">@daod</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=/l?t=L3UvMTU5MDg%3D target="_blank">@daod</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


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

来自:航空航天 / 航天技术航空航天 / 航空技术
3
 
8
已屏蔽 原因:{{ notice.reason }}已屏蔽
{{notice.noticeContent}}
~~空空如也
忆昔长别 作者
4年4个月前 修改于 4年4个月前 IP:黑龙江
866553
被修好了

GIF图片全挂了。。。。。

引用
评论(1)
加载评论中,请稍候...
200字以内,仅用于支线交流,主线讨论请采用回复功能。
折叠评论
星环工贸
3年11个月前 IP:辽宁
880055

封装成exe好了

引用
评论
加载评论中,请稍候...
200字以内,仅用于支线交流,主线讨论请采用回复功能。
折叠评论
yuyue3106
3年11个月前 IP:山东
880122

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


引用
评论(2)
2
加载评论中,请稍候...
200字以内,仅用于支线交流,主线讨论请采用回复功能。
折叠评论

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

所属专业
上级专业
同级专业
忆昔长别
进士 学者 机友 笔友
文章
15
回复
93
学术分
2
2015/12/24注册,5天6时前活动
暂无简介
主体类型:个人
所属领域:无
认证方式:手机号
IP归属地:黑龙江
文件下载
加载中...
{{errorInfo}}
{{downloadWarning}}
你在 {{downloadTime}} 下载过当前文件。
文件名称:{{resource.defaultFile.name}}
下载次数:{{resource.hits}}
上传用户:{{uploader.username}}
所需积分:{{costScores}},{{holdScores}}下载当前附件免费{{description}}
积分不足,去充值
文件已丢失

当前账号的附件下载数量限制如下:
时段 个数
{{f.startingTime}}点 - {{f.endTime}}点 {{f.fileCount}}
视频暂不能访问,请登录试试
仅供内部学术交流或培训使用,请先保存到本地。本内容不代表科创观点,未经原作者同意,请勿转载。
音频暂不能访问,请登录试试
支持的图片格式:jpg, jpeg, png
插入公式
评论控制
加载中...
文号:{{pid}}
投诉或举报
加载中...
{{tip}}
请选择违规类型:
{{reason.type}}

空空如也

加载中...
详情
详情
推送到专栏从专栏移除
设为匿名取消匿名
查看作者
回复
只看作者
加入收藏取消收藏
收藏
取消收藏
折叠回复
置顶取消置顶
评学术分
鼓励
设为精选取消精选
管理提醒
编辑
通过审核
评论控制
退修或删除
历史版本
违规记录
投诉或举报
加入黑名单移除黑名单
查看IP
{{format('YYYY/MM/DD HH:mm:ss', toc)}}