简介：openMotor 是一款面向火箭发动机实验人员的开源内部弹道模拟器，基于 Python 开发，适合具备一定推进剂与喷管基础知识的爱好者、学生及科研人员使用。它依据推进剂特性、颗粒几何形状与喷嘴规格估算燃烧室压力和推力，并采用快速行进法确定推进剂退行方式，从而支持任意芯部几何形状。资源包共 156 个文件，约 351KB，以 77 个 py 源码文件为核心，配合 18 个 ric 图标资源、14 个 ui 界面文件、19 个 yaml 配置及少量 csv 数据、rst 文档与构建脚本，整体结构紧凑、便于二次开发。功能上支持公制与英制切换，内置 BATES、Finocyl、Star 等常见晶粒几何，可从 DXF 文件导入自定义几何，并提供推进剂编辑器、晶粒退行可视化编辑器、ENG 文件导出及 Burnsim 导入导出，UI 支持撤销与重做。目前已有 795 人学习下载，适合用于发动机设计验证、教学演示与内部弹道计算研究。

1. 开源内部弹道模拟器到底在算什么：从一次试车数据对不上说起

如果你做过固体火箭发动机的地面试车，大概率遇到过这种场景：推力曲线和仿真结果差了十几个百分点，压强峰值位置对不上，喉部烧蚀比预期严重，回头翻参数表发现燃速系数用的是手册标称值，而实际批次推进剂的燃速指数偏了零点几。openMotor 这类开源内部弹道模拟器要解决的，就是把这个「对不上」的过程从玄学变成可追溯的计算链。它面向的是火箭发动机实验人员——手里有试车台、有推进剂配方、有喷管几何，需要一个能快速迭代、能改源码、能把自己的燃速模型塞进去的工具，而不是一个黑匣子商业软件。内部弹道（internal ballistics）这个词听起来学术，说白了就是：点火之后，推进剂从表面一层层烧掉，燃气经过喷管加速排出，整个过程里燃烧室压强怎么变、推力怎么变、烧了多久。openMotor 把这一套用 Python 写成开源实现，你可以直接读它的零维内弹道求解器，改燃速公式，换喷管类型，跑参数扫描。适合谁？适合那些不满足于「输入参数出曲线」、想搞清楚每一步假设是什么、并且愿意用代码去验证自己试车数据的人。

2. 零维内弹道求解器怎么跑通：从装环境到第一条推力曲线

2.1 为什么选零维模型而不是 CFD

固体火箭发动机内部弹道的建模粒度大致分三档：零维（集总参数，整个燃烧室一个压强值）、一维（沿轴向分布）、三维 CFD。CFD 当然精细，但一次瞬态计算动辄几十小时，而且推进剂燃速模型、侵蚀燃烧、铝颗粒燃烧这些子模型标定起来极其痛苦。对于实验人员来说，真正高频的需求是：改一个喷喉直径，看压强峰值怎么移；换一版燃速系数，看总冲差多少；调一下点火延迟，看初始压强峰能不能压下去。这些用零维模型几分钟就能跑一轮，精度在工程可接受范围内。openMotor 的核心就是零维内弹道：假设燃烧室内压强空间均匀，用质量守恒把推进剂燃气生成率和喷管质量流率平衡起来，时间推进求解。

零维模型的基本方程是燃烧室质量守恒：

d(ρc Vc)/dt = ρp Ab r - mdot_nozzle

其中 ρc 是燃烧室内气体密度，Vc 是自由容积，ρp 是推进剂密度，Ab 是燃烧面积，r 是燃速，mdot_nozzle 是喷管质量流率。燃速通常用圣罗伯特公式 r = a * P^n，a 和 n 是燃速系数和燃速指数，P 是燃烧室压强。喷管质量流率在亚音速/超音速壅塞条件下用等熵关系算。openMotor 把这套方程离散后做时间积分，同时跟踪燃烧面积随燃去距离的变化——这就是为什么它需要你给出推进剂药柱的几何形状。

2.2 装环境与最小可跑配置

openMotor 是 Python 项目，常见做法是用虚拟环境隔离依赖。我一般会建一个干净的 venv，避免和系统里的 numpy 版本打架。

# 创建虚拟环境，Python 3.9 以上 python -m venv openmotor_env source openmotor_env/bin/activate # Windows 用 openmotor_env\Scripts\activate # 克隆仓库（假设你已经拿到源码包） cd openmotor pip install -r requirements.txt # 启动图形界面 python -m openmotor

如果你只想用它的求解器做批量计算，不走 GUI，可以直接 import 核心模块。下面是一个最小示例，构造一个圆柱形药柱、给定燃速参数、跑一条推力曲线：

