satellitesfactory 0.0.3

Satellite factory — classify, build and catalogue natural satellites for any planetary system: Solar System moons (Moon, Galileans, Titan, Triton…) or custom configurations.
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96

use crate::config::parameters::*;

#[derive(Debug, Clone, Copy, PartialEq, Eq)]
pub enum FormationChannel {
    CircumplanetaryDisk,
    GiantImpact,
    Capture,
    CoAccretion,
    RingSpreadout,
}

pub struct FormationModel {
    pub parent_mass: f64,
    pub parent_radius: f64,
    pub disk_mass_ratio: f64,
    pub disk_outer_edge: f64,
    pub disk_temperature_1r: f64,
}

impl FormationModel {
    pub fn new(
        parent_mass: f64,
        parent_radius: f64,
        disk_mass_ratio: f64,
        disk_outer_edge: f64,
        disk_temperature_1r: f64,
    ) -> Self {
        Self {
            parent_mass,
            parent_radius,
            disk_mass_ratio,
            disk_outer_edge,
            disk_temperature_1r,
        }
    }

    pub fn disk_mass(&self) -> f64 {
        self.parent_mass * self.disk_mass_ratio
    }

    pub fn disk_surface_density(&self, r: f64) -> f64 {
        let r_p = self.parent_radius;
        let sigma_0 = self.disk_mass() / (2.0 * std::f64::consts::PI * r_p * self.disk_outer_edge);
        sigma_0 * (r_p / r)
    }

    pub fn disk_temperature(&self, r: f64) -> f64 {
        self.disk_temperature_1r * (self.parent_radius / r).sqrt()
    }

    pub fn ice_line_radius(&self) -> f64 {
        self.parent_radius * (self.disk_temperature_1r / 170.0).powi(2)
    }

    pub fn isolation_mass(&self, r: f64) -> f64 {
        let sigma = self.disk_surface_density(r);
        let h = hill_radius(r, sigma * r * r, self.parent_mass);
        2.0 * std::f64::consts::PI * r * 2.0 * 3.5 * h * sigma
    }

    pub fn maximum_satellite_mass(&self) -> f64 {
        1e-4 * self.parent_mass
    }

    pub fn formation_timescale(&self, r: f64) -> f64 {
        let sigma = self.disk_surface_density(r);
        let omega = (G * self.parent_mass / r.powi(3)).sqrt();
        if sigma < 1e-30 || omega < 1e-30 {
            return f64::INFINITY;
        }
        let rho_solid = 3000.0;
        let r_body = 100.0e3;
        rho_solid * r_body / (sigma * omega)
    }

    pub fn capture_probability(v_inf: f64, v_esc: f64, gas_drag_factor: f64) -> f64 {
        if v_inf > v_esc {
            return 0.0;
        }
        let x = v_inf / (v_esc + 1e-30);
        gas_drag_factor * (1.0 - x * x)
    }

    pub fn impact_ejecta_fraction(v_impact: f64, v_esc_target: f64) -> f64 {
        let ratio = v_impact / (v_esc_target + 1e-30);
        if ratio < 1.0 {
            0.0
        } else {
            (1.0 - 1.0 / (ratio * ratio)).min(0.5)
        }
    }

    pub fn ring_spreading_timescale(ring_radius: f64, viscosity: f64) -> f64 {
        ring_radius * ring_radius / (viscosity + 1e-30)
    }
}