@ -3,16 +3,18 @@
# include <cmath>
# include <cstdio>
# define MIN_DISTANCE_CLAMP 1.0
static const double SECONDS_PER_DAY = 86400.0 ;
double calculate_kinetic_energy ( CelestialBody * body ) {
double v_squared = body - > global_velocity . x * body - > global_velocity . x +
body - > global_velocity . y * body - > global_velocity . y +
body - > global_velocity . z * body - > global_velocity . z ;
return 0.5 * body - > mass * v_squared ;
Vec3 v = body - > global_velocity ;
return 0.5 * body - > mass * vec3_dot ( v , v ) ;
}
double calculate_potential_energy_pair ( CelestialBody * body1 , CelestialBody * body2 ) {
double distance = vec3_distance ( body1 - > global_position , body2 - > global_position ) ;
if ( distance < 1.0 ) distance = 1.0 ;
if ( distance < MIN_DISTANCE_CLAMP ) distance = MIN_DISTANCE_CLAMP ;
return - G * body1 - > mass * body2 - > mass / distance ;
}
@ -46,28 +48,13 @@ OrbitalMetrics calculate_orbital_metrics(CelestialBody* body, CelestialBody* par
return metrics ;
}
OrbitTracker * create_orbit_tracker ( int body_index ) {
OrbitTracker * tracker = ( OrbitTracker * ) malloc ( sizeof ( OrbitTracker ) ) ;
tracker - > body_index = body_index ;
tracker - > initial_angle = 0.0 ;
tracker - > previous_angle = 0.0 ;
tracker - > quadrant_transitions = 0 ;
tracker - > orbit_completed = false ;
tracker - > time_at_completion = 0.0 ;
tracker - > min_time_days = 100.0 ;
tracker - > inclination = 0.0 ;
tracker - > longitude_of_ascending_node = 0.0 ;
tracker - > argument_of_periapsis = 0.0 ;
tracker - > has_orbital_elements = false ;
return tracker ;
}
OrbitTracker * create_orbit_tracker_with_min_time ( int body_index , double min_time_days ) {
OrbitTracker * tracker = ( OrbitTracker * ) malloc ( sizeof ( OrbitTracker ) ) ;
tracker - > body_index = body_index ;
tracker - > initial_angle = 0.0 ;
tracker - > previous_angle = 0.0 ;
tracker - > quadrant_transitions = 0 ;
tracker - > accumulated_rotation = 0.0 ;
tracker - > wrap_count = 0 ;
tracker - > orbit_completed = false ;
tracker - > time_at_completion = 0.0 ;
tracker - > min_time_days = min_time_days ;
@ -78,6 +65,10 @@ OrbitTracker* create_orbit_tracker_with_min_time(int body_index, double min_time
return tracker ;
}
OrbitTracker * create_orbit_tracker ( int body_index ) {
return create_orbit_tracker_with_min_time ( body_index , 100.0 ) ;
}
OrbitTracker * create_orbit_tracker_3d ( int body_index , double min_time_days ,
double inclination , double lon_ascending_node ,
double argument_of_periapsis ) {
@ -89,14 +80,19 @@ OrbitTracker* create_orbit_tracker_3d(int body_index, double min_time_days,
return tracker ;
}
void reset_orbit_ tracker ( OrbitTracker * tracker ) {
static void reset_tracker_fields ( OrbitTracker * tracker ) {
tracker - > initial_angle = 0.0 ;
tracker - > previous_angle = 0.0 ;
tracker - > quadrant_transitions = 0 ;
tracker - > accumulated_rotation = 0.0 ;
tracker - > wrap_count = 0 ;
tracker - > orbit_completed = false ;
tracker - > time_at_completion = 0.0 ;
}
void reset_orbit_tracker ( OrbitTracker * tracker ) {
reset_tracker_fields ( tracker ) ;
}
void update_orbit_tracker ( OrbitTracker * tracker , CelestialBody * body , CelestialBody * parent , double current_time ) {
if ( tracker - > orbit_completed ) return ;
@ -117,10 +113,11 @@ void update_orbit_tracker(OrbitTracker* tracker, CelestialBody* body, CelestialB
current_angle = atan2 ( relative_pos . y , relative_pos . x ) ;
}
if ( tracker - > quadrant_transitions = = 0 ) {
if ( tracker - > wrap_count = = 0 ) {
tracker - > initial_angle = current_angle ;
tracker - > previous_angle = current_angle ;
tracker - > quadrant_transitions = 1 ;
tracker - > accumulated_rotation = 0.0 ;
tracker - > wrap_count = 1 ;
return ;
}
@ -128,23 +125,20 @@ void update_orbit_tracker(OrbitTracker* tracker, CelestialBody* body, CelestialB
if ( angle_diff > M_PI ) {
angle_diff - = 2.0 * M_PI ;
tracker - > quadrant_transitions + + ;
tracker - > wrap_count + + ;
}
if ( angle_diff < - M_PI ) {
angle_diff + = 2.0 * M_PI ;
tracker - > quadrant_transitions + + ;
tracker - > wrap_count + + ;
}
double total_rotation = current_angle - tracker - > initial_angle ;
if ( total_rotation < - M_PI ) total_rotation + = 2.0 * M_PI ;
if ( total_rotation > M_PI ) total_rotation - = 2.0 * M_PI ;
tracker - > accumulated_rotation + = angle_diff ;
const double SECONDS_PER_DAY = 86400.0 ;
double min_time_seconds = tracker - > min_time_days * SECONDS_PER_DAY ;
if ( tracker - > quadrant_transitions > = 2 & &
if ( tracker - > wrap_count > = 2 & &
current_time > min_time_seconds & &
fabs ( total_rotation ) < 0.05 ) {
tracker - > accumulated_rotation > = 2.0 * M_PI ) {
tracker - > orbit_completed = true ;
tracker - > time_at_completion = current_time ;
}
@ -183,35 +177,35 @@ int dump_simulation_state(SimulationState* sim, const char* label,
offset + = snprintf ( buffer + offset , buffer_size - offset ,
" Spacecraft (%d): \n " , sim - > craft_count ) ;
for ( int i = 0 ; i < sim - > craft_count ; i + + ) {
Spacecraft * s = & sim - > spacecraft [ i ] ;
double r = sqrt ( s - > local_position . x * s - > local_position . x +
s - > local_position . y * s - > local_position . y +
s - > local_position . z * s - > local_position . z ) ;
double v = sqrt ( s - > local_velocity . x * s - > local_velocity . x +
s - > local_velocity . y * s - > local_velocity . y +
s - > local_velocity . z * s - > local_velocity . z ) ;
offset + = snprintf ( buffer + offset , buffer_size - offset ,
" [%d] %s: r=%.1f v=%.1f nu=%.5f a=%.1f e=%.6f, omega=%.6f \n " ,
i , s - > name , r , v ,
s - > orbit . true_anomaly ,
s - > orbit . semi_major_axis ,
s - > orbit . eccentricity ,
s - > orbit . argument_of_periapsis ) ;
offset + = snprintf ( buffer + offset , buffer_size - offset ,
" pos=(%.1f, %.1f, %.1f) vel=(%.1f, %.1f, %.1f) \n " ,
s - > local_position . x , s - > local_position . y , s - > local_position . z ,
s - > local_velocity . x , s - > local_velocity . y , s - > local_velocity . z ) ;
if ( sim - > spacecraft ) {
for ( int i = 0 ; i < sim - > craft_count ; i + + ) {
Spacecraft * s = & sim - > spacecraft [ i ] ;
double r = vec3_magnitude ( s - > local_position ) ;
double v = vec3_magnitude ( s - > local_velocity ) ;
offset + = snprintf ( buffer + offset , buffer_size - offset ,
" [%d] %s: r=%.1f v=%.1f nu=%.5f a=%.1f e=%.6f, omega=%.6f \n " ,
i , s - > name , r , v ,
s - > orbit . true_anomaly ,
s - > orbit . semi_major_axis ,
s - > orbit . eccentricity ,
s - > orbit . argument_of_periapsis ) ;
offset + = snprintf ( buffer + offset , buffer_size - offset ,
" pos=(%.1f, %.1f, %.1f) vel=(%.1f, %.1f, %.1f) \n " ,
s - > local_position . x , s - > local_position . y , s - > local_position . z ,
s - > local_velocity . x , s - > local_velocity . y , s - > local_velocity . z ) ;
}
}
offset + = snprintf ( buffer + offset , buffer_size - offset ,
" Maneuvers (%d): \n " , sim - > maneuver_count ) ;
for ( int i = 0 ; i < sim - > maneuver_count ; i + + ) {
Maneuver * m = & sim - > maneuvers [ i ] ;
offset + = snprintf ( buffer + offset , buffer_size - offset ,
" [%d] %s: craft=%d dir=%d dv=%.4f trigger=%d val=%.2f exec=%d \n " ,
i , m - > name , m - > craft_index , m - > direction , m - > delta_v ,
m - > trigger_type , m - > trigger_value , m - > executed ) ;
if ( sim - > maneuvers ) {
for ( int i = 0 ; i < sim - > maneuver_count ; i + + ) {
Maneuver * m = & sim - > maneuvers [ i ] ;
offset + = snprintf ( buffer + offset , buffer_size - offset ,
" [%d] %s: craft=%d dir=%d dv=%.4f trigger=%d val=%.2f exec=%d \n " ,
i , m - > name , m - > craft_index , m - > direction , m - > delta_v ,
m - > trigger_type , m - > trigger_value , m - > executed ) ;
}
}
return offset ;