from openmotor.motor import Motor from openmotor.grain import Grain from openmotor.nozzle import Nozzle # 定义推进剂燃速：圣罗伯特公式 r = a * P^n # a 单位要匹配：P 用 Pa，r 用 m/s propellant = { "a": 5.0e-5, # 燃速系数，典型值范围 1e-5 ~ 1e-4 "n": 0.35, # 燃速指数，典型值 0.2 ~ 0.6 "density": 1700, # 推进剂密度 kg/m^3 "t_min": 0.0, # 点火延迟 } # 圆柱形内孔燃烧药柱 grain = Grain( outer_diameter=0.08, # 外径 80mm inner_diameter=0.02, # 初始内孔 20mm length=0.15, # 药柱长度 150mm propellant=propellant, ) # 喷管：收敛-扩张，给定喉部直径和出口直径 nozzle = Nozzle( throat_diameter=0.012, # 喉部 12mm exit_diameter=0.030, # 出口 30mm efficiency=0.95, # 喷管效率，考虑摩擦和散热损失 ) motor = Motor(grain=grain, nozzle=nozzle, chamber_diameter=0.09) result = motor.solve(dt=1e-4, t_max=3.0) # 输出关键结果 print(f"峰值压强: {result.max_pressure/1e6:.2f} MPa") print(f"总冲: {result.total_impulse:.1f} N·s") print(f"燃烧时间: {result.burn_time:.2f} s")

这段代码的逻辑是：先定义推进剂燃速模型参数，再给药柱几何，再给喷管几何，最后组装成 Motor 对象做时间积分。dt 是时间步长，1e-4 秒对大多数中小型发动机够用；t_max 是最大仿真时间，设得比预期燃烧时间长一点，保证烧完。参数说明里最关键的是 a 和 n——这两个值直接决定压强曲线形状。a 偏大，整体压强水平上移；n 偏大，压强对燃速的反馈更敏感，容易出振荡。如果你从文献拿到的 a 是基于英制单位的，记得换算，这是新手最常见的翻车点。

2.3 药柱几何与燃烧面积演化

零维模型里，燃烧面积 Ab 不是常数，它随燃去距离变化。openMotor 支持多种药柱类型：内孔燃烧、端面燃烧、内外同时燃烧、星形内孔等。每种几何的 Ab 演化规律不同，直接决定压强曲线的形状。以最简单的圆柱内孔燃烧为例，燃去距离 w 增加时，内孔直径变大，燃烧面积 Ab = π * (d_inner + 2w) * L，是线性增加的。这意味着如果喷喉不变，压强会随燃烧进行逐渐下降——因为燃气生成率随面积增加而增加，但喷管壅塞流率只和压强有关，平衡点会移动。

如果你用的是星形内孔，初始燃烧面积大，燃去一小段后面积可能先增后减，压强曲线就会出现一个「平台」甚至「驼峰」。很多试车曲线上的异常峰，回头查就是药柱几何的 Ab 演化没设对。openMotor 的 Grain 类允许你直接给几何参数，也允许你导入自定义的 Ab-w 曲线。我一般会先用内置几何跑一遍，和试车数据对比，如果形状对不上，再考虑自定义曲线。

3. 燃速模型与喷管参数怎么标定：让仿真曲线贴住试车数据

3.1 从试车数据反推燃速系数的实操步骤

试车数据里最直接可用的是压强-时间曲线。如果你有实测压强曲线，可以用它反推燃速系数。思路是：在每一个时间点上，从压强曲线读出 P，从药柱几何算出当时的 Ab，从喷管几何算出 mdot_nozzle，然后质量守恒给出 ρp * Ab * r = mdot_nozzle + d(ρc Vc)/dt。忽略气体积累项（对大多数发动机，这一项在稳态段很小），得到 r ≈ mdot_nozzle / (ρp * Ab)。有了 r 和对应的 P，就可以拟合 r = a * P^n，取对数后线性回归。

import numpy as np # 假设从试车数据读入压强和推力，这里用模拟数据演示 # P_data: 压强 Pa, t_data: 时间 s # 已知药柱几何和喷管喉部面积 At At = np.pi * (0.012/2)**2 # 喉部面积 m^2 rho_p = 1700 # 推进剂密度 Ab = 0.02 # 该时刻燃烧面积，实际应从几何演化算 # 从压强算喷管质量流率（壅塞条件） def mdot_nozzle(P, At, gamma=1.2, R=300, T=3000): # 简化等熵壅塞流率 return At * P * np.sqrt(gamma / (R * T)) * (2/(gamma+1))**((gamma+1)/(2*(gamma-1))) # 对每个数据点算 r P_data = np.array([2e6, 3e6, 4e6, 5e6]) # 示例压强 r_data = np.array([mdot_nozzle(P, At) / (rho_p * Ab) for P in P_data]) # 对数线性回归：ln(r) = ln(a) + n * ln(P) ln_P = np.log(P_data) ln_r = np.log(r_data) n, ln_a = np.polyfit(ln_P, ln_r, 1) a = np.exp(ln_a) print(f"拟合燃速系数 a = {a:.3e}, 燃速指数 n = {n:.3f}")

