darkmatter 0.0.1

Dark matter simulation engine — gravitational fields, particle dynamics, halo stability, and cosmological constants
Documentation
use crate::constants::physical::{C, G, SOLAR_MASS};
use std::f64::consts::PI;

pub fn capture_rate_sun(rho_local_kg_m3: f64, v_dm: f64, sigma_si_cm2: f64, m_dm_kg: f64) -> f64 {
    let r_sun = 6.957e8;
    let m_sun = SOLAR_MASS;
    let m_h = 1.673e-27;
    let n_h = m_sun / m_h;
    let sigma_si = sigma_si_cm2 * 1e-4;
    let n_dm = rho_local_kg_m3 / m_dm_kg;
    let v_esc_sun = (2.0 * G * m_sun / r_sun).sqrt();
    let geometric = PI * r_sun * r_sun;
    let grav_focus = 1.0 + (v_esc_sun / v_dm).powi(2);
    let prob = (n_h * sigma_si / (PI * r_sun * r_sun)).min(1.0);
    n_dm * v_dm * geometric * grav_focus * prob
}

pub fn capture_rate_planet(
    rho_local_kg_m3: f64,
    v_dm: f64,
    sigma_si_cm2: f64,
    m_dm_kg: f64,
    planet_mass: f64,
    planet_radius: f64,
) -> f64 {
    let m_nucleon = 1.673e-27;
    let n_nucleons = planet_mass / m_nucleon;
    let sigma_si = sigma_si_cm2 * 1e-4;
    let n_dm = rho_local_kg_m3 / m_dm_kg;
    let v_esc = (2.0 * G * planet_mass / planet_radius).sqrt();
    let geometric = PI * planet_radius * planet_radius;
    let grav_focus = 1.0 + (v_esc / v_dm).powi(2);
    let prob = (n_nucleons * sigma_si / (PI * planet_radius * planet_radius)).min(1.0);
    n_dm * v_dm * geometric * grav_focus * prob
}

pub fn equilibrium_annihilation_rate(capture_rate: f64, sigma_v: f64, v_eff: f64) -> f64 {
    let ca = sigma_v / v_eff.powi(3);
    let tau = 1.0 / (capture_rate * ca).sqrt();
    0.5 * capture_rate * (capture_rate * ca * tau * tau / 1.0).tanh()
}

pub fn dm_accumulated_number(
    capture_rate: f64,
    age_seconds: f64,
    sigma_v: f64,
    volume: f64,
) -> f64 {
    let ann_rate = sigma_v / volume;
    let tau_eq = 1.0 / (capture_rate * ann_rate).sqrt();
    let t_ratio = age_seconds / tau_eq;
    (capture_rate / ann_rate).sqrt() * t_ratio.tanh()
}

pub fn dm_density_inside_body(
    n_captured: f64,
    body_radius: f64,
    thermal_radius_fraction: f64,
) -> f64 {
    let r_th = body_radius * thermal_radius_fraction;
    3.0 * n_captured / (4.0 * PI * r_th * r_th * r_th)
}

pub fn heating_from_capture(capture_rate: f64, m_dm_kg: f64) -> f64 {
    capture_rate * m_dm_kg * C * C
}

pub fn neutrino_signal_from_capture(capture_rate: f64, branching_nu: f64, n_nu: f64) -> f64 {
    capture_rate * branching_nu * n_nu
}

pub fn gravitational_focusing_factor(v_dm: f64, v_esc: f64) -> f64 {
    1.0 + (v_esc / v_dm).powi(2)
}