sciforge 0.0.3

A comprehensive scientific computing library in pure Rust with zero dependencies
Documentation
pub fn waddington_landscape(x: f64, param: f64) -> f64 {
    x.powi(4) - param * x.powi(2)
}

pub fn differentiation_potential(x: f64, param: f64) -> f64 {
    4.0 * x.powi(3) - 2.0 * param * x
}

pub fn differentiation_simulate(
    x0: f64,
    param_start: f64,
    param_end: f64,
    noise: f64,
    dt: f64,
    steps: usize,
) -> Vec<(f64, f64)> {
    let mut result = Vec::with_capacity(steps + 1);
    let mut x = x0;
    result.push((param_start, x));
    for step in 0..steps {
        let t_frac = step as f64 / steps as f64;
        let param = param_start + (param_end - param_start) * t_frac;
        let force = -differentiation_potential(x, param);
        let pseudo_noise = noise * ((step as f64 * 0.618).sin() * 2.0 - 1.0);
        x += (force + pseudo_noise) * dt;
        result.push((param, x));
    }
    result
}

pub fn somitogenesis_clock(phase: f64, omega: f64, coupling: f64, neighbor_phase: f64) -> f64 {
    omega + coupling * (neighbor_phase - phase).sin()
}

pub fn somite_clock_simulate(
    phases: &mut [f64],
    omega: f64,
    coupling: f64,
    dt: f64,
    steps: usize,
) -> Vec<Vec<f64>> {
    let n = phases.len();
    let mut result = Vec::with_capacity(steps + 1);
    result.push(phases.to_vec());
    for _ in 0..steps {
        let old = phases.to_vec();
        for i in 0..n {
            let left = if i > 0 { old[i - 1] } else { old[i] };
            let right = if i < n - 1 { old[i + 1] } else { old[i] };
            let avg_neighbor = (left + right) / 2.0;
            let dphase = somitogenesis_clock(old[i], omega, coupling, avg_neighbor);
            phases[i] = old[i] + dphase * dt;
        }
        result.push(phases.to_vec());
    }
    result
}

pub fn lateral_inhibition_step(
    signal: f64,
    neighbor_signal: f64,
    delta: f64,
    notch: f64,
    k: f64,
) -> (f64, f64) {
    let ddelta = -delta + signal;
    let dnotch = -notch + k * neighbor_signal;
    (ddelta, dnotch)
}

pub fn cell_fate_probability(signal: f64, threshold: f64, steepness: f64) -> f64 {
    1.0 / (1.0 + (-(signal - threshold) * steepness).exp())
}