suns 0.0.1

Sun celestial simulation crate for the MilkyWay SolarSystem workspace
Documentation
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)
}