# Moon Colony Logistics（地月运输建模与对比）
本目录包含针对 **MCM 2026 Problem B**（见 `prob/en_B.md`）的一组可运行脚本：用**简化的轨道力学 + 火箭方程**估算把 **1 亿吨（100 million metric tons）**物资送到月球轨道/基地的**时间线**与**能量消耗**，并对比：
- **纯火箭方案**（10 个既有发射场，每站每日 1 发）
- **纯太空电梯方案**（3 个电梯港口持续转运）
- **混合方案**（电梯优先，剩余由火箭补齐，且随着工期变长更多任务转移到低纬发射场以降低纬度惩罚）
> 本项目的“能量”口径：
> - **火箭**：用推进剂化学能近似 \(E \approx m_\text{fuel}\cdot e_s\)（`e_s` 为燃料比能，默认 LH2/LOX 约 15.5 MJ/kg）。
> - **太空电梯**：用“提升能 + 顶端火箭段能量”的**总比能量**近似（默认 157.2 MJ/kg，即 157.2 GJ/吨）。
> 这是一种用于**对比量级与趋势**的工程近似，不等同于严格的热力学效率或完整任务能量收支。
---
### 目录结构与脚本索引
#### 核心脚本
- **`specific_energy_comparison.py`**
地面发射 vs 太空电梯：比能量、燃料比、纬度影响（包含指定发射场列表）。
- **`rocket_launch_scenario.py`**
纯火箭方案（10 发射场 * 每日 1 发）排布与能量-工期曲线；并做发射频率敏感性分析。
- **`combined_scenario.py`**
混合方案（电梯优先 + 火箭补齐），输出组合方案能量-工期曲线，并生成 **100–300 年**区间三方案同图对比。
- **`pareto_optimization.py`**（新增）
**Pareto 前沿 + “膝点（knee point）”**分析：在**组合方案可行区间（约 101–186 年）**内，寻找“时间–能量”双目标权衡的拐点，并输出推荐决策点与图表/CSV。
- **`launch_capacity_analysis.py`**
（可选）对历史火箭年发射量做 Logistic / Gompertz / Richards 拟合与预测，生成 `launch_capacity_*.png` 与 `launch_capacity_results.csv`。
#### 主要输出图表（均在本目录）
| 图表文件 | 含义 |
|---|---|
| `specific_energy_comparison.png` | 地面 vs 电梯的比能量/燃料比对比，含 ΔV 敏感性 |
| `latitude_effects.png` | 月球任务场景下：纬度→自转速度损失→燃料/能量增加趋势（并标注指定发射场） |
| `rocket_launch_timeline.png` | 纯火箭方案：完成年限 vs 总能量（工期变长→更多低纬→能量略降） |
| `rocket_comprehensive.png` | 纯火箭方案综合图：能量趋势 + 分配饼图等 |
| `launch_frequency_analysis.png` | 发射频率敏感性：每站每日发射数提升→完工年限显著下降 |
| `combined_scenario_analysis.png` | 混合方案：总能量/电梯能量/火箭能量随工期变化 + 电梯占比 + 需要启用的火箭站点数 |
| `scenario_comparison.png` | elevator-only / rocket-only / combined 的完工时间与能量柱状对比 |
| **`three_scenarios_comparison.png`** | **100–300 年区间**：纯火箭、纯电梯、混合方案的“年份-能量”折线同图对比 |
| **`pareto_combined_range.png`** | **组合方案（101–186 年）Pareto 前沿 + 膝点检测**（多方法对比、边际收益、构成变化） |
| **`combined_decision.png`** | 组合方案（101–186 年）推荐拐点的“决策图”（更适合直接放报告） |
| `launch_capacity_analysis.png` `launch_capacity_physical.png` `launch_capacity_predictions.png` | （可选）火箭发射能力上限/趋势拟合与预测 |
#### 数据文件（如存在）
- `launch_capacity_results.csv`
- `rocket_launch_counts.csv`
- `launch_sites.csv`
- `rocket launch counts.txt`
---
### 图表总览（本目录所有 PNG）
为便于写报告/插图，这里把当前目录下的所有 PNG 统一列出（按主题分组）。如果你只想看关键结论，优先看 **`three_scenarios_comparison.png`**。
1) 比能量与纬度效应（地面 vs 电梯、发射场纬度）
#### `specific_energy_comparison.png`
#### `latitude_effects.png`
2) 纯火箭方案（排布与敏感性）
#### `rocket_launch_timeline.png`
#### `rocket_comprehensive.png`
#### `launch_distribution_min.png`
#### `launch_distribution_mid.png`
#### `launch_frequency_analysis.png`
3) 混合方案与三方案对比（核心结论图）
#### `combined_scenario_analysis.png`
#### `scenario_comparison.png`
#### `three_scenarios_comparison.png`（100–300 年折线同图对比）
4)（可选）火箭发射能力上限/趋势拟合
#### `launch_capacity_analysis.png`
#### `launch_capacity_physical.png`
#### `launch_capacity_predictions.png`
5) Pareto 前沿与“膝点”（组合方案 101–186 年）
> 说明：`pareto_optimization.py` 与其输出（图片/CSV）均在本目录中。
#### `pareto_combined_range.png`
#### `combined_decision.png`
6) 其它历史遗留图
本目录还包含 `earth_moon_transfer_analysis.png`（早期“地月转移：载荷–燃料/能量关系”图）。当前目录下未发现对应的 `earth_moon_transfer.py` 源码文件；如果你需要我把该脚本补回并与现有口径统一，我可以再加一个文件。
#### `earth_moon_transfer_analysis.png`
### 依赖安装
建议使用 Python 3.10+。
```bash
pip install numpy matplotlib pandas
```
如需运行 `launch_capacity_analysis.py`（含曲线拟合），还需要：
```bash
pip install scipy
```
---
### 数学模型（简化版）概览
#### 1）火箭方程：燃料-载荷的线性比例（在固定 ΔV 与 Isp 下）
火箭方程：
\[
\Delta V = v_e \ln\left(\frac{m_0}{m_f}\right),\quad v_e = I_{sp} g_0
\]
当任务 ΔV、发动机 \(I_{sp}\) 固定时，质量比 \(R=\exp(\Delta V/v_e)\) 为常数。若进一步采用结构系数近似：
\[
m_s=\alpha m_\text{fuel}
\]
则可解得燃料与载荷满足：
\[
\boxed{m_\text{fuel}=k\,m_p}\quad(\text{其中 }k\text{ 为常数，取决于 }\Delta V, I_{sp}, \alpha)
\]
这解释了为什么在同一任务与发动机参数下，“载荷–燃料”是近似**线性**关系。
#### 2）能量口径（用于对比）
- 火箭推进剂能量（近似）：
\[
E_\text{rocket}\approx m_\text{fuel}\,e_s
\]
其中 \(e_s\) 为燃料比能（默认 15.5 MJ/kg）。
- 太空电梯总比能量（本项目用常数近似）：
\[
E_\text{elevator}\approx \epsilon_\text{elevator}\,m_p
\]
默认 \(\epsilon_\text{elevator}=157.2\text{ MJ/kg}\)（可在脚本中改）。
#### 3）纬度导致的自转“损失”近似
地球表面自转线速度：
\[
v(\varphi)=\omega R_E\cos(\varphi)
\]
相对赤道的速度差近似为 ΔV 惩罚：
\[
\Delta v_\text{loss}=v(0^\circ)-v(\varphi)
\]
在月球任务场景下（允许进入与发射场纬度相匹配的倾角、避免昂贵的纯平面变轨），纬度主要通过该项影响火箭能耗。
---
### 如何运行（会生成图片）
#### 1）比能量与纬度效应
```bash
python specific_energy_comparison.py
```
输出（示例）会包括不同发射场的 ΔV 损失、燃料比、能量（MJ/kg）等，并生成：
- `specific_energy_comparison.png`
- `latitude_effects.png`
#### 2）纯火箭方案：工期-能量与排布
```bash
python rocket_launch_scenario.py
```
会生成：
- `rocket_launch_timeline.png`
- `rocket_comprehensive.png`
- `launch_distribution_min.png`、`launch_distribution_mid.png`
- `launch_frequency_analysis.png`
**部分运行输出（示例）**：
```text
Total payload: 100 million metric tons
Payload per launch: 125 metric tons
Total launches needed: 800,000
Number of launch sites: 10
Minimum completion time: 219.2 years
Energy reduction (Extended vs Minimum time): 2.4%
```
解释：
- 由于“每站每日 1 发、载荷 125 吨”的硬约束，纯火箭最短也要约 **219 年**完成。
- 工期拉长后，可把更多发射集中到低纬站点，**能量仅小幅下降（约几个百分点）**——纬度惩罚相对火箭整体能耗并不占主导。
#### 3）混合方案 + 100–300 年三方案同图对比
```bash
python combined_scenario.py
```
会生成：
- `combined_scenario_analysis.png`
- `scenario_comparison.png`
- **`three_scenarios_comparison.png`（100–300 年折线同图对比）**
**部分运行输出（示例）**：
```text
Elevator-only completion: 186.2 years
Minimum completion (all systems): 100.7 years
1. Minimum Time (~101 years): Total energy: 31537 PJ
4. Elevator Only (186 years): Total energy: 15720 PJ
Rocket Only (219 years): Total energy: 50609 PJ
```
解释：
- **混合方案**把最短工期从电梯单独的 ~186 年压缩到 ~101 年（靠火箭补齐），但能量高于纯电梯。
- 当工期 ≥ 186 年时，电梯即可独立完成，混合方案会自动退化为纯电梯曲线（能耗趋于同一水平）。
- 纯火箭能耗显著更高（本模型下约为纯电梯的 **3×** 量级）。
#### 4）Pareto 前沿与“膝点”（只在组合方案 101–186 年区间内找拐点）
> 这个分析回答“在可行区间内，选择多少年最划算”的决策问题：既不想太慢、也不想能耗太高。
> 我们把目标写成双目标：**最小化完工时间 \(T\)** 与 **最小化总能量 \(E\)**，并在 \(T\in[101, 186]\) 年上寻找“边际收益开始明显递减”的拐点（knee point）。
在本目录运行：
```bash
python pareto_optimization.py
```
会生成：
- `pareto_combined_range.png`
- `combined_decision.png`
- `pareto_combined_range.csv`（组合区间的 Pareto 数据）
- `pareto_knee_analysis.png`（完整 100–300 年范围的膝点分析）
- `decision_recommendation.png`（完整范围的决策建议图）
- `pareto_front_data.csv`（完整范围数据）
**组合区间（101–186 年）推荐膝点（示例输出）**：
| 推荐方法 | 年份 | 总能量 (PJ) | 电梯占比 | 火箭占比 |
|---|---:|---:|---:|---:|
| 最大距离法（最稳健） | **139** | **24361** | **74.6%** | **25.4%** |
解读要点：
- **101→139 年**：多给一些时间，能量下降很快（大量火箭任务转移到电梯/低纬火箭站）。
- **139→186 年**：继续拉长工期仍能省能量，但**边际节省明显变小**；186 年达到纯电梯可独立完成的极限点。
---
### 发射场（10 个）与纬度效应（月球任务场景）
本项目按纬度从低到高优先分配火箭发射（“低纬优先”贪心策略），站点列表与纬度如下（与题目一致）：
| 发射场 | 近似纬度 | 备注 |
|---|---:|---|
| French Guiana（Kourou） | 5.2° | 低纬优势明显 |
| Satish Dhawan（India） | 13.7° | 低纬 |
| Texas（Boca Chica） | 26.0° | 美国 |
| Florida（Cape Canaveral） | 28.5° | 美国 |
| California（Vandenberg） | 34.7° | 美国 |
| Virginia（Wallops） | 37.8° | 美国 |
| Taiyuan（China） | 38.8° | 中国 |
| Mahia（New Zealand） | 39.3° | 新西兰（南半球，按 \(|\varphi|\) 处理） |
| Kazakhstan（Baikonur） | 45.6° | 哈萨克斯坦 |
| Alaska（Kodiak） | 57.4° | 高纬，惩罚最大 |
对应曲线可在下图直观看到：
---
### 三方案“年份–能量”折线对比（100–300 年）
下图是你要求的 **100–300 年**区间三条折线同图（不可完成的年份用空缺表示）：
读图要点：
- **100–186 年**：只有混合方案可完成（电梯能力不足、火箭也不够快）。
- **186–219 年**：电梯可完成；火箭仍无法完成；混合与电梯曲线重合（优先用电梯）。
- **≥219 年**：三种方案都可完成，但纯火箭能耗最高。
为了写报告时更快引用，下面给出几组关键年份的能量读数（与脚本运行输出一致，单位 PJ）：
| 年份 | 混合方案 | 纯电梯 | 纯火箭 |
|---:|---:|---:|---:|
| 110 | 29852 | N/A | N/A |
| 150 | 22257 | N/A | N/A |
| 186 | 15720 | 15720 | N/A |
| 220 | 15720 | 15720 | 50605 |
| 300 | 15720 | 15720 | 50187 |
---
### 可调参数（最常改的地方）
#### 火箭相关（`rocket_launch_scenario.py` / `combined_scenario.py`）
- `TOTAL_PAYLOAD`：总物资（吨）
- `PAYLOAD_PER_LAUNCH`：单次发射载荷（吨，默认 125）
- `ISP`、`ALPHA`、`NUM_STAGES`：推进系统与结构参数
- `SPECIFIC_FUEL_ENERGY`：燃料比能（J/kg）
- `DELTA_V_BASE`：赤道基准 ΔV（m/s）
- `LAUNCH_SITES`：发射场纬度列表（可自行替换/增减）
#### 太空电梯相关（`combined_scenario.py`）
- `NUM_ELEVATORS`：电梯数量（默认 3）
- `ELEVATOR_CAPACITY_PER_YEAR`：**每个电梯年运力（吨/年）**
> 题面是 **179,000 t/year**；若你要用 **17,900 t/year**，直接改此值即可。
- `ELEVATOR_SPECIFIC_ENERGY`：电梯总比能量（J/吨），默认 157.2 GJ/吨
---
### 重要假设与局限（写报告时建议明确）
- **脉冲变轨/固定 ΔV**：将复杂轨道机动用固定 ΔV 预算概括。
- **恒定比冲/结构系数**：忽略发动机工况变化、分级差异、载具设计差异。
- **“能量”是工程近似口径**：火箭取推进剂化学能；电梯取提升+顶端段合并比能量。
- **完美运行**：未考虑失败率、维护停机、窗口约束、天气/法规等运营因素。
- **排布算法为贪心**：火箭发射按“低纬优先”分配以降低纬度惩罚；未做更复杂的成本/容量优化（可扩展为线性规划/整数规划）。
---
### 一句话结论（便于写 Summary）
- **纯火箭**在“10 站 * 每日 1 发 * 125 吨/发”的约束下，最短约 **219 年**，且能耗约 **5×10⁴ PJ** 量级。
- **纯电梯**能耗最低（约 **1.57×10⁴ PJ**），但最短约 **186 年**。
- **混合方案**能把工期压到约 **101 年**，能耗介于两者之间，并随着工期增加逐步趋近纯电梯能耗。
- **在混合方案可行区间（101–186 年）内**，用 Pareto/膝点法可给出一个“最划算”的折中点：本模型下约 **139 年**（约 **2.44×10⁴ PJ**，电梯约 **75%**、火箭约 **25%**）。