darkmatter 0.0.2

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

pub fn missing_energy_cross_section(coupling: f64, mediator_mass: f64, sqrt_s: f64) -> f64 {
    let s = sqrt_s * sqrt_s;
    coupling.powi(4) / (16.0 * PI * (s - mediator_mass * mediator_mass).powi(2))
}

pub fn mono_jet_production_rate(
    coupling: f64,
    mediator_mass: f64,
    sqrt_s: f64,
    luminosity: f64,
) -> f64 {
    let sigma = missing_energy_cross_section(coupling, mediator_mass, sqrt_s);
    sigma * luminosity
}

pub fn mediator_width_to_dm(coupling: f64, mediator_mass: f64, m_dm: f64) -> f64 {
    if mediator_mass < 2.0 * m_dm {
        return 0.0;
    }
    let beta = (1.0 - 4.0 * m_dm * m_dm / (mediator_mass * mediator_mass)).sqrt();
    coupling * coupling * mediator_mass * beta.powi(3) / (12.0 * PI)
}

pub fn mediator_width_to_quarks(
    coupling: f64,
    mediator_mass: f64,
    m_quark: f64,
    n_colors: f64,
) -> f64 {
    if mediator_mass < 2.0 * m_quark {
        return 0.0;
    }
    let beta = (1.0 - 4.0 * m_quark * m_quark / (mediator_mass * mediator_mass)).sqrt();
    n_colors * coupling * coupling * mediator_mass * beta / (12.0 * PI)
}

pub fn relic_density_from_sigma_v(sigma_v_cm3_s: f64) -> f64 {
    let thermal = 3e-26;
    0.12 * thermal / sigma_v_cm3_s
}

pub fn dm_nuclear_sigma_from_collider(
    coupling_dm: f64,
    coupling_q: f64,
    mediator_mass: f64,
    m_dm: f64,
    m_nucleon: f64,
) -> f64 {
    let mu = m_dm * m_nucleon / (m_dm + m_nucleon);
    coupling_dm * coupling_dm * coupling_q * coupling_q * mu * mu / (PI * mediator_mass.powi(4))
}

pub fn eft_validity_check(mediator_mass: f64, sqrt_s: f64) -> bool {
    mediator_mass > sqrt_s
}

#[derive(Clone, Copy, Debug, PartialEq)]
pub struct SimplifiedModel {
    pub mediator_mass: f64,
    pub coupling_dm: f64,
    pub coupling_sm: f64,
    pub m_dm: f64,
    pub spin: &'static str,
}

impl SimplifiedModel {
    pub fn total_width(&self, m_quark: f64) -> f64 {
        let w_dm = mediator_width_to_dm(self.coupling_dm, self.mediator_mass, self.m_dm);
        let w_qq = mediator_width_to_quarks(self.coupling_sm, self.mediator_mass, m_quark, 3.0);
        w_dm + 6.0 * w_qq
    }

    pub fn branching_ratio_dm(&self, m_quark: f64) -> f64 {
        let w_dm = mediator_width_to_dm(self.coupling_dm, self.mediator_mass, self.m_dm);
        let total = self.total_width(m_quark);
        if total > 0.0 { w_dm / total } else { 0.0 }
    }

    pub fn branching_ratio_visible(&self, m_quark: f64) -> f64 {
        1.0 - self.branching_ratio_dm(m_quark)
    }
}