Skip to main content

sc_neurocore_engine/neurons/rate/
fractional_lif.rs

1// SPDX-License-Identifier: AGPL-3.0-or-later
2// Commercial license available
3// © Concepts 1996–2026 Miroslav Šotek. All rights reserved.
4// © Code 2020–2026 Miroslav Šotek. All rights reserved.
5// ORCID: 0009-0009-3560-0851
6// Contact: www.anulum.li | protoscience@anulum.li
7// SC-NeuroCore — Fractional LIF neuron model
8
9/// Fractional-order LIF — Grünwald-Letnikov approximation. Teka et al. 2014.
10#[derive(Clone, Debug)]
11pub struct FractionalLIFNeuron {
12    pub v: f64,
13    pub v_rest: f64,
14    pub v_reset: f64,
15    pub v_threshold: f64,
16    pub alpha: f64,
17    pub resistance: f64,
18    pub dt: f64,
19    history: Vec<f64>,
20    gl_coeffs: Vec<f64>,
21    _max_hist: usize,
22}
23
24impl FractionalLIFNeuron {
25    pub fn new(alpha: f64, max_hist: usize) -> Self {
26        let mut coeffs = vec![0.0; max_hist + 1];
27        coeffs[0] = 1.0;
28        for j in 1..=max_hist {
29            coeffs[j] = coeffs[j - 1] * (1.0 - (alpha + 1.0) / j as f64);
30        }
31        Self {
32            v: 0.0,
33            v_rest: 0.0,
34            v_reset: 0.0,
35            v_threshold: 1.0,
36            alpha,
37            resistance: 1.0,
38            dt: 1.0,
39            history: vec![0.0; max_hist],
40            gl_coeffs: coeffs,
41            _max_hist: max_hist,
42        }
43    }
44    pub fn step(&mut self, current: f64) -> i32 {
45        // Grünwald-Letnikov: D^α v ≈ (1/dt^α) Σ_j c_j v(t-j·dt)
46        let mut gl_sum = 0.0;
47        let n = self.history.len().min(self.gl_coeffs.len() - 1);
48        for j in 0..n {
49            gl_sum += self.gl_coeffs[j + 1] * self.history[n - 1 - j];
50        }
51        let rhs = -(self.v - self.v_rest) + self.resistance * current;
52        self.v = rhs * self.dt.powf(self.alpha) - gl_sum;
53        // Shift history
54        let len = self.history.len();
55        if len > 0 {
56            for i in 0..len - 1 {
57                self.history[i] = self.history[i + 1];
58            }
59            self.history[len - 1] = self.v;
60        }
61        if self.v >= self.v_threshold {
62            self.v = self.v_reset;
63            1
64        } else {
65            0
66        }
67    }
68    pub fn reset(&mut self) {
69        self.v = self.v_rest;
70        self.history.fill(0.0);
71    }
72}
73
74#[cfg(test)]
75mod tests {
76    use super::*;
77
78    #[test]
79    fn frac_lif_fires() {
80        let mut n = FractionalLIFNeuron::new(0.8, 50);
81        let t: i32 = (0..200).map(|_| n.step(2.0)).sum();
82        assert!(t > 0);
83    }
84
85    #[test]
86    fn frac_lif_reset() {
87        let mut n = FractionalLIFNeuron::new(0.8, 50);
88        for _ in 0..100 {
89            n.step(2.0);
90        }
91        n.reset();
92        assert!((n.v - n.v_rest).abs() < 1e-10);
93    }
94
95    #[test]
96    fn frac_lif_nan_no_panic() {
97        FractionalLIFNeuron::new(0.8, 50).step(f64::NAN);
98    }
99}