企业官网建设流程全解析

简介：面向飞行器设计与气动学方向的Matlab仿真资料，聚焦三级火箭发射弹道主动段的三次点火过程，适合航空宇航专业学生、课程设计及弹道仿真初学者使用。压缩包共9个文件，包含8个.m源码文件与1张运行结果图，代码以主程序main.m为核心，拆分为大气密度、声速、质量、推进剂及脉冲点火等独立函数模块，便于逐段理解、调试与重组。包体仅28KB，轻量紧凑，适合快速下载与本地运行。目前已有188人学习，可用于验证三级火箭从点火到入轨的完整主动段弹道。通过该资源可快速理解主动段弹道建模思路、三次点火时序安排、目标轨道参数注入条件，并借助运行结果图对照检查仿真输出；函数模块化设计也便于扩展制导律或推进剂模型，是开展火箭弹道仿真课题的高效起点。

1. 三级火箭主动段弹道仿真：三次点火入轨的Matlab工程到底值不值得拆

每次做运载火箭方案论证，最烦的不是写动力学方程，而是反复调点火时刻。手里这个【气动学】三级火箭发射弹道主动段仿真（三次点火达到目标轨道）的Matlab源码包，就是把这件事打包成一个能直接改参数复现的工程：从起飞到三次点火入轨，气动阻力、重力转弯、级间分离和关机点判断都在同一套代码里闭环。你不需要从零搭模型，只需要改初始质量、比冲和目标轨道根数，就能得到主动段的高度-速度曲线、入轨点参数和三级各自的工作段数据。适合做运载火箭总体方案、弹道设计课程作业或初步入轨能力评估的从业者；对新手，这份源码最大的价值是把“三次点火为什么要分成三段时间”变成看得见的数值结果，而不是论文里的一句结论。

2. 主动段弹道建模：坐标系、力学方程与三级点火逻辑

2.1 从发射坐标系到地心惯性系：程序转弯怎么建模

常见的误解是用地面高度直接代替地心距。主动段仿真必须把发射点的地理坐标先换成地心惯性系位置，否则引力方向、科里奥利项和入轨后的轨道根数全部对不上。拆这个资源包时，第一件事是去查它的坐标转换函数。常见做法是取发射点经度、地心纬度和海拔，先得到发射点在地固系中的位置，再根据格林尼治时角旋转到地心惯性系。

对应的Matlab代码一般是这个结构：

% init_params.m 发射点坐标转换到地心惯性系 lat = deg2rad(28.5); % 发射点地心纬度，单位rad lon = deg2rad(102.0); % 发射点经度，单位rad alt = 0.05e3; % 发射点海拔，单位m Re = 6378.137e3; % 地球赤道半径 f = 1/298.257; % 地球扁率 % 地心纬度和地理纬度不完全相同，这一步做扁率修正 R0 = Re / sqrt(1 - f*(2-f)*sin(lat)^2); x_earth_fixed = (R0 + alt) * cos(lat) * cos(lon); y_earth_fixed = (R0 + alt) * cos(lat) * sin(lon); z_earth_fixed = ((1 - f)^2 * R0 + alt) * sin(lat); % 主动段只有十来分钟，格林尼治时角这里按固定值给 GAST = deg2rad(120.0); r0_ECI = [cos(GAST) sin(GAST) 0; -sin(GAST) cos(GAST) 0; 0 0 1] * [x_earth_fixed; y_earth_fixed; z_earth_fixed];

这里最关键的是扁率修正和自转速度初值。发射台本身随地球自转，初始速度不是零，而是地球自转角速度叉乘位置矢量。不少源码包会漏掉这一项，导致入轨点东西向速度偏差几百米每秒。拿到源代码后，我会先在命令行里手动算一遍发射点的r0_ECI和v0_ECI，再去看主脚本里有没有给初始速度赋值，这一眼就能看出包的质量。

2.2 推力、气动阻力与质量变化：三级火箭的动力学方程

主动段仿真的状态量一般是三组：位置矢量、速度矢量、当前质量。动力学方程把推力、引力、气动力和地球自转相关的哥氏项全部加在加速度上。

