suns 0.0.4

Sun celestial simulation crate for the MilkyWay SolarSystem workspace
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
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140

use sciforge::hub::domain::common::constants::K_B;

use crate::{MU_0, PROTON_MASS, SOLAR_SURFACE_GRAVITY};

pub const PROMINENCE_TEMPERATURE: f64 = 8000.0;
pub const PROMINENCE_DENSITY: f64 = 1e-10;
pub const PROMINENCE_LENGTH_M: f64 = 2e8;
pub const PROMINENCE_HEIGHT_M: f64 = 5e7;
pub const PROMINENCE_WIDTH_M: f64 = 5e6;
pub const PROMINENCE_MAGNETIC_FIELD_T: f64 = 1e-3;
pub const PROMINENCE_LIFETIME_S: f64 = 86400.0 * 30.0;
pub const FILAMENT_CHANNEL_LENGTH_M: f64 = 3e8;

pub struct Prominence {
    pub height_m: f64,
    pub length_m: f64,
    pub width_m: f64,
    pub temperature_k: f64,
    pub density_kgm3: f64,
    pub magnetic_field_t: f64,
}

impl Default for Prominence {
    fn default() -> Self {
        Self::quiescent()
    }
}

impl Prominence {
    pub fn quiescent() -> Self {
        Self {
            height_m: PROMINENCE_HEIGHT_M,
            length_m: PROMINENCE_LENGTH_M,
            width_m: PROMINENCE_WIDTH_M,
            temperature_k: PROMINENCE_TEMPERATURE,
            density_kgm3: PROMINENCE_DENSITY,
            magnetic_field_t: PROMINENCE_MAGNETIC_FIELD_T,
        }
    }

    pub fn active() -> Self {
        Self {
            height_m: 3e7,
            length_m: 1e8,
            width_m: 3e6,
            temperature_k: 1.5e4,
            density_kgm3: 5e-10,
            magnetic_field_t: 5e-3,
        }
    }

    pub fn volume_m3(&self) -> f64 {
        self.height_m * self.length_m * self.width_m
    }

    pub fn mass_kg(&self) -> f64 {
        self.density_kgm3 * self.volume_m3()
    }

    pub fn thermal_energy(&self) -> f64 {
        let n = self.density_kgm3 / PROTON_MASS;
        1.5 * n * K_B * self.temperature_k * self.volume_m3()
    }

    pub fn magnetic_energy(&self) -> f64 {
        self.magnetic_field_t.powi(2) * self.volume_m3() / (2.0 * MU_0)
    }

    pub fn gravitational_potential_energy(&self) -> f64 {
        self.mass_kg() * SOLAR_SURFACE_GRAVITY * self.height_m
    }

    pub fn magnetic_tension_force(&self) -> f64 {
        self.magnetic_field_t.powi(2) / (MU_0 * self.length_m)
    }

    pub fn magnetic_pressure_gradient(&self) -> f64 {
        self.magnetic_field_t.powi(2) / (2.0 * MU_0 * self.width_m)
    }

    pub fn kippenhahn_schluter_support(&self) -> f64 {
        let weight = self.density_kgm3 * SOLAR_SURFACE_GRAVITY;
        let magnetic_support = self.magnetic_field_t.powi(2) / (MU_0 * self.height_m);
        magnetic_support / weight
    }

    pub fn thermal_instability_criterion(&self) -> bool {
        let corona_temp = 1e6;
        self.temperature_k < corona_temp * 0.1
    }

    pub fn alfven_speed(&self) -> f64 {
        self.magnetic_field_t / (MU_0 * self.density_kgm3).sqrt()
    }

    pub fn sound_speed(&self) -> f64 {
        let gamma = 5.0 / 3.0;
        let mu = 0.6;
        (gamma * K_B * self.temperature_k / (mu * PROTON_MASS)).sqrt()
    }

    pub fn free_fall_time(&self) -> f64 {
        (2.0 * self.height_m / SOLAR_SURFACE_GRAVITY).sqrt()
    }

    pub fn eruption_speed(&self) -> f64 {
        let e_mag = self.magnetic_energy();
        let m = self.mass_kg();
        (2.0 * e_mag / m).sqrt().min(1e6)
    }

    pub fn cooling_time(&self) -> f64 {
        let n = self.density_kgm3 / PROTON_MASS;
        let lambda = 1e-23;
        3.0 * K_B * self.temperature_k / (n * lambda)
    }

    pub fn emission_measure(&self) -> f64 {
        let n_e = self.density_kgm3 / PROTON_MASS;
        n_e * n_e * self.volume_m3()
    }
}

pub fn drain_rate(density: f64, area: f64, gravity: f64, height: f64) -> f64 {
    let v = (2.0 * gravity * height).sqrt();
    density * area * v
}

pub fn prominence_oscillation_period(length: f64, field: f64, density: f64) -> f64 {
    let va = field / (MU_0 * density).sqrt();
    2.0 * length / va
}

pub fn coronal_rain_velocity(height: f64) -> f64 {
    (2.0 * SOLAR_SURFACE_GRAVITY * height).sqrt()
}

pub fn magnetic_dip_angle(bx: f64, bz: f64) -> f64 {
    bz.atan2(bx)
}