这段代码的关键假设是忽略气体积累项和侵蚀燃烧。如果你的发动机 L/D 比较大（长径比超过 10），侵蚀燃烧会让实际燃速比静态标定值高，拟合出来的 n 会偏大。这时候要么分段拟合，要么在模型里加侵蚀燃烧修正。参数说明：gamma 是燃气比热比，典型值 1.2 左右；R 是气体常数，T 是燃烧温度，这两个值影响 mdot_nozzle 的绝对值，但对 n 的拟合影响不大，因为取对数后是常数偏移。

3.2 喷管效率与喉部烧蚀的修正

喷管效率 efficiency 这个参数，在 openMotor 里用来修正实际推力与理论推力的偏差。理论推力来自等熵膨胀，实际有摩擦损失、散热损失、非轴向流动损失。对于小型发动机，efficiency 取 0.90~0.95 比较常见；大型发动机可以到 0.96~0.98。如果你发现仿真推力比试车推力系统性偏高，先别急着改燃速，检查 efficiency 是不是设成了 1.0。

喉部烧蚀是另一个容易被忽略的因素。固体发动机工作时，喉部材料在高温燃气冲刷下会逐渐扩大，导致压强下降、推力下降。openMotor 支持给喉部直径一个随时间变化的曲线。如果你有试车后的喉部测量数据，可以把它拟合成线性或幂律扩大，导入模型。没有实测数据的话，保守做法是给一个 2%~5% 的喉部扩大率，看曲线尾部是否更贴合。

3.3 参数扫描：一次跑完 20 组喷喉直径

实验人员经常需要回答「如果喉部大 1mm，压强峰值降多少」。手动改参数跑 20 次太慢，openMotor 的求解器可以嵌在循环里做扫描。

import numpy as np from openmotor.motor import Motor from openmotor.grain import Grain from openmotor.nozzle import Nozzle throat_diameters = np.linspace(0.010, 0.015, 20) # 10mm 到 15mm results = [] for dt_throat in throat_diameters: nozzle = Nozzle(throat_diameter=dt_throat, exit_diameter=0.030, efficiency=0.95) grain = Grain(outer_diameter=0.08, inner_diameter=0.02, length=0.15, propellant={"a": 5.0e-5, "n": 0.35, "density": 1700, "t_min": 0.0}) motor = Motor(grain=grain, nozzle=nozzle, chamber_diameter=0.09) res = motor.solve(dt=1e-4, t_max=3.0) results.append((dt_throat, res.max_pressure, res.total_impulse, res.burn_time)) for d, p, i, t in results: print(f"喉部 {d*1000:.1f}mm: 峰值压强 {p/1e6:.2f} MPa, 总冲 {i:.1f} N·s, 燃烧时间 {t:.2f} s")

这个扫描跑完大概几分钟，能直接给出「喉部直径-峰值压强」的敏感度曲线。注意每次循环都要重新构造 Grain 和 Nozzle 对象，不要复用——openMotor 的求解过程会修改对象内部状态，复用会导致结果错乱。这是我自己踩过的坑，第一次跑扫描时结果完全不对，查了半天发现是对象状态污染。

4. 避坑与排查：试车对不上时先查这五件事

4.1 压强曲线整体偏高或偏低

现象：仿真峰值压强比试车高 15% 以上，曲线形状相似。原因：最常见的是燃速系数 a 偏大，或者推进剂密度给低了。另一个可能是喷管喉部面积算错——有人把直径当半径代入面积公式，差 4 倍。解决：先核对单位，a 的国际单位是 m/s/Pa^n，如果你从文献拿到的是 mm/s 和 MPa，换算关系是 a_SI = a_lit * 1e-3 / (1e6)^n。再核对喉部面积，At = π * (d/2)^2，不是 π * d^2。

4.2 初始压强峰对不上

现象：试车曲线开头有一个尖锐的压强峰，仿真里没有或者太小。原因：点火延迟和初始自由容积没设对。点火药产生的燃气先填充燃烧室自由容积，压强快速上升，然后推进剂才被点燃。openMotor 里 t_min 控制点火延迟，但自由容积还和药柱与壳体之间的间隙有关。解决：把点火延迟调到和试车压强上升沿一致，同时检查 chamber_diameter 和 grain outer_diameter 的差值，这个环形间隙就是初始自由容积的主要来源。

