sc_neurocore_engine/neurons/trivial/
quadratic_if.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 — Quadratic Integrate-and-Fire Neuron
8
9/// Quadratic Integrate-and-Fire — canonical Type-I excitability.
10/// dv/dt = v² + I, reset at v_peak.
11#[derive(Clone, Debug)]
12pub struct QuadraticIFNeuron {
13    pub v: f64,
14    pub v_reset: f64,
15    pub v_peak: f64,
16    pub dt: f64,
17}
18
19impl QuadraticIFNeuron {
20    pub fn new(v_reset: f64, v_peak: f64, dt: f64) -> Self {
21        Self {
22            v: v_reset,
23            v_reset,
24            v_peak,
25            dt,
26        }
27    }
28
29    fn valid_numeric_contract(&self) -> bool {
30        self.v.is_finite()
31            && self.v_reset.is_finite()
32            && self.v_peak.is_finite()
33            && self.dt.is_finite()
34            && self.v < self.v_peak
35            && self.v_reset < self.v_peak
36            && self.dt > 0.0
37    }
38
39    pub fn step(&mut self, current: f64) -> i32 {
40        if !self.valid_numeric_contract() || !current.is_finite() {
41            return 0;
42        }
43        let (next_v, spiked) = self.exact_candidate(current);
44        if !next_v.is_finite() {
45            return 0;
46        }
47        self.v = next_v;
48        if spiked {
49            1
50        } else {
51            0
52        }
53    }
54
55    pub fn reset(&mut self) {
56        self.v = self.v_reset;
57    }
58
59    fn exact_candidate(&self, current: f64) -> (f64, bool) {
60        if current > 0.0 {
61            let root_i = current.sqrt();
62            let phase = (self.v / root_i).atan();
63            let peak_phase = (self.v_peak / root_i).atan();
64            let next_phase = phase + root_i * self.dt;
65            if next_phase >= peak_phase || next_phase >= std::f64::consts::FRAC_PI_2 {
66                return (self.v_reset, true);
67            }
68            return (root_i * next_phase.tan(), false);
69        }
70        if current == 0.0 {
71            let denominator = 1.0 - self.v * self.dt;
72            if denominator <= 0.0 {
73                return (self.v_reset, true);
74            }
75            let next_v = self.v / denominator;
76            if next_v >= self.v_peak {
77                return (self.v_reset, true);
78            }
79            return (next_v, false);
80        }
81
82        let root_i = (-current).sqrt();
83        if (self.v + root_i).abs() <= 1e-15 {
84            return (self.v, false);
85        }
86        let numerator_ratio = (self.v - root_i) / (self.v + root_i);
87        let evolved_ratio = numerator_ratio * (2.0 * root_i * self.dt).exp();
88        let denominator = 1.0 - evolved_ratio;
89        if (numerator_ratio < 1.0 && evolved_ratio >= 1.0) || denominator.abs() <= 1e-15 {
90            return (self.v_reset, true);
91        }
92        let next_v = root_i * (1.0 + evolved_ratio) / denominator;
93        if next_v >= self.v_peak {
94            (self.v_reset, true)
95        } else {
96            (next_v, false)
97        }
98    }
99}
100
101impl Default for QuadraticIFNeuron {
102    fn default() -> Self {
103        Self::new(-1.0, 1.0, 0.01)
104    }
105}
106
107#[cfg(test)]
108mod tests {
109    use super::*;
110
111    #[test]
112    fn qif_fires_with_positive_input() {
113        let mut n = QuadraticIFNeuron::default();
114        let total: i32 = (0..1000).map(|_| n.step(0.5)).sum();
115        assert!(total > 0);
116    }
117    #[test]
118    fn qif_silent_without_input() {
119        let mut n = QuadraticIFNeuron::default();
120        let t: i32 = (0..1000).map(|_| n.step(0.0)).sum();
121        assert_eq!(t, 0);
122    }
123    #[test]
124    fn qif_reset_clears_state() {
125        let mut n = QuadraticIFNeuron::default();
126        for _ in 0..100 {
127            n.step(0.5);
128        }
129        n.reset();
130        assert!((n.v - n.v_reset).abs() < 1e-10);
131    }
132    #[test]
133    fn qif_bounded() {
134        let mut n = QuadraticIFNeuron::default();
135        for _ in 0..1000 {
136            n.step(10.0);
137        }
138        assert!(n.v.is_finite());
139    }
140    #[test]
141    fn qif_nan_no_panic() {
142        let mut n = QuadraticIFNeuron::default();
143        let before = n.v;
144        assert_eq!(n.step(f64::NAN), 0);
145        assert_eq!(n.v, before);
146    }
147    #[test]
148    fn qif_nonfinite_increment_preserves_state() {
149        let mut n = QuadraticIFNeuron {
150            v: -0.25,
151            ..Default::default()
152        };
153        let before = n.v;
154        assert_eq!(n.step(-1.0e308), 0);
155        assert_eq!(n.v, before);
156    }
157    #[test]
158    fn qif_matches_exact_positive_current_flow() {
159        let mut n = QuadraticIFNeuron::default();
160        let root_i = 0.5_f64.sqrt();
161        let expected = root_i * ((n.v / root_i).atan() + root_i * n.dt).tan();
162        assert_eq!(n.step(0.5), 0);
163        assert!((n.v - expected).abs() < 1e-12);
164    }
165    #[test]
166    fn qif_preserves_negative_current_fixed_point() {
167        let mut n = QuadraticIFNeuron::default();
168        assert_eq!(n.step(-1.0), 0);
169        assert_eq!(n.v, -1.0);
170    }
171    #[test]
172    fn qif_exact_flow_resets_on_peak_crossing() {
173        let mut n = QuadraticIFNeuron {
174            v: 0.95,
175            dt: 0.5,
176            ..Default::default()
177        };
178        assert_eq!(n.step(1.0), 1);
179        assert_eq!(n.v, n.v_reset);
180    }
181    #[test]
182    fn qif_negative_no_crash() {
183        let mut n = QuadraticIFNeuron::default();
184        for _ in 0..500 {
185            n.step(-5.0);
186        }
187        assert!(n.v.is_finite());
188    }
189}
sc_neurocore_engine/neurons/trivial/quadratic_if.rs

sc_neurocore_engine/neurons/trivial/
quadratic_if.rs
