darkmatter 0.0.2

Dark matter simulation engine — gravitational fields, particle dynamics, halo stability, and cosmological constants
Documentation
use crate::constants::cosmic::OMEGA_BARYON;
use crate::constants::physical::HBAR;
use std::f64::consts::PI;

#[derive(Clone, Copy, Debug, PartialEq)]
pub struct AsymmetricDarkMatter {
    pub mass_gev: f64,
    pub asymmetry_ratio: f64,
}

impl AsymmetricDarkMatter {
    pub fn mass_kg(&self) -> f64 {
        self.mass_gev * 1.782_661_92e-27
    }

    pub fn number_density_cm3(&self, rho_gev_cm3: f64) -> f64 {
        rho_gev_cm3 / self.mass_gev
    }

    pub fn relic_density(&self) -> f64 {
        OMEGA_BARYON * self.asymmetry_ratio * self.mass_gev / 0.938_272
    }

    pub fn annihilation_freeze_out_temperature_gev(&self) -> f64 {
        self.mass_gev / 20.0
    }

    pub fn symmetric_component_depletion(&self, sigma_v_cm3_s: f64) -> f64 {
        let x_fo = self.mass_gev / self.annihilation_freeze_out_temperature_gev();
        let n_eq_ratio = (-x_fo).exp();
        n_eq_ratio * sigma_v_cm3_s / 3e-26
    }
}

pub fn adm_mass_estimate_gev(omega_dm: f64, omega_b: f64) -> f64 {
    0.938_272 * omega_dm / omega_b
}

pub fn baryon_dm_ratio(omega_dm: f64, omega_b: f64) -> f64 {
    omega_dm / omega_b
}

pub fn nugget_binding_energy_gev(n_quarks: f64, lambda_qcd_gev: f64) -> f64 {
    n_quarks.powf(2.0 / 3.0) * lambda_qcd_gev
}

pub fn composite_adm_radius_fm(n_constituents: f64, constituent_mass_gev: f64) -> f64 {
    let r_bohr = HBAR / (constituent_mass_gev * 1.782_661_92e-27 * 137.036 * 3e8);
    r_bohr * n_constituents.powf(1.0 / 3.0) * 1e15
}

pub fn adm_self_interaction_cross_section_cm2(radius_fm: f64) -> f64 {
    let r_cm = radius_fm * 1e-13;
    4.0 * PI * r_cm * r_cm
}

pub fn bbn_constraint_mass_lower_gev(omega_dm: f64, omega_b: f64) -> f64 {
    let ratio = omega_dm / omega_b;
    0.938_272 * ratio * 0.1
}

pub fn direct_detection_sigma_cm2(mass_gev: f64, mediator_mass_gev: f64, coupling: f64) -> f64 {
    let mu = mass_gev * 0.938_272 / (mass_gev + 0.938_272);
    coupling.powi(4) * mu * mu * 1.782_661_92e-27 * 1.782_661_92e-27
        / (PI * mediator_mass_gev.powi(4) * 1.782_661_92e-27 * 1.782_661_92e-27)
        * 1e-36
}

pub fn mirror_dm_temperature_ratio(epsilon: f64) -> f64 {
    epsilon.powf(1.0 / 3.0)
}

pub fn dark_atom_binding_energy_ev(alpha_dark: f64, mass_dark_electron_gev: f64) -> f64 {
    0.5 * alpha_dark * alpha_dark * mass_dark_electron_gev * 1e9
}

pub fn dark_recombination_redshift(alpha_dark: f64, mass_dark_electron_gev: f64, xi: f64) -> f64 {
    let binding_ev = dark_atom_binding_energy_ev(alpha_dark, mass_dark_electron_gev);
    let binding_k = binding_ev * 1.160_451e4;
    binding_k / (xi * 25.0)
}

pub fn adm_capture_rate_in_star(
    rho_dm_gev_cm3: f64,
    sigma_cm2: f64,
    v_star_escape_km_s: f64,
    v_dm_km_s: f64,
    mass_star_solar: f64,
    dm_mass_gev: f64,
) -> f64 {
    let v_esc = v_star_escape_km_s * 1e5;
    let v_dm = v_dm_km_s * 1e5;
    let n_dm = rho_dm_gev_cm3 / dm_mass_gev;
    let r_star = 6.957e10 * (mass_star_solar).powf(0.8);
    let geometric = PI * r_star * r_star;
    let grav_focus = 1.0 + (v_esc / v_dm).powi(2);
    n_dm * sigma_cm2 * v_dm * geometric * grav_focus
}