4.3 曲线尾部下降太快

现象：仿真在燃烧后期压强掉得比试车快。原因：喉部烧蚀没考虑，或者药柱尾部有未燃尽推进剂。解决：给喉部直径加一个随时间扩大的曲线，扩大率从试车后喉部测量反推。如果是药柱几何问题，检查 Ab-w 曲线尾部是否下降过陡——实际药柱可能有包覆层或者残药，不会完全按理想几何烧完。

4.4 求解器报错或结果振荡

现象：时间步长调小后结果反而振荡，或者求解器直接报数值错误。原因：零维模型在压强快速变化时刚性较强，显式积分步长太大会发散。解决：dt 从 1e-4 降到 1e-5 试试，或者换用自适应步长。如果还振荡，检查燃速指数 n 是不是设得太大（超过 0.7 时系统容易不稳定），以及喷管质量流率公式在低压下是否出现了负值——低压下壅塞条件不成立，需要切换到亚音速流公式。

4.5 总冲对不上但压强曲线很贴

现象：压强曲线几乎重合，但总冲差 10%。原因：推力系数算错，或者喷管出口直径给错。总冲是推力对时间积分，推力 = 推力系数 * 喉部面积 * 压强。如果压强对但总冲不对，问题一定在推力系数或喉部面积。解决：核对 exit_diameter 和 throat_diameter 的比值，膨胀比不对会直接影响推力系数。另外检查 efficiency 是不是设成了 1.0，实际发动机达不到理想膨胀。

5. 把 openMotor 嵌进你的试车数据处理流程

5.1 用实测压强曲线自动标定燃速参数

前面 3.1 给了手动反推的思路，实际用的时候我会写一个自动标定脚本：读入试车压强曲线，用最小二乘拟合 a 和 n，然后把拟合值写回 openMotor 跑一遍，对比仿真和实测的压强曲线残差。如果残差在稳态段小于 5%，说明标定可用；如果残差集中在某一段，说明那一段有侵蚀燃烧或者几何异常，需要分段处理。

import numpy as np from scipy.optimize import minimize # 假设有实测压强曲线 P_meas(t) 和对应时间 t_meas # 以及一个函数 simulate(a, n) 返回仿真压强曲线 P_sim(t) # 目标：最小化 P_sim 和 P_meas 在重叠时间段的均方误差 def objective(params): a, n = params P_sim = simulate(a, n) # 调用 openMotor 求解器 # 插值到相同时间点后算 MSE P_sim_interp = np.interp(t_meas, t_sim, P_sim) return np.mean((P_sim_interp - P_meas)**2) # 初始猜测 x0 = [5e-5, 0.35] # 对数空间优化，因为 a 跨数量级 res = minimize(lambda x: objective([np.exp(x[0]), x[1]]), [np.log(x0[0]), x0[1]], method='Nelder-Mead') a_opt, n_opt = np.exp(res.x[0]), res.x[1] print(f"标定结果: a = {a_opt:.3e}, n = {n_opt:.3f}")

这个流程的关键是：优化变量取 a 的对数，因为 a 的典型值跨 1~2 个数量级，直接优化容易陷入局部极小。另外，如果试车数据噪声大，先做滑动平均再拟合，否则拟合出的 n 会偏大——噪声让压强和燃速的相关性变强，这是统计上的假象。

5.2 批量对比不同推进剂配方的压强曲线

如果你在做推进剂配方筛选，可以把每个配方的燃速参数作为一组输入，批量跑 openMotor，输出峰值压强、总冲、燃烧时间、最大推力等指标，做成对比表。我一般会固定药柱几何和喷管，只变燃速参数，这样对比的是配方本身的差异。注意不同配方的燃速指数 n 可能差很多，n 大的配方在压强升高时燃速增长快，容易出压强峰，需要配合更大的喷喉来压平曲线。

5.3 一个容易被忽略的验证习惯

每次标定完参数，不要只看拟合残差，还要做一件事：把标定后的参数代回一个「留出」的试车工况——比如同一配方但不同喷喉直径的试车数据——看仿真能不能预测准。如果只在标定工况上贴得好，换个工况就偏，说明你的模型里缺了某个物理效应（通常是侵蚀燃烧或喉部烧蚀）。这个留出验证的习惯，比任何拟合优度都更能说明模型是否可信。我自己吃过亏，一组参数在标定工况上残差 2%，换了个喷喉直接偏 20%，回头查是侵蚀燃烧没考虑。希望帮到你。

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