273 lines
8.2 KiB
Python
273 lines
8.2 KiB
Python
"""
|
||
地月转移轨道：有效载荷与燃料消耗（能量）关系分析
|
||
|
||
基于齐奥尔科夫斯基火箭方程和霍曼转移轨道模型
|
||
"""
|
||
|
||
import numpy as np
|
||
import matplotlib
|
||
matplotlib.use('Agg') # 使用非交互式后端
|
||
import matplotlib.pyplot as plt
|
||
from matplotlib import rcParams
|
||
|
||
# 设置中文字体支持
|
||
rcParams['font.sans-serif'] = ['Arial Unicode MS', 'SimHei', 'DejaVu Sans']
|
||
rcParams['axes.unicode_minus'] = False
|
||
|
||
# ============== 物理常数 ==============
|
||
G0 = 9.81 # 地表重力加速度 (m/s²)
|
||
|
||
# 地月转移轨道典型ΔV预算 (m/s)
|
||
DELTA_V_TLI = 3100 # 地月转移注入 (Trans-Lunar Injection)
|
||
DELTA_V_LOI = 800 # 月球轨道捕获 (Lunar Orbit Insertion)
|
||
DELTA_V_TOTAL = DELTA_V_TLI + DELTA_V_LOI # 总ΔV
|
||
|
||
# 不同发动机类型的比冲 (秒)
|
||
ENGINES = {
|
||
'固体火箭': 280,
|
||
'液氧煤油': 350,
|
||
'液氧液氢': 450,
|
||
'离子推进 (仅供参考)': 3000,
|
||
}
|
||
|
||
# 结构系数 (结构质量/燃料质量)
|
||
ALPHA = 0
|
||
|
||
|
||
def exhaust_velocity(isp):
|
||
"""
|
||
计算排气速度
|
||
|
||
参数:
|
||
isp: 比冲 (秒)
|
||
返回:
|
||
排气速度 (m/s)
|
||
"""
|
||
return isp * G0
|
||
|
||
|
||
def mass_ratio(delta_v, isp):
|
||
"""
|
||
计算质量比 R = m0/mf
|
||
|
||
基于齐奥尔科夫斯基方程: ΔV = ve * ln(m0/mf)
|
||
|
||
参数:
|
||
delta_v: 速度增量 (m/s)
|
||
isp: 比冲 (秒)
|
||
返回:
|
||
质量比 R
|
||
"""
|
||
ve = exhaust_velocity(isp)
|
||
return np.exp(delta_v / ve)
|
||
|
||
|
||
def fuel_mass(payload_mass, delta_v, isp, alpha=ALPHA):
|
||
"""
|
||
计算所需燃料质量
|
||
|
||
公式推导:
|
||
m0 = mp + ms + mf (初始质量 = 有效载荷 + 结构 + 燃料)
|
||
mf_final = mp + ms (最终质量)
|
||
ms = alpha * m_fuel (结构质量与燃料成比例)
|
||
|
||
由火箭方程: R = m0 / mf_final
|
||
解得: m_fuel = (mp + ms) * (R - 1)
|
||
m_fuel = mp * (R - 1) / (1 - alpha * (R - 1))
|
||
|
||
参数:
|
||
payload_mass: 有效载荷质量 (kg)
|
||
delta_v: 速度增量 (m/s)
|
||
isp: 比冲 (秒)
|
||
alpha: 结构系数
|
||
返回:
|
||
燃料质量 (kg)
|
||
"""
|
||
R = mass_ratio(delta_v, isp)
|
||
denominator = 1 - alpha * (R - 1)
|
||
|
||
if denominator <= 0:
|
||
# 物理上不可能: 结构质量太大
|
||
return np.inf
|
||
|
||
return payload_mass * (R - 1) / denominator
|
||
|
||
|
||
def total_energy(fuel_mass, isp):
|
||
"""
|
||
计算燃料的动能（近似能量消耗）
|
||
|
||
能量 E = 0.5 * m_fuel * ve²
|
||
|
||
参数:
|
||
fuel_mass: 燃料质量 (kg)
|
||
isp: 比冲 (秒)
|
||
返回:
|
||
能量 (焦耳)
|
||
"""
|
||
ve = exhaust_velocity(isp)
|
||
return 0.5 * fuel_mass * ve ** 2
|
||
|
||
|
||
def characteristic_energy(delta_v):
|
||
"""
|
||
计算特征能量 C3 (用于轨道能量分析)
|
||
|
||
C3 = V_infinity² (相对于地球的双曲线超速平方)
|
||
对于地月转移，C3 ≈ -2 km²/s² (束缚轨道)
|
||
|
||
简化计算: C3 ≈ (ΔV)² - V_escape²
|
||
"""
|
||
# 近地轨道逃逸速度约 10.9 km/s，轨道速度约 7.8 km/s
|
||
# 从LEO出发需要的额外速度约 3.1 km/s
|
||
return (delta_v / 1000) ** 2 # km²/s²
|
||
|
||
|
||
def plot_payload_fuel_relationship():
|
||
"""
|
||
绘制有效载荷与燃料消耗的关系图
|
||
"""
|
||
# 有效载荷范围: 100 kg 到 50000 kg
|
||
payloads = np.linspace(100, 50000, 500)
|
||
|
||
fig, axes = plt.subplots(2, 2, figsize=(14, 12))
|
||
|
||
# ========== 图1: 不同发动机的燃料消耗 ==========
|
||
ax1 = axes[0, 0]
|
||
for engine_name, isp in ENGINES.items():
|
||
fuels = [fuel_mass(p, DELTA_V_TOTAL, isp) for p in payloads]
|
||
fuels = np.array(fuels)
|
||
# 过滤掉无穷大值
|
||
valid = fuels < 1e9
|
||
ax1.plot(payloads[valid]/1000, fuels[valid]/1000,
|
||
label=f'{engine_name} (Isp={isp}s)', linewidth=2)
|
||
|
||
ax1.set_xlabel('有效载荷质量 (吨)', fontsize=12)
|
||
ax1.set_ylabel('燃料质量 (吨)', fontsize=12)
|
||
ax1.set_title('地月转移：有效载荷 vs 燃料消耗\n(ΔV = 3.9 km/s)', fontsize=14)
|
||
ax1.legend(loc='upper left')
|
||
ax1.grid(True, alpha=0.3)
|
||
ax1.set_xlim(0, 50)
|
||
|
||
# ========== 图2: 燃料/载荷比 ==========
|
||
ax2 = axes[0, 1]
|
||
for engine_name, isp in list(ENGINES.items())[:3]: # 排除离子推进
|
||
ratios = [fuel_mass(p, DELTA_V_TOTAL, isp) / p for p in payloads]
|
||
ax2.axhline(y=ratios[0], label=f'{engine_name} (Isp={isp}s)', linewidth=2)
|
||
|
||
ax2.set_xlabel('有效载荷质量 (吨)', fontsize=12)
|
||
ax2.set_ylabel('燃料/载荷 质量比', fontsize=12)
|
||
ax2.set_title('燃料效率比较\n(质量比与载荷无关，仅取决于Isp和ΔV)', fontsize=14)
|
||
ax2.legend()
|
||
ax2.grid(True, alpha=0.3)
|
||
ax2.set_ylim(0, 10)
|
||
|
||
# ========== 图3: 能量消耗 ==========
|
||
ax3 = axes[1, 0]
|
||
for engine_name, isp in list(ENGINES.items())[:3]:
|
||
fuels = np.array([fuel_mass(p, DELTA_V_TOTAL, isp) for p in payloads])
|
||
energies = total_energy(fuels, isp) / 1e12 # 转换为TJ
|
||
valid = energies < 1e6
|
||
ax3.plot(payloads[valid]/1000, energies[valid],
|
||
label=f'{engine_name}', linewidth=2)
|
||
|
||
ax3.set_xlabel('有效载荷质量 (吨)', fontsize=12)
|
||
ax3.set_ylabel('能量消耗 (TJ)', fontsize=12)
|
||
ax3.set_title('能量消耗 vs 有效载荷\nE = ½ × m_fuel × ve²', fontsize=14)
|
||
ax3.legend()
|
||
ax3.grid(True, alpha=0.3)
|
||
|
||
# ========== 图4: ΔV 敏感性分析 ==========
|
||
ax4 = axes[1, 1]
|
||
delta_v_range = np.linspace(3000, 5000, 100) # 3-5 km/s
|
||
payload_fixed = 10000 # 10吨有效载荷
|
||
isp_reference = 450 # 液氧液氢
|
||
|
||
fuel_vs_dv = [fuel_mass(payload_fixed, dv, isp_reference) for dv in delta_v_range]
|
||
ax4.plot(delta_v_range/1000, np.array(fuel_vs_dv)/1000,
|
||
'b-', linewidth=2, label='燃料质量')
|
||
|
||
# 标记关键点
|
||
ax4.axvline(x=DELTA_V_TOTAL/1000, color='r', linestyle='--',
|
||
label=f'标准任务 ΔV={DELTA_V_TOTAL/1000} km/s')
|
||
|
||
ax4.set_xlabel('ΔV (km/s)', fontsize=12)
|
||
ax4.set_ylabel('燃料质量 (吨)', fontsize=12)
|
||
ax4.set_title(f'ΔV敏感性分析\n(有效载荷=10吨, Isp=450s)', fontsize=14)
|
||
ax4.legend()
|
||
ax4.grid(True, alpha=0.3)
|
||
|
||
plt.tight_layout()
|
||
plt.savefig('/Volumes/Files/code/mm/20260130_b/earth_moon_transfer_analysis.png',
|
||
dpi=150, bbox_inches='tight')
|
||
print("图表已生成!")
|
||
# plt.show() # 非交互式环境下注释掉
|
||
|
||
return fig
|
||
|
||
|
||
def print_summary():
|
||
"""
|
||
打印关键数据摘要
|
||
"""
|
||
print("=" * 60)
|
||
print("地月转移轨道分析摘要")
|
||
print("=" * 60)
|
||
print(f"\n任务参数:")
|
||
print(f" - 总ΔV需求: {DELTA_V_TOTAL/1000:.1f} km/s")
|
||
print(f" - TLI (地月转移注入): {DELTA_V_TLI/1000:.1f} km/s")
|
||
print(f" - LOI (月球轨道捕获): {DELTA_V_LOI/1000:.1f} km/s")
|
||
print(f" - 结构系数 α: {ALPHA}")
|
||
|
||
print(f"\n不同发动机的燃料效率 (每吨有效载荷所需燃料):")
|
||
print("-" * 50)
|
||
|
||
payload_ref = 1000 # 1吨参考载荷
|
||
|
||
for engine_name, isp in ENGINES.items():
|
||
R = mass_ratio(DELTA_V_TOTAL, isp)
|
||
fuel = fuel_mass(payload_ref, DELTA_V_TOTAL, isp)
|
||
ve = exhaust_velocity(isp)
|
||
|
||
if fuel < 1e9:
|
||
print(f"\n {engine_name}:")
|
||
print(f" 比冲 Isp = {isp} s")
|
||
print(f" 排气速度 ve = {ve:.0f} m/s")
|
||
print(f" 质量比 R = {R:.2f}")
|
||
print(f" 燃料/载荷比 = {fuel/payload_ref:.2f}")
|
||
print(f" 1吨载荷需燃料: {fuel/1000:.2f} 吨")
|
||
else:
|
||
print(f"\n {engine_name}:")
|
||
print(f" 比冲 Isp = {isp} s (化学火箭无法实现)")
|
||
|
||
print("\n" + "=" * 60)
|
||
print("关键数学关系:")
|
||
print("=" * 60)
|
||
print("""
|
||
齐奥尔科夫斯基方程:
|
||
ΔV = ve × ln(m₀/mf)
|
||
|
||
质量比:
|
||
R = exp(ΔV / ve)
|
||
|
||
燃料-载荷关系:
|
||
m_fuel / m_payload = (R - 1) / (1 - α(R - 1))
|
||
|
||
能量消耗:
|
||
E = ½ × m_fuel × ve²
|
||
|
||
其中:
|
||
ve = Isp × g₀ (排气速度)
|
||
α = 结构系数 (约0.1)
|
||
""")
|
||
|
||
|
||
if __name__ == "__main__":
|
||
# 打印分析摘要
|
||
print_summary()
|
||
|
||
# 绘制关系图
|
||
print("\n正在生成图表...")
|
||
fig = plot_payload_fuel_relationship()
|
||
print("图表已保存至: earth_moon_transfer_analysis.png")
|