% rocket_dyn.m 三级火箭主动段动力学 function dydt = rocket_dyn(t, y, p) r = y(1:3); v = y(4:6); m = y(7); r_norm = norm(r); g = -p.mu * r / r_norm^3; stage = p.current_stage; is_burn = p.current_burn; mdot = 0; thrust = 0; if is_burn mdot = p.mdot_list(stage); thrust = p.Isp(stage) * p.g0 * mdot; end % 推力方向：速度方向叠加程序转弯角 thrust_dir = v / norm(v); acc_thrust = thrust / m * thrust_dir; % 大气密度与气动阻力 altitude = r_norm - p.Re; rho = p.atmo_density(altitude); v_rel = v - cross(p.omega_earth, r); D = 0.5 * rho * norm(v_rel)^2 * p.Cd * p.A_ref(stage); acc_drag = -D / m * v_rel / norm(v_rel); dydt = [v; acc_thrust + g + acc_drag ... - 2*cross(p.omega_earth, v) ... - cross(p.omega_earth, cross(p.omega_earth, r)); -mdot]; end

参数上要看三个地方。第一是比冲Isp，它决定同样推进剂能换多少速度增量，三级火箭通常一级用煤油，三级用氢氧，Isp差别会直接影响状态机里的推力切换。第二是参考面积A_ref，一级最大，三级最小；如果源码包把A_ref写成常数，高空段的阻力误差会小，低空段偏大。第三是pitch_profile，它不是死的角度而是时间的函数，一二级要尽量沿速度方向，三级在入轨前会有小幅俯仰调整。如果动力学的dydt里少了- 2*cross(omega_earth, v)这一项，你会看到高度曲线整体升高，但那不是火箭真的更强，而是坐标系没转对。

2.3 三级火箭的质量结构：干重、推进剂和级间比

三级火箭和单级最大的区别是“丢质量”。级间分离本质上是把已经用完的壳体、发动机和剩余结构重量抛弃掉，让下一级以更小的初始质量继续加速。源码包里通常会有一张参数表，我复现时会把它改成这样：

一级初始质量包含二级和三级以及有效载荷，所以二三级在下行过程中也在持续消耗一级的推进剂。这个表如果某一列比例失调，比如三级比冲特别低，那么三次点火即使全部按时发生也到不了目标轨道，最后你会看到速度曲线在天顶方向缓慢增长，这不是数值问题，是方案问题。

光看推力不够，还要算推重比和重力损失。一级推重比通常在1.3到1.5之间，太低则起飞后掉高度，太高则气动加热和动压过大。二级和三级因为已经进入稀薄大气，推重比可以逐渐降低。三次点火的价值在于：一级只负责把火箭送出稠密大气，二级做加速和高度提升，三级做最终入轨修正。每一次点火都对应一个明确任务段，而不是简单地“多开一次机”。

3. 把三点火跑起来：源码结构、主循环和关机条件

3.1 源码文件划分与解压后怎么找到入口

这类课程源码包解压后的结构一般是：一个主脚本、一个动力学函数、一个事件函数、一个绘图脚本、一个参数初始化脚本。文件清单通常是：

我一般直接在main.m所在目录右键运行。如果源码包把文件放在多个子目录，注意看有没有addpath('./utils')，没有的话要手动加一次。另外，Matlab的函数名不能和脚本名冲突，如果里面既有plot.m又有plot_traj.m，运行时调用plot会被用户脚本覆盖，这个坑经常出现在课程包打包不干净的时候。

% main.m 主入口 clear; close all; clc; p = init_params(); % 读取全部参数 tspan = [0 p.burnout_time(3) + 300]; % 总仿真时间略大于三级关机点 y0 = [p.r0_ECI; p.v0_ECI; p.m0]; % 初始状态 options = odeset('RelTol',1e-9,'AbsTol',1e-9, ... 'Events', @(t,y) rocket_event(t,y,p)); [t, y, te, ye, ie] = ode45(@(t,y) rocket_dyn(t,y,p), tspan, y0, options); % 取关机点状态 r_final = ye(end,1:3); v_final = ye(end,4:6); m_final = ye(end,7); fprintf('关机点地心距 %.3f km, 速度 %.3f km/s\n', ... norm(r_final)/1e3, norm(v_final)/1e3);

