Course: ASEN 6014 — Space Vehicle Guidance and Control, CU Boulder
Topics: Relative orbital elements (ROE), eccentric orbit formation control, impulsive maneuver design, delta-v lower bound, J2/SRP perturbation analysis, Monte Carlo robustness, boundedness assessment
This project implements and evaluates the closed-form impulsive controller for spacecraft formation flying in eccentric orbits derived by Chernick & D'Amico (2018). The controller analytically computes the locations and magnitudes of tangential and normal burns that drive a deputy spacecraft from an initial ROE state to a desired target, using the Eccentric Relative Orbital Elements (EROE) framework.
The implementation is validated in three layers of fidelity:
| Layer | Model | Purpose |
|---|---|---|
| Linear | ROE STM only | Verify controller algebra |
| Two-body truth | Cartesian ECI + COE conversions | Check nonlinear fidelity |
| EROE truth | Linear EROE plant with J2 + SRP | Perturbation sensitivity + Monte Carlo |
Chief orbit: a=10,000 km, e=0.2, i=30°, Ω=40°, ω=25°, M₀=10°
Burn schedule (one-orbit control window, δaF = 0 for boundedness):
| Burn | Type | Mean anomaly [rad] | Time [s] | ΔvT [m/s] | ΔvN [m/s] |
|---|---|---|---|---|---|
| 1 | Tangential | 1.745 | 2,488 | −2.90 | 0 |
| 2 | Tangential | 3.316 | 4,976 | −2.52 | 0 |
| 3 | Tangential | 4.887 | 7,464 | +4.83 | 0 |
| 4 | Normal | 6.361 | 9,799 | 0 | +3.75 |
Total ΔvΣ = 14.0 m/s vs theoretical lower bound Δv_LB = 7.7 m/s → ratio ≈ 1.81 (expected for a fixed 3-burn triplet vs a continuously-distributed optimum)
Figure 1: Per-burn and cumulative control effort. Burns 1-3 are tangential; burn 4 is cross-track.
Figure 2: Relative trajectory in the RT plane during the one-orbit control window. Red markers show the four burn locations.
Four environment models were compared: Kepler only, Kepler+J2, Kepler+SRP, and Kepler+J2+SRP. The controller drives all cases to nearly the same final ROE; perturbations affect only the long-term free-drift behavior.
Figure 3: ROE time histories across all four perturbation models. Burns occur during the first ~2.5 hr. J2 introduces bounded deviations in δλ; SRP drives secular drift in δe.
Figure 4: RT-plane trajectories over 8 free-drift orbits after the control window. Kepler produces a nearly closed ellipse; J2 and SRP slowly deform it.
Key finding: J2 induces bounded periodic modulation in δλ and the eccentricity vector. SRP produces secular growth in δλ due to the constant differential solar forcing. Combined J2+SRP shows the largest long-term deviation but the final controlled state differs from target only at the 10⁻⁵–10⁻⁶ level.
Phase 3A: Navigation noise (1σ ROE estimation errors applied to both δα₀ and δαF)
| Component | RMS error |
|---|---|
| δa | 7.1 × 10⁻⁵ |
| δλ | 6.7 × 10⁻³ |
| δex | 9.7 × 10⁻⁴ |
| δey | 1.5 × 10⁻³ |
| δix | 1.4 × 10⁻⁴ |
| δiy | 1.4 × 10⁻⁴ |
Mean total Δv: 14.9 m/s (RMS deviation: 3.6 m/s from nominal)
Figure 5: Phase 3A — distribution of final ROE error norm over 5,000 navigation-noise runs. Broad spread dominated by δλ and δe components.
Phase 3B: Thruster magnitude error (±5% 1σ on each burn)
| Component | RMS error |
|---|---|
| δex | 5.5 × 10⁻⁴ |
| δey | 1.15 × 10⁻³ |
| δix | 2.1 × 10⁻⁵ |
| δiy | 1.3 × 10⁻⁵ |
Mean total Δv: 14.0 m/s (RMS deviation: 0.36 m/s — nearly unchanged)
Figure 6: Phase 3B — distribution of final ROE error norm under 5% thruster errors. Distribution is tighter than Phase 3A, confirming navigation is the dominant error source.
Navigation uncertainty is the dominant robustness concern: it forces the controller to recompute all burns from noisy data, producing 10× larger final-state errors and Δv scatter than thruster errors alone.
The nominal final state (δaF ≈ 0) stays bounded; a worst-case 3σ bias in δa drives secular growth exceeding 2,000 km relative distance within 50 orbits.
Figure 7: 50-orbit free drift. Nominal (blue solid): δa ≈ 0, bounded relative distance. Worst-case (red dashed): 3σ δa bias drives monotonic drift in δλ and unbounded relative distance.
This confirms the theoretical requirement: δaF = 0 is essential for long-term bounded formation flight under J2+SRP.
The 6-element state vector is:
δα = [δa, δλe, δex, δey, δix, δiy]
where δλe = (Md − Mc) + η(Δω + ΔΩ cos i) is the eccentric mean longitude (differs from the EROE dlam by an η factor).
Key properties:
- δa = 0 → bounded (no secular along-track drift)
- (δex, δey) → in-plane eccentricity vector; controls relative ellipse shape
- (δix, δiy) → cross-track inclination vector; controlled by a single normal burn
The Γ matrix maps an RTN impulse Δv to an instantaneous ROE jump:
Δδα = Γ(ν) · Δv
Scaled by 1/(n·a), so Γ · dv [km/s] = δα [dimensionless]. At perigee (ν=0), the tangential column has its maximum effect on δa (G₁₂ = 2(1+e)/η).
Notable structure: Tangential burns have zero direct effect on δλe (G₂₂ = 0). The three-burn in-plane solution therefore leaves δλe unconstrained — this is expected and explains the large δλe mismatch in the truth model.
Three tangential burns at mean anomalies M₁, M₂, M₃ (evenly spaced in the control window) solve a 3×3 linear system for Δv₁, Δv₂, Δv₃ targeting (δa, δex, δey). The validity condition (eccentricity-dominated regime) requires:
||Δδe|| ≥ max(|Δδa|·Ce / [2(1+e)], |Δδλe|·Ce / [3(1+e)·ΔM])
where Ce = √(3e⁴ − 7e² + 4) → 2 as e → 0.
A single normal burn at argument of latitude θN = atan2(Δδiy, Δδix):
|ΔvN| = n·a · (1 + e·cos(νN))/η · ||Δδi||
The (1 + e·cos(ν)) factor captures how orbital speed modulates normal thrust efficiency in eccentric orbits.
Δv_LB / (n·a·η) = max(|Δδa| / [2(1+e)], |Δδλe| / [3(1+e)·ΔM], ||Δδe|| / Ce)
For this test case (Δδa = 0, Δδλe = 0), the eccentricity term dominates and Δv_LB = 7.7 m/s. The 3-burn solution achieves ~1.81× this bound — within the expected range for a fixed-location triplet.
eccentric-orbit-formation-control/
├── README.md
├── .gitignore
├── src/
│ ├── test_cd_with_eroe_truth.m # Main script: all 4 phases (run this first)
│ ├── test_ecc_controller_truth.m # Phase 2: two-body Cartesian truth validation
│ ├── step1_test_script.m # Sanity check: STM consistency
│ ├── step2_test_script.m # Linear-model controller validation
│ │
│ ├── cd_eccentric_controller.m # Top-level controller wrapper
│ ├── cd_eccentric_controller_refined.m # Refined controller (Newton hook, see note)
│ ├── cd_ecc_inplane_controller.m # In-plane: regime check + burn placement
│ ├── cd_ecc_inplane_triplet.m # 3x3 linear solve for tangential burns
│ ├── cd_ecc_outofplane.m # Single-burn normal maneuver
│ ├── cd_eq44_residual_and_jac.m # Placeholder for Eq. 44 Newton refinement
│ ├── cd_dv_lower_bound.m # Analytical delta-v lower bound (Eq. 20)
│ ├── build_maneuver_list_from_cd_controller.m # Assemble sorted burn list
│ ├── ecc_eccentric_dominated_regime.m # Regime validity check
│ │
│ ├── gamma_eccentric_cd.m # ROE control input matrix Γ(ν)
│ ├── roe_stm_eccentric.m # Eccentric ROE state transition matrix
│ ├── roe_propagate_eccentric.m # ROE propagation via STM
│ ├── roe_from_coe_cd.m # COE -> CD ROE conversion
│ ├── eroe_propagate_with_maneuvers.m # EROE dynamics + impulsive jumps + J2/SRP
│ │
│ ├── coe2rv_eci.m # COE -> ECI Cartesian
│ ├── rv2coe_eci.m # ECI Cartesian -> COE
│ ├── apply_impulse_rtn.m # Apply RTN impulse to ECI state
│ ├── rtn_frame.m # Compute RTN unit vectors from r, v
│ ├── true_anomaly_from_mean.m # Kepler solver: M -> ν
│ └── mean_from_true.m # Inverse: ν -> M
│
└── figures/
├── fig01_control_effort.png # Per-burn and cumulative dv
├── fig02_rt_plane_burns.png # RT plane with burn markers
├── fig03_roe_time_histories.png # ROE components vs time (4 pert. models)
├── fig04_rt_plane_perturbations.png # RT plane, all perturbation models
├── fig05_monte_carlo_nav_noise.png # Phase 3A histogram
├── fig06_monte_carlo_thruster_error.png # Phase 3B histogram
└── fig07_boundedness.png # Phase 4: 50-orbit drift
Add src/ to your MATLAB path, then call:
addpath('src')
% Full simulation (all phases, ~2 min for Monte Carlo)
test_cd_with_eroe_truth()
% Two-body Cartesian truth validation
test_ecc_controller_truth()
% Sanity check: STM consistency
step1_test_script()
% Linear-model controller validation with plots
step2_test_script()Perturbation toggles (set at the top of test_cd_with_eroe_truth.m):
cases(1).opts = struct('includeJ2', false, 'includeSRP', false); % Kepler
cases(2).opts = struct('includeJ2', true, 'includeSRP', false); % +J2
cases(3).opts = struct('includeJ2', false, 'includeSRP', true, 'srp', srp_opts); % +SRP
cases(4).opts = struct('includeJ2', true, 'includeSRP', true, 'srp', srp_opts); % +J2+SRPMonte Carlo sample count (at top of Phase 3A/3B sections):
N_mc_nav = 5000; % navigation noise runs
N_mc_thr = 5000; % thruster error runs| Parameter | Value |
|---|---|
| Semi-major axis a | 10,000 km |
| Eccentricity e | 0.20 |
| Inclination i | 30° |
| RAAN Ω | 40° |
| Arg. of perigee ω | 25° |
| Epoch mean anomaly M₀ | 10° |
| μ | 398,600 km³/s² |
| J₂ | 1.08262668 × 10⁻³ |
| SRP ΔB | 3 × 10⁻⁴ m²/kg |
| Control window | 2π rad (one full orbit) |
| Nav noise 1σ | [5×10⁻⁵, 5×10⁻³, 5×10⁻⁴, 5×10⁻⁴, 10⁻⁴, 10⁻⁴] |
| Thruster error 1σ | 5% (multiplicative, per burn) |
cd_eq44_residual_and_jac is a documented placeholder that returns F = 0, J = I. cd_eccentric_controller_refined therefore produces the same result as the basic controller for this implementation. The placeholder structure is preserved so the Newton refinement from Eq. (44) of Chernick & D'Amico can be added without restructuring the codebase.
ROE convention: roe_from_coe_cd uses the Chernick–D'Amico δλe definition, which differs from the Delurgio–D'Amico EROE dlam by an η factor. For bounded formations (δa = 0), the coupling term in the STM vanishes and the two conventions produce identical trajectories in the ROE time-series plots.
eroe_propagate_with_maneuvers bridges the two conventions: it initializes the state from CD ROE, uses the EROE Kepler plant matrix for continuous propagation, and applies impulsive jumps via gamma_eccentric_cd. This is an approximation that is accurate for the bounded (δa ≈ 0) regime explored in Phases 3–4.
Chernick, M. & D'Amico, S. (2018). New Closed-Form Solutions for Optimal Impulsive Control of Spacecraft Relative Motion. Journal of Guidance, Control, and Dynamics, 41(2), 301–319. https://doi.org/10.2514/1.G002814
- MATLAB R2020b or later
- No toolboxes required — uses only built-in
ode45/expm/norm/cross eroe_propagate_with_maneuverscallseroe_A_J2_numeric,eroe_A_SRP, anderoe_series_to_rtnfrom Project 1 (formation-flying-eccentric-orbits/src/) if J2 or SRP are enabled. Copy those three files intosrc/or add the Project 1src/folder to your path.