sciforge 0.0.3

A comprehensive scientific computing library in pure Rust with zero dependencies
Documentation
pub fn metabolic_network_steady_state(
    stoich: &[Vec<f64>],
    rates: impl Fn(&[f64]) -> Vec<f64>,
    initial: &[f64],
    dt: f64,
    steps: usize,
    tolerance: f64,
) -> Vec<f64> {
    let mut conc = initial.to_vec();
    for _ in 0..steps {
        let v = rates(&conc);
        let mut max_change = 0.0_f64;
        for (i, stoich_row) in stoich.iter().enumerate() {
            let mut dc = 0.0;
            for (j, &vj) in v.iter().enumerate() {
                dc += stoich_row[j] * vj;
            }
            conc[i] += dc * dt;
            conc[i] = conc[i].max(0.0);
            max_change = max_change.max(dc.abs() * dt);
        }
        if max_change < tolerance {
            break;
        }
    }
    conc
}

pub fn arrhenius(a: f64, ea: f64, t: f64) -> f64 {
    use crate::constants::K_B;
    a * (-ea / (K_B * t)).exp()
}

pub fn q10_factor(rate1: f64, rate2: f64, t1: f64, t2: f64) -> f64 {
    (rate2 / rate1).powf(10.0 / (t2 - t1))
}

pub fn metabolic_control_coefficient(
    flux_perturbed: f64,
    flux_original: f64,
    enzyme_perturbed: f64,
    enzyme_original: f64,
) -> f64 {
    let df = (flux_perturbed - flux_original) / flux_original;
    let de = (enzyme_perturbed - enzyme_original) / enzyme_original;
    if de.abs() < 1e-30 {
        return 0.0;
    }
    df / de
}

pub fn gibbs_free_energy(delta_g0: f64, t: f64, q: f64) -> f64 {
    use crate::constants::{K_B, N_A};
    let r = K_B * N_A;
    delta_g0 + r * t * q.ln()
}

pub fn mass_action_ratio(
    products: &[f64],
    reactants: &[f64],
    stoich_p: &[f64],
    stoich_r: &[f64],
) -> f64 {
    let num: f64 = products
        .iter()
        .zip(stoich_p.iter())
        .map(|(&c, &s)| c.powf(s))
        .product();
    let den: f64 = reactants
        .iter()
        .zip(stoich_r.iter())
        .map(|(&c, &s)| c.powf(s))
        .product();
    if den < 1e-30 {
        return f64::INFINITY;
    }
    num / den
}

pub fn reaction_quotient_vs_keq(q: f64, keq: f64) -> f64 {
    (keq / q.max(1e-30)).ln()
}

pub fn flux_control_summation(coefficients: &[f64]) -> f64 {
    coefficients.iter().sum()
}

pub fn elasticity_coefficient(
    rate: f64,
    metabolite: f64,
    delta_rate: f64,
    delta_metabolite: f64,
) -> f64 {
    if metabolite.abs() < 1e-30 || rate.abs() < 1e-30 {
        return 0.0;
    }
    (delta_rate / rate) / (delta_metabolite / metabolite)
}

pub fn supply_demand_modular(
    supply_flux: f64,
    demand_flux: f64,
    linking_metabolite: f64,
    elasticity_supply: f64,
    elasticity_demand: f64,
) -> f64 {
    let deviation = supply_flux - demand_flux;
    linking_metabolite + deviation / (elasticity_demand - elasticity_supply).abs().max(1e-30)
}