这里tspan上界不能拍脑袋填，一般先猜一个总时长，再根据事件返回的关机时间修正。RelTol用1e-9看起来高，但主动段有多个质量突变点，容差不收紧会在级间分离附近返回错误的时间步。如果你只想快速看结果，先跑通main.m，再回来读参数，这是最快的路径。

3.2 三级点火状态机：怎么在Matlab里切换级和关机

三级火箭的推力不是一个连续函数，而是一段一段的。Matlab里最常见的做法是用事件函数和persistent状态变量配合，把整个主动段拆成“一级点火、一级关机、级间分离、二级点火、二级关机、级间分离、三级点火、三级关机”八个事件。事件函数只负责检测零点和终止条件，真正的状态切换由主脚本在每次事件停止后完成。

% rocket_event.m 事件函数：检测关机点与级间分离 function [value, isterminal, direction] = rocket_event(t, y, p) r = y(1:3); v = y(4:6); m = y(7); alt = norm(r) - p.Re; % 一级推进剂耗尽：剩余质量等于一级干重+二级整箭质量+三级整箭质量+载荷 value1 = m - (p.m1_dry + p.m2_init + p.m3_init + p.m_payload); % 二级推进剂耗尽：剩余质量等于二级干重+三级整箭质量+载荷 value2 = m - (p.m2_dry + p.m3_init + p.m_payload); % 三级关机：近地点高度首次达到目标值 value3 = compute_perigee(r, v) - p.target_perigee; value = [value1; value2; value3]; isterminal = [1; 1; 1]; direction = [0; 0; 0]; end

每次ode45停在事件点后，主脚本检查ie返回的事件编号，更新p.current_stage和p.current_burn。一个容易翻车的细节：value1的判据必须减去本级分离后仍然跟着火箭走的全部结构质量，不能只减当前级的推进剂。比如一级关机时，一级壳体和二级、三级、载荷都还在箭上，如果只减一级推进剂，事件触发会偏晚，一级关机时速度偏高，后面的二级和三级就得跟着补偿，最后入轨点差一大截。

3.3 主动段结束判据：半长轴、偏心率和入轨时刻

“三次点火达到目标轨道”不能只看高度到了就算成功。轨道是否圆、近地点是否抬起来，要看半长轴和偏心率。仿真结束那一刻，如果偏心率太大，说明推力方向没沿当地水平面充分加速；如果半长轴小于地球半径加目标近地点高度，说明能量不够。

% orbit_elements.m 由状态向量计算经典轨道根数 function [a, ecc, inc, raan, argp, nu] = orbit_elements(r, v, mu) r_norm = norm(r); v_norm = norm(v); h = cross(r, v); h_norm = norm(h); e_vec = cross(v, h) / mu - r / r_norm; ecc = norm(e_vec); energy = v_norm^2 / 2 - mu / r_norm; a = -mu / (2 * energy); inc = acos(h(3) / h_norm); raan = atan2(h(1), -h(2)); argp = atan2(norm(cross([0;0;1], h)) * e_vec(3), ... dot(cross([0;0;1], h), e_vec)); nu = acos(dot(e_vec, r) / (ecc * r_norm)); end

这个脚本里inc算的是地心惯性系倾角，不是发射方位角。如果目标轨道倾角和发射点纬度相差很大，三次点火顺序再怎么调也救不回来，因为它们只改变面内能量分配，不改轨道面。源码包最后给出的入轨时刻，我一般会在te向量里找倒数第二个事件点，那才是三级真正关机时刻，之后进入无动力滑行段，不能继续算作主动段。

4. 调参与复现：把默认目标轨道改成自己的入轨条件

4.1 初始参数与推力矢量配置

