20 KiB
Patched Conics and SOI Transition Implementation Plan
Date: January 14, 2026 Status: Planning Phase Branch: To be created
Overview
This plan implements support for patched conics trajectory simulation, enabling satellites to transition between multiple spheres of influence (SOI) in complex orbital scenarios:
- Planet → Star → Planet transfers
- Star → Planet → Moon rendezvous
- Multi-leg interplanetary missions
Current State Analysis
✅ What's Already Working
-
SOI Transitions Are Already Implemented
- Lines 103-121 in
simulation.cppshow coordinate transformation logic - Converts local→global using old parent, then global→local using new parent
- This is essentially Phase 3 implementation from hierarchical_frames_plan.md
- Lines 103-121 in
-
Local Frame Integration (Phase 2) ✅
- All bodies integrate in local coordinates
- Global coordinates computed after each timestep
- Improved numerical precision for nested orbits
-
Parent-First Update Order (Phase 4) ✅
- Root bodies skipped in loop (fixed at origin)
- Child bodies integrate using parent coordinates
- Hierarchical update order implemented
🔴 Critical Issues for Patched Conics
Issue 1: Cannot Transition to Root Bodies
if (new_parent != body->parent_index && new_parent != -1) {
// ... transition logic
}
The new_parent != -1 condition prevents switching to Sun (parent_index = -1). This breaks the scenario: Planet→Sun→Moon rendezvous is impossible.
Issue 2: Hysteresis Barrier
The 0.5x hysteresis factor in find_dominant_body() (line 71) creates one-way barriers:
- Planet→Sun: Possible (easy entry)
- Sun→Planet: Impossible (can't exit due to hysteresis)
Issue 3: Integration After Transition Transition happens before integration in the same timestep, using coordinates from the end of the previous timestep. This may cause velocity discontinuities.
⚠️ Potential Issues
Issue 4: Numerical Precision Satellites crossing between star/planet/moon scales will see position magnitude changes of 10⁸ to 10¹¹ meters, potentially losing precision.
Issue 5: Fixed Timestep 60s timestep may be too coarse for fast orbital phases (moon capture) and too slow for deep-space phases.
Implementation Phases
Phase 1: Fix Root Body Transitions (Critical)
Goal: Allow satellites to switch to/from Sun (parent_index = -1)
Changes:
-
Remove
new_parent != -1check insimulation.cppline 104 -
Add special handling for root body transitions:
if (new_parent != body->parent_index) { // Convert local → global using old parent if (body->parent_index >= 0) { // old_parent is a real body CelestialBody* old_parent = &sim->bodies[body->parent_index]; body->position = vec3_add(body->local_position, old_parent->position); body->velocity = vec3_add(body->local_velocity, old_parent->velocity); } else { // old_parent is root (Sun): local = global body->position = body->local_position; body->velocity = body->local_velocity; } body->parent_index = new_parent; // Convert global → local using new parent if (new_parent >= 0) { // new_parent is a real body CelestialBody* new_parent_body = &sim->bodies[new_parent]; body->local_position = vec3_sub(body->position, new_parent_body->position); body->local_velocity = vec3_sub(body->velocity, new_parent_body->velocity); } else { // new_parent is root (Sun): global = local body->local_position = body->position; body->local_velocity = body->velocity; } } -
Update
find_dominant_body()to properly handle -1 returns
Files to modify:
src/simulation.cpp(lines 103-121)src/simulation.h(no changes needed)
Estimated complexity: Low Risk: Medium (affects core transition logic)
Expected outcome:
- ✅ Satellites can transition to/from Sun
- ✅ Enables Planet→Sun→Planet transfers
- ✅ Enables Star→Planet→Moon rendezvous
Tests to add:
- Test satellite transitioning from Earth to Sun
- Test satellite transitioning from Sun to Mars
- Test full round-trip: Earth→Sun→Earth
Phase 2: Remove or Modify Hysteresis (Critical for Round-Trips)
Current Problem: The 0.5x hysteresis factor prevents oscillation but creates one-way barriers:
if (distance < dist_to_current * 0.5) {
dominant = i;
}
Option A: Remove Hysteresis
- Remove 0.5x factor (line 71 in
simulation.cpp) - Allow switching to closest body at all times
- Pros: Simple, enables all transitions
- Cons: May cause oscillation at SOI boundaries
Option B: Adaptive Hysteresis (Recommended)
- Keep hysteresis but only apply when already in SOI
- Allow free switching when outside current SOI
- Pros: Prevents oscillation while enabling round-trips
- Cons: More complex logic
Option B Implementation:
if (can_switch && i != dominant) {
if (dominant == -1) {
dominant = i;
} else {
CelestialBody* current_parent = &sim->bodies[dominant];
double dist_to_current = vec3_distance(body->position, current_parent->position);
if (outside_current_soi) {
// Outside current SOI: switch to closest body (no hysteresis)
if (distance < dist_to_current) {
dominant = i;
}
} else {
// Inside current SOI: apply hysteresis to prevent oscillations
if (distance < dist_to_current * 0.5) {
dominant = i;
}
}
}
}
Files to modify:
src/simulation.cpp(line 70-75)
Estimated complexity: Low-Medium Risk: Medium (may affect transition behavior)
Expected outcome:
- ✅ Enables round-trip transitions (Earth→Sun→Earth)
- ✅ Maintains stability by preventing oscillation when inside SOI
- ✅ Allows free switching when outside current SOI
Tests to add:
- Validate Satellite→Planet→Satellite round-trip
- Validate Earth→Sun→Mars→Sun→Earth full round-trip
Phase 3: Refactor Transition to Separate Function
Goal: Cleaner code, easier to test, better separation of concerns
Current state: Transition logic is inline in update_simulation() (lines 105-120)
Proposed refactoring:
Add to simulation.h:
void transition_to_new_parent(SimulationState* sim, CelestialBody* body,
int old_parent_idx, int new_parent_idx);
Add to simulation.cpp:
void transition_to_new_parent(SimulationState* sim, CelestialBody* body,
int old_parent_idx, int new_parent_idx) {
// Current state is in old parent's frame
Vec3 old_local_pos = body->local_position;
Vec3 old_local_vel = body->local_velocity;
// Convert to global coordinates using old parent
if (old_parent_idx >= 0 && old_parent_idx < sim->body_count) {
CelestialBody* old_parent = &sim->bodies[old_parent_idx];
body->position = vec3_add(old_local_pos, old_parent->position);
body->velocity = vec3_add(old_local_vel, old_parent->velocity);
} else {
// old_parent is root (Sun): local = global
body->position = old_local_pos;
body->velocity = old_local_vel;
}
// Update parent index
body->parent_index = new_parent_idx;
// Convert to local coordinates using new parent
if (new_parent_idx >= 0 && new_parent_idx < sim->body_count) {
CelestialBody* new_parent_body = &sim->bodies[new_parent_idx];
body->local_position = vec3_sub(body->position, new_parent_body->position);
body->local_velocity = vec3_sub(body->velocity, new_parent_body->velocity);
} else {
// new_parent is root (Sun): global = local
body->local_position = body->position;
body->local_velocity = body->velocity;
}
}
Update update_simulation():
int new_parent = find_dominant_body(sim, i);
if (new_parent != body->parent_index) {
transition_to_new_parent(sim, body, body->parent_index, new_parent);
}
Files to modify:
src/simulation.h(add function declaration)src/simulation.cpp(extract and refactor transition logic)
Estimated complexity: Low Risk: Low (pure refactor, no behavior change)
Expected outcome:
- ✅ Cleaner code with better separation of concerns
- ✅ Easier to unit test transition logic
- ✅ Follows Phase 3 plan from hierarchical_frames_plan.md
Tests to add:
- Unit tests for
transition_to_new_parent()with all scenarios:- Body→Body transition
- Body→Root transition
- Root→Body transition
- Root→Root transition (edge case)
Phase 4: Multi-Body Transition Test Configs
Goal: Create realistic test scenarios for patched conics
Test Config 1: Satellite Rendezvous
File: tests/configs/satellite_rendezvous.toml
[[bodies]]
name = "Sun"
mass = 1.989e30
radius = 6.96e8
position = { x = 0.0, y = 0.0, z = 0.0 }
parent_index = -1
color = { r = 1.0, g = 1.0, b = 0.0 }
eccentricity = 0.0
semi_major_axis = 0.0
[[bodies]]
name = "Earth"
mass = 5.972e24
radius = 6.371e6
position = { x = 1.496e11, y = 0.0, z = 0.0 }
parent_index = 0
color = { r = 0.0, g = 0.5, b = 1.0 }
eccentricity = 0.0
semi_major_axis = 1.496e11
[[bodies]]
name = "Moon"
mass = 7.342e22
radius = 1.737e6
position = { x = 1.49984e11, y = 0.0, z = 0.0 }
parent_index = 1
color = { r = 0.7, g = 0.7, b = 0.7 }
eccentricity = 0.0
semi_major_axis = 3.844e8
[[bodies]]
name = "Satellite"
mass = 1.0e3
radius = 1.0e1
position = { x = 1.500e11, y = 1.0e8, z = 0.0 }
parent_index = 1
color = { r = 1.0, g = 0.0, b = 1.0 }
eccentricity = 0.2
semi_major_axis = 4.0e8
Scenario: Satellite launches from Earth, transfers to Moon, returns
Test Config 2: Interplanetary Transfer
File: tests/configs/interplanetary_transfer.toml
[[bodies]]
name = "Sun"
mass = 1.989e30
radius = 6.96e8
position = { x = 0.0, y = 0.0, z = 0.0 }
parent_index = -1
color = { r = 1.0, g = 1.0, b = 0.0 }
eccentricity = 0.0
semi_major_axis = 0.0
[[bodies]]
name = "Earth"
mass = 5.972e24
radius = 6.371e6
position = { x = 1.496e11, y = 0.0, z = 0.0 }
parent_index = 0
color = { r = 0.0, g = 0.5, b = 1.0 }
eccentricity = 0.0
semi_major_axis = 1.496e11
[[bodies]]
name = "Mars"
mass = 6.39e23
radius = 3.3895e6
position = { x = 2.279e11, y = 0.0, z = 0.0 }
parent_index = 0
color = { r = 0.8, g = 0.3, b = 0.1 }
eccentricity = 0.0
semi_major_axis = 2.279e11
[[bodies]]
name = "Probe"
mass = 1.0e3
radius = 1.0e1
position = { x = 1.496e11, y = 0.0, z = 0.0 }
parent_index = 1
color = { r = 0.0, g = 1.0, b = 0.0 }
eccentricity = 0.5
semi_major_axis = 1.888e11
Scenario: Probe: Earth→Sun→Mars
Test Config 3: Moon Capture
File: tests/configs/moon_capture.toml
[[bodies]]
name = "Sun"
mass = 1.989e30
radius = 6.96e8
position = { x = 0.0, y = 0.0, z = 0.0 }
parent_index = -1
color = { r = 1.0, g = 1.0, b = 0.0 }
eccentricity = 0.0
semi_major_axis = 0.0
[[bodies]]
name = "Jupiter"
mass = 1.898e27
radius = 6.9911e7
position = { x = 7.785e11, y = 0.0, z = 0.0 }
parent_index = 0
color = { r = 0.9, g = 0.7, b = 0.5 }
eccentricity = 0.0
semi_major_axis = 7.785e11
[[bodies]]
name = "Ganymede"
mass = 1.48e23
radius = 2.634e6
position = { x = 7.796e11, y = 0.0, z = 0.0 }
parent_index = 1
color = { r = 0.6, g = 0.6, b = 0.5 }
eccentricity = 0.0
semi_major_axis = 1.070e9
[[bodies]]
name = "Comet"
mass = 1.0e14
radius = 5.0e3
position = { x = 1.0e12, y = 0.0, z = 0.0 }
parent_index = 0
color = { r = 0.5, g = 0.8, b = 1.0 }
eccentricity = 2.0
semi_major_axis = -5.0e11
Scenario: Comet: Sun→Jupiter→Ganymede→Sun
Files to create:
tests/configs/satellite_rendezvous.tomltests/configs/interplanetary_transfer.tomltests/configs/moon_capture.toml
Estimated complexity: Low Risk: Low (configs only, no code changes)
Expected outcome:
- ✅ Realistic test scenarios for patched conics
- ✅ Covers multi-leg missions
- ✅ Tests root body transitions
Phase 5: Adaptive Timestepping (Performance + Accuracy)
Goal: Different timesteps for different orbital phases
Problem: Fixed 60s timestep is:
- Too coarse for fast orbital phases (moon capture)
- Too slow for deep-space phases
Proposed solution: Adaptive timestep based on orbital period
Implementation:
double calculate_adaptive_timestep(CelestialBody* body, CelestialBody* parent) {
if (parent == NULL || body->semi_major_axis <= 0.0) {
return 60.0; // Default timestep
}
// Calculate orbital period using Kepler's third law
double T = 2.0 * M_PI * sqrt(pow(body->semi_major_axis, 3) / (G * parent->mass));
// Use 1/1000 of orbital period as timestep
double adaptive_dt = T / 1000.0;
// Clamp to reasonable bounds
adaptive_dt = fmax(adaptive_dt, 10.0); // Minimum 10s
adaptive_dt = fmin(adaptive_dt, 600.0); // Maximum 600s
return adaptive_dt;
}
Changes required:
- Add per-body timesteps to
SimulationState - Update
update_simulation()to use adaptive timesteps - Add synchronization mechanism (multiple timesteps)
Complexity: High Risk: High (may affect energy conservation, requires careful testing)
Expected outcome:
- ✅ Better accuracy for fast orbits (moon capture)
- ✅ Faster simulation for deep-space phases
- ✅ Energy conserved across transitions
Tests to add:
- Verify energy drift with adaptive timesteps
- Verify orbital period accuracy with adaptive timesteps
- Test stability across SOI transitions
Phase 6: Enhanced Debugging & Visualization
Goal: Better tools for debugging transitions and planning missions
Features to add:
6.1 Transition Logging
void log_transition(SimulationState* sim, CelestialBody* body,
int old_parent, int new_parent, double time) {
printf("[%.2f days] %s: parent %d → %d\n",
time / 86400.0, body->name, old_parent, new_parent);
}
6.2 Visual SOI Boundaries
- Add wireframe spheres in renderer for SOI radii
- Toggle with keyboard key
- Use different colors for different bodies
6.3 Trajectory Prediction
- Predict future trajectory for given time span
- Show predicted SOI crossings
- Assist with mission planning
Files to modify:
src/simulation.cpp(logging)src/renderer.cpp(SOI visualization)src/renderer.cpp(trajectory prediction - new function)
Estimated complexity: Medium Risk: Low (mostly UI/visualization)
Expected outcome:
- ✅ Better debugging tools for transitions
- ✅ Visual SOI boundaries
- ✅ Trajectory prediction for mission planning
Open Questions
Question 1: Hysteresis Strategy
Context: The 0.5x hysteresis factor prevents oscillation but creates one-way barriers.
Options:
-
Option A: Remove hysteresis entirely
- Pros: Simple, enables all transitions
- Cons: May cause oscillation at SOI boundaries
-
Option B: Adaptive hysteresis (recommended)
- Pros: Prevents oscillation while enabling round-trips
- Cons: More complex logic
Recommendation: Option B - keep stability while enabling your use case
Decision needed: Which approach should we implement?
Question 2: Integration Timing
Context: Transition happens before integration in the same timestep, using coordinates from the end of the previous timestep. This may cause velocity discontinuities.
Options:
-
Option A: Keep current approach (transition at start of timestep)
- Pros: Simple, works for most cases
- Cons: May have slight velocity discontinuities
-
Option B: Transition in middle of timestep (half-step approach)
- Pros: Better accuracy, smoother transitions
- Cons: More complex, requires half-step RK4
-
Option C: Transition after integration, then integrate again
- Pros: Ensures continuity
- Cons: Doubles computation, complex
Recommendation: Start with Option A, move to Option B if needed
Decision needed: Should we implement half-step transitions for better accuracy?
Question 3: Test Priorities
Context: We have three test scenarios to validate patched conics.
Options:
-
Option A: Start with "Satellite Rendezvous" (Earth→Moon→Earth)
- Pros: Simpler, 2-body transition, quick to implement
- Cons: Doesn't test root body transitions
-
Option B: Start with "Interplanetary Transfer" (Earth→Sun→Mars)
- Pros: Tests root body transitions, realistic scenario
- Cons: More complex, longer simulation time
-
Option C: Implement all three in parallel
- Pros: Comprehensive coverage
- Cons: More work upfront
Recommendation: Option A first, then Option B, then Option C
Decision needed: Which test scenario should we start with?
Question 4: Adaptive Timesteps Priority
Context: Fixed 60s timestep may be suboptimal for different orbital phases.
Options:
-
Option A: Implement adaptive timesteps now
- Pros: Better accuracy and performance from the start
- Cons: Adds complexity, may delay other features
-
Option B: Use fixed timesteps for now, optimize later
- Pros: Simpler implementation, focus on transitions first
- Cons: May need to revisit later
Recommendation: Option B - get transitions working first, then optimize
Decision needed: Is adaptive timestepping critical for your use case?
Success Criteria
Phase 1 Success
- Satellite can transition to/from Sun
- Tests pass for Earth→Sun→Mars
- No energy conservation issues during transitions
Phase 2 Success
- Round-trip transitions work (Earth→Sun→Earth)
- No oscillation at SOI boundaries
- Tests validate stability
Phase 3 Success
- Transition logic extracted to separate function
- Unit tests cover all transition scenarios
- Code is cleaner and more maintainable
Phase 4 Success
- Three test configs created
- All test scenarios pass
- Configs are realistic and well-documented
Phase 5 Success (Optional)
- Adaptive timesteps implemented
- Energy drift < 1% with adaptive timesteps
- Performance improved for deep-space phases
Phase 6 Success (Optional)
- Transition logging works
- SOI boundaries visualized
- Trajectory prediction functional
Timeline Estimate
- Phase 1: 1-2 days (fix root body transitions)
- Phase 2: 1 day (adaptive hysteresis)
- Phase 3: 0.5 days (refactoring)
- Phase 4: 1 day (test configs)
- Phase 5: 2-3 days (adaptive timesteps - optional)
- Phase 6: 1-2 days (debugging/visualization - optional)
Total for core features (Phases 1-4): 4-5 days Total with all features: 8-10 days
Dependencies
Required
- ✅ Hierarchical coordinate frames (complete)
- ✅ Local frame integration (complete)
- ✅ Parent-first update order (complete)
- ⏸️ Root body transitions (Phase 1 - pending)
Optional
- ⏸️ Adaptive timesteps (Phase 5 - optional)
- ⏸️ Debugging tools (Phase 6 - optional)
Risks and Mitigations
High Risk
-
Energy conservation during transitions
- Mitigation: Careful testing, energy drift checks
- Backup: Keep old code for rollback
-
Numerical precision across scales
- Mitigation: Use local frames (already implemented)
- Backup: Double precision (already used)
Medium Risk
-
Oscillation at SOI boundaries
- Mitigation: Adaptive hysteresis (Phase 2)
- Backup: Increase hysteresis factor if needed
-
Timestep too coarse for fast orbits
- Mitigation: Adaptive timesteps (Phase 5 - optional)
- Backup: Reduce fixed timestep if needed
Low Risk
- Code complexity increases
- Mitigation: Good unit tests, refactoring (Phase 3)
- Backup: Keep functions small and focused
References
docs/hierarchical_frames_plan.md- Phase 3: SOI Transition with Frame Transformdocs/implementation_plan.md- SOI Transition Mechanics sectionsrc/simulation.cpp- Current SOI transition implementation (lines 103-121)tests/test_soi_transition.cpp- Current SOI transition tests