cometsfactory 0.0.3

use crate::config::parameters::{AU, G, SOLAR_MASS};

/// Gravitational binding energy of a spherical body (J).
/// E = 3 G M² / (5 R)
pub fn binding_energy(mass: f64, radius: f64) -> f64 {
    if radius < 1.0 {
        return 0.0;
    }
    3.0 * G * mass * mass / (5.0 * radius)
}

/// Tidal acceleration from the Sun at distance r (m) on a body of radius R.
/// a_tidal ≈ 2 G M_sun R / r³
pub fn solar_tidal_acceleration(body_radius: f64, distance_au: f64) -> f64 {
    let r = distance_au * AU;
    if r < 1.0 {
        return 0.0;
    }
    2.0 * G * SOLAR_MASS * body_radius / (r * r * r)
}

/// Hill sphere radius (m) at heliocentric distance a (AU) for body of mass m.
pub fn hill_radius(semi_major_au: f64, body_mass: f64) -> f64 {
    let a = semi_major_au * AU;
    a * (body_mass / (3.0 * SOLAR_MASS)).powf(1.0 / 3.0)
}

/// Roche limit distance (m) for a fluid satellite of given density
/// around the Sun.
/// d = R_sun * (2 ρ_sun / ρ_body)^(1/3)
/// We approximate using solar radius = 6.957e8 m, solar mean density = 1408 kg/m³
pub fn roche_limit_solar(body_density: f64) -> f64 {
    let r_sun = 6.957e8;
    let rho_sun = 1408.0;
    r_sun * (2.0 * rho_sun / body_density.max(1.0)).powf(1.0 / 3.0)
}

/// Orbital velocity at distance r_au from the Sun (m/s).
pub fn orbital_velocity(r_au: f64) -> f64 {
    let r = r_au * AU;
    if r < 1.0 {
        return 0.0;
    }
    (G * SOLAR_MASS / r).sqrt()
}

/// Specific orbital energy (J/kg) for a bound orbit with semi-major axis a (AU).
pub fn specific_orbital_energy(semi_major_au: f64) -> f64 {
    let a = semi_major_au * AU;
    if a < 1.0 {
        return 0.0;
    }
    -G * SOLAR_MASS / (2.0 * a)
}