From b4f820c4497fe19b33f34f7403e7fb3d857c660b Mon Sep 17 00:00:00 2001 From: cinnaboot Date: Thu, 16 Apr 2026 22:13:11 +0000 Subject: [PATCH] add initial hohmann rendezvous implementation --- src/rendezvous_hohmann.cpp | 132 +++++++++++++ src/rendezvous_hohmann.h | 76 ++++++++ tests/test_rendezvous_hohmann.cpp | 292 +++++++++++++++++++++++++++++ tests/test_rendezvous_hohmann.toml | 70 +++++++ 4 files changed, 570 insertions(+) create mode 100644 src/rendezvous_hohmann.cpp create mode 100644 src/rendezvous_hohmann.h create mode 100644 tests/test_rendezvous_hohmann.cpp create mode 100644 tests/test_rendezvous_hohmann.toml diff --git a/src/rendezvous_hohmann.cpp b/src/rendezvous_hohmann.cpp new file mode 100644 index 0000000..f8a6fea --- /dev/null +++ b/src/rendezvous_hohmann.cpp @@ -0,0 +1,132 @@ +#include "rendezvous_hohmann.h" +#include +#include + +// Mean motion: n = sqrt(mu / a^3) +static double calc_mean_motion(double radius, double mass) { + double mu = G * mass; + return sqrt(mu / pow(radius, 3)); +} + +// Hohmann transfer time (half orbit of transfer ellipse) +static double hohmann_transfer_time(double r1, double r2, double mass) { + double mu = G * mass; + double a_transfer = (r1 + r2) / 2.0; + double T_transfer = 2.0 * M_PI * sqrt(pow(a_transfer, 3) / mu); + return T_transfer / 2.0; +} + +// Calculate required angular separation at first burn +// For Hohmann transfer: target should be at specific angle when chaser burns +// Returns: required angular separation in radians [0, 2π) +static double required_separation(double r1, double r2, double mass) { + double transfer_time = hohmann_transfer_time(r1, r2, mass); + double n2 = calc_mean_motion(r2, mass); + double target_angle = n2 * transfer_time; + return M_PI - target_angle; +} + + + +// Normalize angle to [0, 2π) +static double normalize_angle_2pi(double angle) { + while (angle < 0.0) { + angle += 2.0 * M_PI; + } + while (angle >= 2.0 * M_PI) { + angle -= 2.0 * M_PI; + } + return angle; +} + +// Normalize angle to [-π, π] for shortest path +static double normalize_angle_pi(double angle) { + angle = normalize_angle_2pi(angle); + while (angle > M_PI) { + angle -= 2.0 * M_PI; + } + while (angle < -M_PI) { + angle += 2.0 * M_PI; + } + return angle; +} + +// Calculate wait time before starting Hohmann transfer +// Determines how long to wait before executing the first burn so that +// both chaser and target arrive at the interception point simultaneously. +// Returns: wait time in seconds. Positive = wait, negative = transfer already late +double calculate_wait_time_for_hohmann( + double initial_orbit_radius, + double target_orbit_radius, + double angular_separation, + double central_mass +) { + double required_sep = required_separation(initial_orbit_radius, target_orbit_radius, central_mass); + double n1 = calc_mean_motion(initial_orbit_radius, central_mass); + double n2 = calc_mean_motion(target_orbit_radius, central_mass); + double rel_angular_vel = n1 - n2; + + // Normalize required separation to [0, 2*pi) + required_sep = normalize_angle_2pi(required_sep); + + // Normalize current separation to [0, 2*pi) + double sep = normalize_angle_2pi(angular_separation); + + // How much more angle needs to close + double angle_to_close = required_sep - sep; + + // Normalize to [-pi, pi] for shortest path + angle_to_close = normalize_angle_pi(angle_to_close); + + // Wait time = angle_to_close / relative_angular_velocity + return angle_to_close / rel_angular_vel; +} + +// Calculate required angular separation for Hohmann transfer +// Computes the ideal angle between chaser and target at the moment +// of first burn to ensure simultaneous arrival at target orbit. +// Returns: required angular separation in radians (-2π, 2π) +double calculate_required_separation_for_hohmann( + double initial_orbit_radius, + double target_orbit_radius, + double central_mass +) { + double required_sep = required_separation(initial_orbit_radius, target_orbit_radius, central_mass); + return normalize_angle_pi(required_sep); +} + +// Verify spacecraft is on correct Hohmann transfer orbit +// Checks if current orbit matches expected Hohmann transfer parameters. +// Returns: true if orbit is on Hohmann transfer, false otherwise +bool verify_hohmann_transfer_orbit( + const OrbitalElements* orbit, + double r1, + double r2, + double tolerance +) { + double expected_a = (r1 + r2) / 2.0; + double actual_a = orbit->semi_major_axis; + double diff = fabs(actual_a - expected_a); + return diff < tolerance; +} + +// Check if Hohmann transfer is complete +// Determines if transfer time has elapsed and spacecraft is at target radius. +// Returns: true if transfer is complete, false otherwise +bool hohmann_transfer_complete( + double transfer_start_time, + double current_time, + double transfer_time, + double target_radius, + double current_radius, + double tolerance +) { + // Check if enough time has elapsed + if (current_time < transfer_start_time + transfer_time - tolerance) { + return false; + } + + // Check if at target radius + double radius_diff = fabs(current_radius - target_radius); + return radius_diff < tolerance; +} diff --git a/src/rendezvous_hohmann.h b/src/rendezvous_hohmann.h new file mode 100644 index 0000000..a48f594 --- /dev/null +++ b/src/rendezvous_hohmann.h @@ -0,0 +1,76 @@ +#ifndef RENDEZVOUS_HOHMANN_H +#define RENDEZVOUS_HOHMANN_H + +#include "orbital_mechanics.h" +#include "simulation.h" +#include "maneuver.h" + +// Hohmann transfer calculations for coplanar circular orbits +// Provides functions to plan and execute Hohmann transfers between +// circular coplanar orbits, including phasing calculations to +// intercept moving targets. + +// Calculate wait time before starting Hohmann transfer +// Determines how long to wait before executing the first burn so that +// both chaser and target arrive at the interception point simultaneously. +// Returns: wait time in seconds. Positive = wait, negative = transfer already late +double calculate_wait_time_for_hohmann( + double initial_orbit_radius, // Current circular orbit radius (meters) + double target_orbit_radius, // Target circular orbit radius (meters) + double angular_separation, // Current angle from chaser to target (radians, [0, 2π)) + double central_mass // Mass of central body (kg) +); + +// Calculate required angular separation for Hohmann transfer +// Computes the ideal angle between chaser and target at the moment +// of first burn to ensure simultaneous arrival at target orbit. +// Returns: required angular separation in radians (-2π, 2π) +double calculate_required_separation_for_hohmann( + double initial_orbit_radius, // Current circular orbit radius (meters) + double target_orbit_radius, // Target circular orbit radius (meters) + double central_mass // Mass of central body (kg) +); + +// Create first burn maneuver for Hohmann transfer +// Adds a prograde burn maneuver to enter the Hohmann transfer orbit. +// Returns: maneuver index on success, -1 on failure +int create_hohmann_departure_maneuver( + SimulationState* sim, // Simulation state to add maneuver to + int chaser_index, // Index of chaser spacecraft + double delta_v, // Delta-v magnitude from calculate_hohmann_transfer (m/s) + double trigger_true_anomaly // True anomaly for trigger (-1 for immediate) +); + +// Create second burn maneuver for Hohmann transfer +// Adds a circularization burn maneuver to match target orbit. +// Returns: maneuver index on success, -1 on failure +int create_hohmann_arrival_maneuver( + SimulationState* sim, // Simulation state to add maneuver to + int chaser_index, // Index of chaser spacecraft + double delta_v, // Delta-v magnitude from calculate_hohmann_transfer (m/s) + double trigger_true_anomaly // True anomaly at target orbit for trigger +); + +// Verify spacecraft is on correct Hohmann transfer orbit +// Checks if current orbit matches expected Hohmann transfer parameters. +// Returns: true if orbit is on Hohmann transfer, false otherwise +bool verify_hohmann_transfer_orbit( + const OrbitalElements* orbit, // Current orbit after first burn + double r1, // Initial orbit radius (meters) + double r2, // Target orbit radius (meters) + double tolerance // Acceptable deviation in semi-major axis (meters) +); + +// Check if Hohmann transfer is complete +// Determines if transfer time has elapsed and spacecraft is at target radius. +// Returns: true if transfer is complete, false otherwise +bool hohmann_transfer_complete( + double transfer_start_time, // When first burn was executed + double current_time, // Current simulation time + double transfer_time, // Expected transfer duration (from calculate_hohmann_transfer) + double target_radius, // Target orbit radius (meters) + double current_radius, // Current orbital radius (meters) + double tolerance // Radius tolerance for completion check (meters) +); + +#endif // RENDEZVOUS_HOHMANN_H diff --git a/tests/test_rendezvous_hohmann.cpp b/tests/test_rendezvous_hohmann.cpp new file mode 100644 index 0000000..f29b7df --- /dev/null +++ b/tests/test_rendezvous_hohmann.cpp @@ -0,0 +1,292 @@ +#include +#include +#include "../src/physics.h" +#include "../src/orbital_mechanics.h" +#include "../src/simulation.h" +#include "../src/orbital_objects.h" +#include "../src/rendezvous_hohmann.h" +#include "../src/config_loader.h" +#include +#include + + +using Catch::Matchers::WithinAbs; + + +// ============================================================================ +// Helper Functions +// ============================================================================ + +// TODO: Add find_spacecraft_by_name() to simulation.h interface +// FIXME: This helper should be part of the public simulation API for testability +static int find_spacecraft_by_name(SimulationState* sim, const char* name) { + for (int i = 0; i < sim->craft_count; i++) { + if (strcmp(sim->spacecraft[i].name, name) == 0) { + return i; + } + } + return -1; +} + + +TEST_CASE("Config loading for Hohmann transfer", "[rendezvous_hohmann][config]") { + const double TIME_STEP = 30.0; + + SimulationState* sim = create_simulation(3, 5, 10, TIME_STEP); + + REQUIRE(load_system_config(sim, "tests/test_rendezvous_hohmann.toml")); + + REQUIRE(sim->body_count == 1); + REQUIRE(std::string(sim->bodies[0].name) == "Earth"); + + REQUIRE(sim->craft_count == 3); + REQUIRE(std::string(sim->spacecraft[0].name) == "Target_Satellite"); + REQUIRE(std::string(sim->spacecraft[1].name) == "Chaser_Lower"); + REQUIRE(std::string(sim->spacecraft[2].name) == "Chaser_Higher"); + + REQUIRE(sim->spacecraft[0].parent_index == 0); + REQUIRE(sim->spacecraft[1].parent_index == 0); + REQUIRE(sim->spacecraft[2].parent_index == 0); + + // Verify initial orbits + REQUIRE_THAT(sim->spacecraft[0].orbit.semi_major_axis, + WithinAbs(6.771e6, 1.0)); // 400 km altitude + REQUIRE_THAT(sim->spacecraft[1].orbit.semi_major_axis, + WithinAbs(6.671e6, 1.0)); // 300 km altitude + REQUIRE_THAT(sim->spacecraft[2].orbit.semi_major_axis, + WithinAbs(6.871e6, 1.0)); // 500 km altitude + + // Verify initial true anomalies + REQUIRE_THAT(sim->spacecraft[0].orbit.true_anomaly, + WithinAbs(0.0, 0.001)); + REQUIRE_THAT(sim->spacecraft[1].orbit.true_anomaly, + WithinAbs(4.71238898038469, 0.001)); // 270° = 3π/2 + REQUIRE_THAT(sim->spacecraft[2].orbit.true_anomaly, + WithinAbs(1.5707963267948966, 0.001)); // 90° = π/2 + + destroy_simulation(sim); +} + +TEST_CASE("Calculate wait time for Hohmann transfer (lower to higher)", "[rendezvous_hohmann][phasing]") { + const double TIME_STEP = 30.0; + + SimulationState* sim = create_simulation(3, 5, 10, TIME_STEP); + + REQUIRE(load_system_config(sim, "tests/test_rendezvous_hohmann.toml")); + + int target_idx = find_spacecraft_by_name(sim, "Target_Satellite"); + int chaser_lower_idx = find_spacecraft_by_name(sim, "Chaser_Lower"); + + REQUIRE(target_idx >= 0); + REQUIRE(chaser_lower_idx >= 0); + + Spacecraft* target = &sim->spacecraft[target_idx]; + Spacecraft* chaser = &sim->spacecraft[chaser_lower_idx]; + CelestialBody* earth = &sim->bodies[0]; + + initialize_orbital_objects(sim); + + double r1 = vec3_magnitude(chaser->local_position); + double r2 = vec3_magnitude(target->local_position); + + SECTION("Zero angular separation - immediate transfer not possible") { + // If chaser is directly behind target, need to wait for target to move ahead + double angular_separation = 0.0; + double wait_time = calculate_wait_time_for_hohmann(r1, r2, angular_separation, earth->mass); + + INFO("Wait time: " << wait_time << " s"); + // Since lower orbit is faster, chaser will catch up, so wait time should be positive + REQUIRE_THAT(wait_time, WithinAbs(1358.16, 1.0)); + } + + SECTION("Small angular separation") { + double angular_separation = 0.5; // ~29 degrees + double wait_time = calculate_wait_time_for_hohmann(r1, r2, angular_separation, earth->mass); + + INFO("Angular separation: " << angular_separation << " rad"); + INFO("Wait time: " << wait_time << " s"); + + REQUIRE_THAT(wait_time, WithinAbs(-18192.7, 1.0)); + } + + SECTION("Large angular separation (near 2π)") { + double angular_separation = 6.0; // ~344 degrees + double wait_time = calculate_wait_time_for_hohmann(r1, r2, angular_separation, earth->mass); + + INFO("Angular separation: " << angular_separation << " rad"); + INFO("Wait time: " << wait_time << " s"); + + REQUIRE_THAT(wait_time, WithinAbs(12431.2, 1.0)); + } + + destroy_simulation(sim); +} + +TEST_CASE("Calculate wait time for Hohmann transfer (higher to lower)", "[rendezvous_hohmann][phasing]") { + const double TIME_STEP = 30.0; + + SimulationState* sim = create_simulation(3, 5, 10, TIME_STEP); + + REQUIRE(load_system_config(sim, "tests/test_rendezvous_hohmann.toml")); + + int target_idx = find_spacecraft_by_name(sim, "Target_Satellite"); + int chaser_higher_idx = find_spacecraft_by_name(sim, "Chaser_Higher"); + + REQUIRE(target_idx >= 0); + REQUIRE(chaser_higher_idx >= 0); + + Spacecraft* target = &sim->spacecraft[target_idx]; + Spacecraft* chaser = &sim->spacecraft[chaser_higher_idx]; + CelestialBody* earth = &sim->bodies[0]; + + initialize_orbital_objects(sim); + + double r1 = vec3_magnitude(chaser->local_position); + double r2 = vec3_magnitude(target->local_position); + + SECTION("Zero angular separation - target must catch up") { + // Higher orbit is slower, so target must catch up + double angular_separation = 0.0; + double wait_time = calculate_wait_time_for_hohmann(r1, r2, angular_separation, earth->mass); + + INFO("Wait time: " << wait_time << " s"); + REQUIRE_THAT(wait_time, WithinAbs(1414.46, 1.0)); + } + + SECTION("Small angular separation") { + double angular_separation = 0.3; // ~17 degrees + double wait_time = calculate_wait_time_for_hohmann(r1, r2, angular_separation, earth->mass); + + INFO("Angular separation: " << angular_separation << " rad"); + INFO("Wait time: " << wait_time << " s"); + + REQUIRE_THAT(wait_time, WithinAbs(13586.2, 1.0)); + } + + destroy_simulation(sim); +} + +TEST_CASE("Calculate required separation for Hohmann transfer", "[rendezvous_hohmann][phasing]") { + const double TIME_STEP = 30.0; + + SimulationState* sim = create_simulation(3, 5, 10, TIME_STEP); + + REQUIRE(load_system_config(sim, "tests/test_rendezvous_hohmann.toml")); + + int target_idx = find_spacecraft_by_name(sim, "Target_Satellite"); + int chaser_lower_idx = find_spacecraft_by_name(sim, "Chaser_Lower"); + int chaser_higher_idx = find_spacecraft_by_name(sim, "Chaser_Higher"); + + REQUIRE(target_idx >= 0); + REQUIRE(chaser_lower_idx >= 0); + REQUIRE(chaser_higher_idx >= 0); + + Spacecraft* target = &sim->spacecraft[target_idx]; + Spacecraft* chaser_lower = &sim->spacecraft[chaser_lower_idx]; + Spacecraft* chaser_higher = &sim->spacecraft[chaser_higher_idx]; + CelestialBody* earth = &sim->bodies[0]; + + initialize_orbital_objects(sim); + + double r_lower = vec3_magnitude(chaser_lower->local_position); + double r_target = vec3_magnitude(target->local_position); + double r_higher = vec3_magnitude(chaser_higher->local_position); + + SECTION("Lower to higher transfer") { + double required_separation = calculate_required_separation_for_hohmann(r_lower, r_target, earth->mass); + + INFO("Required separation: " << required_separation << " rad"); + INFO("Required separation (deg): " << required_separation * 180.0 / M_PI << "°"); + + REQUIRE_THAT(required_separation, WithinAbs(0.034734, 0.001)); + } + + SECTION("Higher to lower transfer") { + double required_separation = calculate_required_separation_for_hohmann(r_higher, r_target, earth->mass); + + INFO("Required separation: " << required_separation << " rad"); + INFO("Required separation (deg): " << required_separation * 180.0 / M_PI << "°"); + + REQUIRE_THAT(required_separation, WithinAbs(-0.0348625, 0.001)); + } + + SECTION("Equal radii - no transfer needed") { + double required_separation = calculate_required_separation_for_hohmann(r_target, r_target, earth->mass); + + INFO("Required separation: " << required_separation << " rad"); + + REQUIRE_THAT(required_separation, WithinAbs(0.0, 0.001)); + } + + destroy_simulation(sim); +} + +SCENARIO("Complete Hohmann transfer phasing workflow", "[rendezvous_hohmann][workflow]") { + const double TIME_STEP = 10.0; + + SimulationState* sim = create_simulation(3, 5, 10, TIME_STEP); + + REQUIRE(load_system_config(sim, "tests/test_rendezvous_hohmann.toml")); + + int target_idx = find_spacecraft_by_name(sim, "Target_Satellite"); + int chaser_lower_idx = find_spacecraft_by_name(sim, "Chaser_Lower"); + + REQUIRE(target_idx >= 0); + REQUIRE(chaser_lower_idx >= 0); + + Spacecraft* target = &sim->spacecraft[target_idx]; + Spacecraft* chaser = &sim->spacecraft[chaser_lower_idx]; + CelestialBody* earth = &sim->bodies[0]; + + initialize_orbital_objects(sim); + + SECTION("Calculate and verify phasing for lower-to-higher transfer") { + double r1 = vec3_magnitude(chaser->local_position); + double r2 = vec3_magnitude(target->local_position); + + INFO("Chaser orbit radius: " << r1 << " m"); + INFO("Target orbit radius: " << r2 << " m"); + + // Calculate current angular separation + double current_sep = angular_distance(chaser->orbit.true_anomaly, target->orbit.true_anomaly); + + INFO("Current angular separation: " << current_sep << " rad"); + INFO("Current angular separation (deg): " << current_sep * 180.0 / M_PI << "°"); + + // Calculate required separation for Hohmann + double required_sep = calculate_required_separation_for_hohmann(r1, r2, earth->mass); + + INFO("Required separation: " << required_sep << " rad"); + INFO("Required separation (deg): " << required_sep * 180.0 / M_PI << "°"); + + // Calculate wait time + double wait_time = calculate_wait_time_for_hohmann(r1, r2, current_sep, earth->mass); + + INFO("Wait time: " << wait_time << " s"); + + REQUIRE_THAT(wait_time, WithinAbs(-60062.651728, 0.1)); + + // Verify orbits are circular + REQUIRE_THAT(chaser->orbit.eccentricity, WithinAbs(0.0, 0.001)); + REQUIRE_THAT(target->orbit.eccentricity, WithinAbs(0.0, 0.001)); + } + + SECTION("Calculate Hohmann transfer parameters") { + double r1 = vec3_magnitude(chaser->local_position); + double r2 = vec3_magnitude(target->local_position); + + // Use existing calculate_hohmann_transfer from maneuver.h + HohmannTransfer hohmann = calculate_hohmann_transfer(r1, r2, earth->mass); + + INFO("First burn delta-v: " << hohmann.dv1 << " m/s"); + INFO("Second burn delta-v: " << hohmann.dv2 << " m/s"); + INFO("Transfer time: " << hohmann.transfer_time << " s"); + INFO("Target true anomaly: " << hohmann.true_anomaly_2 << " rad"); + + REQUIRE_THAT(hohmann.dv1, WithinAbs(28.699077, 0.01)); + REQUIRE_THAT(hohmann.dv2, WithinAbs(28.592521, 0.01)); + REQUIRE_THAT(hohmann.transfer_time, WithinAbs(2741.813778, 0.1)); + } + + destroy_simulation(sim); +} diff --git a/tests/test_rendezvous_hohmann.toml b/tests/test_rendezvous_hohmann.toml new file mode 100644 index 0000000..e07c651 --- /dev/null +++ b/tests/test_rendezvous_hohmann.toml @@ -0,0 +1,70 @@ +# Test Configuration: Hohmann Transfer Between Coplanar Circular Orbits +# Three spacecraft in circular, coplanar LEO orbits around Earth +# Tests Hohmann transfer phasing calculations in both directions +# +# Configuration: +# - Target_Satellite: middle orbit (reference) +# - Chaser_Lower: lower orbit, will transfer UP to target +# - Chaser_Higher: higher orbit, will transfer DOWN to target + +[[bodies]] +name = "Earth" +mass = 5.972e24 +radius = 6.371e6 +parent_index = -1 +color = { r = 0.0, g = 0.5, b = 1.0 } +orbit = { + semi_major_axis = 0.0, + eccentricity = 0.0, + true_anomaly = 0.0 +} + +# ========== TARGET SPACECRAFT ========== +# Circular LEO orbit at 400 km altitude +# Reference orbit for Hohmann transfers +[[spacecraft]] +name = "Target_Satellite" +mass = 500.0 +parent_index = 0 +orbit = { + semi_major_axis = 6.771e6, + eccentricity = 0.0, + true_anomaly = 0.0, + inclination = 0.0, + longitude_of_ascending_node = 0.0, + argument_of_periapsis = 0.0 +} + +# ========== CHASER LOWER ========== +# Circular LEO orbit at 300 km altitude (lower than target) +# Will perform Hohmann transfer UP to target orbit +# Starts 90 degrees behind target to test phasing +[[spacecraft]] +name = "Chaser_Lower" +mass = 500.0 +parent_index = 0 +orbit = { + semi_major_axis = 6.671e6, + eccentricity = 0.0, + true_anomaly = 4.71238898038469, + inclination = 0.0, + longitude_of_ascending_node = 0.0, + argument_of_periapsis = 0.0 +} + +# ========== CHASER HIGHER ========== +# Circular LEO orbit at 500 km altitude (higher than target) +# Will perform Hohmann transfer DOWN to target orbit +# Starts 90 degrees ahead of target to test phasing +[[spacecraft]] +name = "Chaser_Higher" +mass = 500.0 +parent_index = 0 +orbit = { + semi_major_axis = 6.871e6, + eccentricity = 0.0, + true_anomaly = 1.5707963267948966, + inclination = 0.0, + longitude_of_ascending_node = 0.0, + argument_of_periapsis = 0.0 +}