当你想把默认的300km圆轨道改成600km太阳同步轨道，不能只改target_perigee。三级点火时刻、三级比冲和有效载荷质量都会影响最终结果。先列几个我调参时优先动的地方：

% init_params.m 部分参数 p.target_perigee = 600e3; % 目标近地点高度 p.target_apogee = 620e3; % 目标远地点高度 p.m_payload = 4500; % 有效载荷质量 p.Isp(3) = 438; % 三级比冲，单位秒 p.thrust(3) = 110e3; % 三级推力，单位N p.burnout_time = [146.2; 295.8; 1180.3]; % 预置关机时间

注意burnout_time在加了事件函数后只作为初值，真正关机时间由事件决定。我给第三级换推力时，会把三级点火时刻从默认的970秒改成1010秒附近，让二级有更长的滑行段去调整近地点幅角。常见做法是先用一组粗略参数跑通，再看每级分离时的高度和速度，决定是提前还是推后下一次点火。

推力矢量方向也要一起改。三级轨控段的俯仰角不能继续沿用一二级的重力转弯角，否则推力多出来的竖直分量会让远地点抬高但近地点不变，看上去像入轨了，实际上离心率偏大。源码包里如果只给固定pitch_profile，我会改成多项式插值：

% pitch_profile 三级按时间推进方式 pitch = deg2rad(12.0) - 0.0021*(t - 970.0); pitch = max(pitch, deg2rad(1.5));

三个数的含义依次是：初始俯仰角、角速度斜率、最小俯仰角限制。三级点火时火箭已经接近水平，所以俯仰角通常在5°到15°之间缓慢减小。不能设成负数，负数意味着推力指向地面，轨道会立即坠落。如果你发现远地点增量远大于近地点增量，第一反应不是去调点火时间，而是去看俯仰角。

4.2 用二分法确定第三次点火时刻

三级点火时刻不是从参数表里读出来的，而是算出来的。最笨也最可靠的办法是二分法：固定一二级参数，让第三级点火时刻在一定范围内变化，观察关机点的近地点高度是否达到目标值。

% bisection_ignition.m 二分法搜索三级点火时刻 t_low = 960; t_high = 1100; for k = 1:20 t_mid = 0.5 * (t_low + t_high); p.ignition_time(3) = t_mid; [~, ~, ~, ye, ~] = run_main(p); % 封装后的主仿真函数 perigee = compute_perigee(ye(end,1:3), ye(end,4:6)); if perigee > 600e3 t_high = t_mid; else t_low = t_mid; end end

这里判断逻辑用的是“结果高于目标就压低点火时刻”，因为点火时刻直接影响三级入轨时的当地速度方向；具体单调方向要看火箭当前轨道能量，实际项目里我会先用三五个样本点画出近地点随点火时刻的曲线，再决定把t_low还是t_high往哪边挪。20次迭代虽然慢，但在主动段仿真规模下足够稳定。课程级别做到这里已经能说明三次点火的物理意义，工程上则会用序列二次规划或遗传算法去替代。

4.3 多级迭代与轨道面修正：一个简单的参数收敛表

调参过程中我习惯把每一轮的输出记下来，判断趋势是否合理。比如某一轮结果：

可以看出点火时间每推迟一秒，近地点抬升约0.15km，远地点下降约0.2km。如果数据出现近地点和远地点同增，大概率是推力方向里有过多竖直分量，再调时间也压不回目标圆轨道。这个时候要回头修俯仰角剖面，而不是继续二分法。很多仿真包跑出来结果不合理，不是代码错了，是参数之间不满足火箭方程的质量约束。

5. 避坑：初看能跑但结果离谱的五个常见问题

5.1 现象：高度曲线正常，偏心率算出来大于1

这是最典型的“能跑但没意义”。现象是绘图上的高度和速度都在合理范围内，但用轨道根数函数算出的ecc大于1，仿真对象处于双曲线逃逸状态。原因：关机条件用了“速度达到一个固定值”，而这个速度是按圆轨道公式算出来的，却忽略了当时的高度和当地引力，导致速度值偏高。解决方法：关机条件改用能量或半长轴，用a = -mu/(2*epsilon)判断，只有半长轴落在目标范围内才算入轨。

