Skip to main content

sc_neurocore_engine/neurons/biophysical/
plant_r15.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 — Plant R15 Neuron Model
8
9//! Plant R15 Aplysia parabolic-burster dynamics.
10
11use super::safe_rate;
12
13/// Plant R15 — Aplysia parabolic burster. Plant & Kim 1976.
14#[derive(Clone, Debug)]
15pub struct PlantR15Neuron {
16    pub v: f64,
17    pub m: f64,
18    pub h: f64,
19    pub n: f64,
20    pub ca: f64,
21    pub g_na: f64,
22    pub g_k: f64,
23    pub g_ca: f64,
24    pub g_l: f64,
25    pub g_kca: f64,
26    pub e_na: f64,
27    pub e_k: f64,
28    pub e_ca: f64,
29    pub e_l: f64,
30    pub c_m: f64,
31    pub k_ca: f64,
32    pub tau_ca: f64,
33    pub dt: f64,
34    pub v_threshold: f64,
35}
36
37impl PlantR15Neuron {
38    pub fn new() -> Self {
39        Self {
40            v: -50.0,
41            m: 0.05,
42            h: 0.6,
43            n: 0.3,
44            ca: 0.1,
45            g_na: 4.0,
46            g_k: 0.3,
47            g_ca: 0.004,
48            g_l: 0.003,
49            g_kca: 0.03,
50            e_na: 30.0,
51            e_k: -75.0,
52            e_ca: 140.0,
53            e_l: -40.0,
54            c_m: 1.0,
55            k_ca: 0.0085,
56            tau_ca: 500.0,
57            dt: 0.05,
58            v_threshold: -10.0,
59        }
60    }
61    pub fn step(&mut self, current: f64) -> i32 {
62        let v_prev = self.v;
63        for _ in 0..5 {
64            let am = safe_rate(0.1, 50.0, self.v, 10.0, 1.0);
65            let bm = 4.0 * (-(self.v + 75.0) / 18.0).exp();
66            let ah = 0.07 * (-(self.v + 50.0) / 20.0).exp();
67            let bh = 1.0 / (1.0 + (-(self.v + 20.0) / 10.0).exp());
68            let an = safe_rate(0.01, 55.0, self.v, 10.0, 0.1);
69            let bn = 0.125 * (-(self.v + 65.0) / 80.0).exp();
70            self.m += (am * (1.0 - self.m) - bm * self.m) * self.dt;
71            self.h += (ah * (1.0 - self.h) - bh * self.h) * self.dt;
72            self.n += (an * (1.0 - self.n) - bn * self.n) * self.dt;
73            let m_ca = 1.0 / (1.0 + (-(self.v + 25.0) / 5.0).exp());
74            let kca_act = self.ca / (0.5 + self.ca);
75            let i_na = self.g_na * self.m.powi(3) * self.h * (self.v - self.e_na);
76            let i_k = self.g_k * self.n.powi(4) * (self.v - self.e_k);
77            let i_ca = self.g_ca * m_ca.powi(2) * (self.v - self.e_ca);
78            let i_kca = self.g_kca * kca_act * (self.v - self.e_k);
79            let i_l = self.g_l * (self.v - self.e_l);
80            self.v += (-i_na - i_k - i_ca - i_kca - i_l + current) / self.c_m * self.dt;
81            self.ca = (self.ca + (-self.k_ca * i_ca - self.ca / self.tau_ca) * self.dt).max(0.0);
82        }
83        if self.v >= self.v_threshold && v_prev < self.v_threshold {
84            1
85        } else {
86            0
87        }
88    }
89    pub fn reset(&mut self) {
90        self.v = -50.0;
91        self.m = 0.05;
92        self.h = 0.6;
93        self.n = 0.3;
94        self.ca = 0.1;
95    }
96}
97impl Default for PlantR15Neuron {
98    fn default() -> Self {
99        Self::new()
100    }
101}
102
103#[cfg(test)]
104mod tests {
105    use super::*;
106
107    #[test]
108    fn default_matches_constructor_state() {
109        let default = PlantR15Neuron::default();
110        let constructed = PlantR15Neuron::new();
111        assert_eq!(default.v, constructed.v);
112    }
113
114    #[test]
115    fn plant_r15_fires() {
116        let mut n = PlantR15Neuron::new();
117        let t: i32 = (0..500).map(|_| n.step(2.0)).sum();
118        assert!(t > 0);
119    }
120
121    // -- PlantR15 --
122    #[test]
123    fn plant_r15_silent_without_input() {
124        let mut n = PlantR15Neuron::new();
125        // R15 is a parabolic burster — may burst spontaneously
126        for _ in 0..500 {
127            n.step(0.0);
128        }
129        assert!(n.v.is_finite());
130    }
131    #[test]
132    fn plant_r15_reset_clears_state() {
133        let mut n = PlantR15Neuron::new();
134        for _ in 0..100 {
135            n.step(2.0);
136        }
137        n.reset();
138        assert!((n.v - (-50.0)).abs() < 1e-10);
139        assert!((n.ca - 0.1).abs() < 1e-10);
140    }
141    #[test]
142    fn plant_r15_moderate_input_stable() {
143        // Plant R15 is a parabolic burster — moderate input stability
144        let mut n = PlantR15Neuron::new();
145        for _ in 0..500 {
146            n.step(2.0);
147        }
148        assert!(n.v.is_finite());
149    }
150    #[test]
151    fn plant_r15_ca_dynamics() {
152        let mut n = PlantR15Neuron::new();
153        for _ in 0..500 {
154            n.step(2.0);
155        }
156        assert!(n.ca >= 0.0, "Ca²⁺ must be non-negative");
157        assert!(n.ca.is_finite());
158    }
159    #[test]
160    fn plant_r15_weak_negative_no_crash() {
161        let mut n = PlantR15Neuron::new();
162        for _ in 0..200 {
163            n.step(-1.0);
164        }
165        assert!(n.v.is_finite());
166    }
167    #[test]
168    fn plant_r15_nan_no_panic() {
169        let mut n = PlantR15Neuron::new();
170        n.step(f64::NAN);
171    }
172}