5.2 现象：三次点火只触发了两次，第二次或第三次没发生

现象是te向量里只有一级关机事件，第二级的推进剂消耗值一直停在初始值附近。原因：事件函数里本级推进剂耗尽的判据用了全局变量，但ode45在积分推进过程中会把状态变量复制到多个临时副本，全局变量更新滞后，导致事件值永远不为零。解决：不要用全局变量切换级，把stage和burn以字段形式放进结构体p，在事件触发后由主脚本统一修改p，再继续下一次ode45；或者用persistent，并在每次事件后显式清理。

5.3 现象：低空段速度一直比理论值低，阻力怎么调都没用

现象：初始段速度曲线平缓，超出了气动阻力能解释的范围。原因：动力学里缺少地球自转产生的牵连速度，发射台随地球自转的初速度没有被加入到y0里。解决：初始速度至少要加上cross(omega_earth, r0_ECI)，否则起飞后相对惯性系的速度少了一圈自转分量，在赤道附近可差约465m/s。这个数字对入轨速度7.7km/s来说不小，直接决定一级关机点能不能选对。

5.4 现象：级间分离前后报积分警告，步长反复缩小

现象：ode45在分离时刻附近报“步长必须小于某些点的距离”的警告，甚至中断。原因：级间分离时质量突变和推力突降同时发生，是一个强不连续点，事件函数没有把分离严格切成两段积分。解决：把每个分离点和关机点都定义为terminal事件，在分离点结束积分，更新参数后从该点重新启动下一次ode45；不要试图在同一个tspan内跨过不连续点。这个做法看似多写了几个循环，实际比提高容差靠谱得多。

5.5 现象：近地点和远地点目标都对，但倾角不对

现象：入轨瞬间载荷应该进600km太阳同步轨道，算出的轨道面却偏离目标几度。原因：发射方位角、发射点纬度和目标倾角的关系没满足，或者初始速度的东向分量被方向设置弄丢了。解决：先计算该纬度能直接到达的最小轨道倾角，如果目标倾角小于发射点纬度，必然要做平面机动，主动段仿真无法覆盖；在init_params里将发射方位角设置为beta = asin(cos(inc_target)/cos(lat))，并回头验证y0速度方向确保东向分量不为零。

6. 让仿真结果更可信：主动段关机点的六个轨道根数校验

先说一个我每次拿到类似源码包都会做的小动作：把ode45返回的关机点状态，分别用自编六根数函数和一个独立的解析公式算一遍，如果两者对不上，说明代码里混用了坐标系。六根数校验不需要额外工具箱，我用的是能量守恒和角动量守恒交叉验证。

% 校验1：机械能守恒 specific_energy = norm(v)^2/2 - mu/norm(r); semi_major = -mu/(2*specific_energy); % 校验2：角动量守恒 h = cross(r,v); e_vec = cross(v,h)/mu - r/norm(r); ecc = norm(e_vec); rp = semi_major*(1-ecc); ra = semi_major*(1+ecc);

这里的巧妙之处在于：半长轴和偏心率是两个不依赖坐标系旋转的物理量。只要主动段仿真没有引入额外误差，关机点的rp和ra应当与目标轨道根数一致。如果rp刚好、ra偏大，说明推力方向太陡；如果两者都偏小，说明点火总时长不够或比冲偏低。我还会把事件函数中的“三级关机”和“目标轨道”改成两条件，比如同时要求半长轴误差小于10km和偏心率误差小于0.005，才让仿真停止，这样能避免“高度到了但形状不对”的假阳性。

验证之后，我习惯把三级点火时刻作为横轴，画出近地点高度变化曲线。这条曲线通常有一个近似单调区间，我会在二分法收敛后把最终参数写回init_params.m，然后跑一次完整仿真存档。从那以后，我每次调火箭弹道主动段，都强制先跑一遍六根数校验再去看高度曲线，否则默认拿到的结果只是“能画出来”而不是“能入轨”。希望帮到你。

本文还有配套的精品资源，点击获取