Skip to main content

sc_neurocore_engine/
network_runner.rs

1// SPDX-License-Identifier: AGPL-3.0-or-later | Commercial license available
2// © Concepts 1996–2026 Miroslav Šotek. All rights reserved.
3// © Code 2020–2026 Miroslav Šotek. All rights reserved.
4// ORCID: 0009-0009-3560-0851
5// Contact: www.anulum.li | protoscience@anulum.li
6// SC-NeuroCore — Network runner: high-performance Rust simulation backend
7
8//! High-performance network simulation backend.
9//!
10//! Replaces the Python per-neuron loop with Rayon-parallel Rust execution
11//! over CSR-stored projections and heterogeneous neuron populations.
12
13use rayon::prelude::*;
14
15use crate::neuron::*;
16use crate::neurons::*;
17
18// ── NeuronVariant ───────────────────────────────────────────────────
19
20/// Enum dispatch across all neuron models that accept `step(f64) -> i32`.
21///
22/// Models with non-standard signatures (multi-input, integer-input,
23/// rate-output) are excluded — they require specialized projection wiring.
24#[allow(clippy::large_enum_variant)]
25pub enum NeuronVariant {
26    // neuron.rs
27    Izhikevich(Izhikevich),
28    AdEx(AdExNeuron),
29    ExpIF(ExpIfNeuron),
30    Lapicque(LapicqueNeuron),
31    HomeostaticLif(HomeostaticLif),
32
33    // biophysical.rs
34    HodgkinHuxley(HodgkinHuxleyNeuron),
35    TraubMiles(TraubMilesNeuron),
36    WangBuzsaki(WangBuzsakiNeuron),
37    ConnorStevens(ConnorStevensNeuron),
38    DestexheThalamic(DestexheThalamicNeuron),
39    HuberBraun(HuberBraunNeuron),
40    GolombFS(GolombFSNeuron),
41    Pospischil(PospischilNeuron),
42    MainenSejnowski(MainenSejnowskiNeuron),
43    DeSchutterPurkinje(DeSchutterPurkinjeNeuron),
44    PlantR15(PlantR15Neuron),
45    Prescott(PrescottNeuron),
46    MihalasNiebur(MihalasNieburNeuron),
47    GLIF(GLIFNeuron),
48    GIFPopulation(GIFPopulationNeuron),
49    AvRonCardiac(AvRonCardiacNeuron),
50    DurstewitzDopamine(DurstewitzDopamineNeuron),
51    HillTononi(HillTononiNeuron),
52    BertramPhantom(BertramPhantomBurster),
53    Yamada(YamadaNeuron),
54
55    // simple_spiking.rs
56    FitzHughNagumo(FitzHughNagumoNeuron),
57    MorrisLecar(MorrisLecarNeuron),
58    HindmarshRose(HindmarshRoseNeuron),
59    ResonateAndFire(ResonateAndFireNeuron),
60    FitzHughRinzel(FitzHughRinzelNeuron),
61    McKean(McKeanNeuron),
62    TermanWang(TermanWangOscillator),
63    GutkinErmentrout(GutkinErmentroutNeuron),
64    WilsonHR(WilsonHRNeuron),
65    Chay(ChayNeuron),
66    ChayKeizer(ChayKeizerNeuron),
67    ShermanRinzelKeizer(ShermanRinzelKeizerNeuron),
68    ButeraRespiratory(ButeraRespiratoryNeuron),
69    EPropALIF(EPropALIFNeuron),
70    SuperSpike(SuperSpikeNeuron),
71    LearnableNeuron(LearnableNeuronModel),
72    Pernarowski(PernarowskiNeuron),
73
74    // trivial.rs (simple IF variants)
75    QuadraticIF(QuadraticIFNeuron),
76    Theta(ThetaNeuron),
77    PerfectIntegrator(PerfectIntegratorNeuron),
78    GatedLIF(GatedLIFNeuron),
79    NonlinearLIF(NonlinearLIFNeuron),
80    SFA(SFANeuron),
81    MAT(MATNeuron),
82    KLIF(KLIFNeuron),
83    InhibitoryLIF(InhibitoryLIFNeuron),
84    ComplementaryLIF(ComplementaryLIFNeuron),
85    ParametricLIF(ParametricLIFNeuron),
86    NonResettingLIF(NonResettingLIFNeuron),
87    AdaptiveThresholdIF(AdaptiveThresholdIFNeuron),
88    SigmaDelta(SigmaDeltaNeuron),
89    EnergyLIF(EnergyLIFNeuron),
90    ClosedFormContinuous(ClosedFormContinuousNeuron),
91
92    // maps.rs
93    ChialvoMap(ChialvoMapNeuron),
94    RulkovMap(RulkovMapNeuron),
95    IbarzTanakaMap(IbarzTanakaMapNeuron),
96    MedvedevMap(MedvedevMapNeuron),
97    CazellesMap(CazellesMapNeuron),
98    CourageNekorkinMap(CourageNekorkinMapNeuron),
99
100    // hardware.rs (f64 input subset)
101    BrainScaleSAdEx(BrainScaleSAdExNeuron),
102    SpiNNakerLIF(SpiNNakerLIFNeuron),
103    NeuroGrid(NeuroGridNeuron),
104    DPI(DPINeuron),
105
106    // multi_compartment.rs (single-f64-input subset)
107    MarderSTG(MarderSTGNeuron),
108    RallCable(RallCableNeuron),
109    BoothRinzel(BoothRinzelNeuron),
110    Dendrify(DendrifyNeuron),
111
112    // rate.rs (spiking subset with step(f64)->i32)
113    LiquidTimeConstant(LiquidTimeConstantNeuron),
114    ParallelSpiking(ParallelSpikingNeuron),
115    FractionalLIF(FractionalLIFNeuron),
116
117    // special.rs (step(f64)->i32 subset)
118    StochasticIF(StochasticIFNeuron),
119    GalvesLocherbach(GalvesLocherbachNeuron),
120    SpikeResponse(SpikeResponseNeuron),
121    GLM(GLMNeuron),
122    Arcane(ArcaneNeuron),
123}
124
125macro_rules! dispatch_step {
126    ($self:expr, $current:expr,
127     $($variant:ident),* $(,)?) => {
128        match $self {
129            $(NeuronVariant::$variant(n) => n.step($current),)*
130        }
131    };
132}
133
134macro_rules! dispatch_reset {
135    ($self:expr,
136     $($variant:ident),* $(,)?) => {
137        match $self {
138            $(NeuronVariant::$variant(n) => n.reset(),)*
139        }
140    };
141}
142
143/// All variant names in one place for the dispatch macros.
144macro_rules! all_variants {
145    ($mac:ident, $($args:tt)*) => {
146        $mac!($($args)*
147            Izhikevich, AdEx, ExpIF, Lapicque, HomeostaticLif,
148            HodgkinHuxley, TraubMiles, WangBuzsaki, ConnorStevens,
149            DestexheThalamic, HuberBraun, GolombFS,
150            Pospischil, MainenSejnowski, DeSchutterPurkinje,
151            PlantR15, Prescott, MihalasNiebur, GLIF, GIFPopulation,
152            AvRonCardiac, DurstewitzDopamine, HillTononi, BertramPhantom, Yamada,
153            FitzHughNagumo, MorrisLecar, HindmarshRose, ResonateAndFire,
154            FitzHughRinzel, McKean, TermanWang, GutkinErmentrout, WilsonHR,
155            Chay, ChayKeizer, ShermanRinzelKeizer, ButeraRespiratory,
156            EPropALIF, SuperSpike, LearnableNeuron, Pernarowski,
157            QuadraticIF, Theta, PerfectIntegrator, GatedLIF, NonlinearLIF,
158            SFA, MAT, KLIF, InhibitoryLIF, ComplementaryLIF, ParametricLIF,
159            NonResettingLIF, AdaptiveThresholdIF, SigmaDelta, EnergyLIF,
160            ClosedFormContinuous,
161            ChialvoMap, RulkovMap, IbarzTanakaMap, MedvedevMap,
162            CazellesMap, CourageNekorkinMap,
163            BrainScaleSAdEx, SpiNNakerLIF, NeuroGrid, DPI,
164            MarderSTG, RallCable, BoothRinzel, Dendrify,
165            LiquidTimeConstant, ParallelSpiking, FractionalLIF,
166            StochasticIF, GalvesLocherbach, SpikeResponse, GLM,
167            Arcane,
168        )
169    };
170}
171
172impl NeuronVariant {
173    pub fn step(&mut self, current: f64) -> i32 {
174        all_variants!(dispatch_step, self, current,)
175    }
176
177    pub fn reset(&mut self) {
178        all_variants!(dispatch_reset, self,)
179    }
180
181    pub fn soma_voltage(&self) -> f64 {
182        match self {
183            NeuronVariant::Izhikevich(n) => n.v,
184            NeuronVariant::AdEx(n) => n.v,
185            NeuronVariant::ExpIF(n) => n.v,
186            NeuronVariant::Lapicque(n) => n.v,
187            NeuronVariant::HomeostaticLif(n) => n.v,
188            NeuronVariant::HodgkinHuxley(n) => n.v,
189            NeuronVariant::TraubMiles(n) => n.v,
190            NeuronVariant::WangBuzsaki(n) => n.v,
191            NeuronVariant::ConnorStevens(n) => n.v,
192            NeuronVariant::DestexheThalamic(n) => n.v,
193            NeuronVariant::HuberBraun(n) => n.v,
194            NeuronVariant::GolombFS(n) => n.v,
195            NeuronVariant::Pospischil(n) => n.v,
196            NeuronVariant::MainenSejnowski(n) => n.vs,
197            NeuronVariant::DeSchutterPurkinje(n) => n.v,
198            NeuronVariant::PlantR15(n) => n.v,
199            NeuronVariant::Prescott(n) => n.v,
200            NeuronVariant::MihalasNiebur(n) => n.v,
201            NeuronVariant::GLIF(n) => n.v,
202            NeuronVariant::GIFPopulation(n) => n.v,
203            NeuronVariant::AvRonCardiac(n) => n.v,
204            NeuronVariant::DurstewitzDopamine(n) => n.v,
205            NeuronVariant::HillTononi(n) => n.v,
206            NeuronVariant::BertramPhantom(n) => n.v,
207            NeuronVariant::Yamada(n) => n.v,
208            NeuronVariant::FitzHughNagumo(n) => n.v,
209            NeuronVariant::MorrisLecar(n) => n.v,
210            NeuronVariant::HindmarshRose(n) => n.x,
211            NeuronVariant::ResonateAndFire(n) => n.x,
212            NeuronVariant::FitzHughRinzel(n) => n.v,
213            NeuronVariant::McKean(n) => n.v,
214            NeuronVariant::TermanWang(n) => n.v,
215            NeuronVariant::GutkinErmentrout(n) => n.v,
216            NeuronVariant::WilsonHR(n) => n.v,
217            NeuronVariant::Chay(n) => n.v,
218            NeuronVariant::ChayKeizer(n) => n.v,
219            NeuronVariant::ShermanRinzelKeizer(n) => n.v,
220            NeuronVariant::ButeraRespiratory(n) => n.v,
221            NeuronVariant::EPropALIF(n) => n.v,
222            NeuronVariant::SuperSpike(n) => n.v,
223            NeuronVariant::LearnableNeuron(n) => n.v,
224            NeuronVariant::Pernarowski(n) => n.v,
225            NeuronVariant::QuadraticIF(n) => n.v,
226            NeuronVariant::Theta(n) => n.theta,
227            NeuronVariant::PerfectIntegrator(n) => n.v,
228            NeuronVariant::GatedLIF(n) => n.v,
229            NeuronVariant::NonlinearLIF(n) => n.v,
230            NeuronVariant::SFA(n) => n.v,
231            NeuronVariant::MAT(n) => n.v,
232            NeuronVariant::KLIF(n) => n.v,
233            NeuronVariant::InhibitoryLIF(n) => n.v,
234            NeuronVariant::ComplementaryLIF(n) => n.v_pos,
235            NeuronVariant::ParametricLIF(n) => n.v,
236            NeuronVariant::NonResettingLIF(n) => n.v,
237            NeuronVariant::AdaptiveThresholdIF(n) => n.v,
238            NeuronVariant::SigmaDelta(n) => n.sigma,
239            NeuronVariant::EnergyLIF(n) => n.v,
240            NeuronVariant::ClosedFormContinuous(n) => n.x,
241            NeuronVariant::ChialvoMap(n) => n.x,
242            NeuronVariant::RulkovMap(n) => n.x,
243            NeuronVariant::IbarzTanakaMap(n) => n.x,
244            NeuronVariant::MedvedevMap(n) => n.x,
245            NeuronVariant::CazellesMap(n) => n.x,
246            NeuronVariant::CourageNekorkinMap(n) => n.x,
247            NeuronVariant::BrainScaleSAdEx(n) => n.v,
248            NeuronVariant::SpiNNakerLIF(n) => n.v,
249            NeuronVariant::NeuroGrid(n) => n.v_s,
250            NeuronVariant::DPI(n) => n.i_mem,
251            NeuronVariant::MarderSTG(n) => n.v,
252            NeuronVariant::RallCable(n) => n.v.first().copied().unwrap_or(0.0),
253            NeuronVariant::BoothRinzel(n) => n.vs,
254            NeuronVariant::Dendrify(n) => n.v_s,
255            NeuronVariant::LiquidTimeConstant(n) => n.x,
256            NeuronVariant::ParallelSpiking(_) => 0.0,
257            NeuronVariant::FractionalLIF(n) => n.v,
258            NeuronVariant::StochasticIF(n) => n.v,
259            NeuronVariant::GalvesLocherbach(n) => n.v,
260            NeuronVariant::SpikeResponse(n) => n.v,
261            NeuronVariant::GLM(n) => n.mu,
262            NeuronVariant::Arcane(n) => n.v_fast,
263        }
264    }
265}
266
267// ── PopulationRunner ────────────────────────────────────────────────
268
269pub struct PopulationRunner {
270    neurons: Vec<NeuronVariant>,
271    spikes: Vec<u8>,
272    currents: Vec<f64>,
273}
274
275const CHUNK_SIZE: usize = 256;
276
277impl PopulationRunner {
278    pub fn new(neurons: Vec<NeuronVariant>) -> Self {
279        let n = neurons.len();
280        Self {
281            neurons,
282            spikes: vec![0u8; n],
283            currents: vec![0.0; n],
284        }
285    }
286
287    pub fn len(&self) -> usize {
288        self.neurons.len()
289    }
290
291    pub fn is_empty(&self) -> bool {
292        self.neurons.is_empty()
293    }
294
295    pub fn step_all(&mut self) {
296        let neurons = &mut self.neurons;
297        let spikes = &mut self.spikes;
298        let currents = &self.currents;
299
300        neurons
301            .par_chunks_mut(CHUNK_SIZE)
302            .zip(spikes.par_chunks_mut(CHUNK_SIZE))
303            .zip(currents.par_chunks(CHUNK_SIZE))
304            .for_each(|((n_chunk, s_chunk), c_chunk)| {
305                for i in 0..n_chunk.len() {
306                    let fired = n_chunk[i].step(c_chunk[i]);
307                    s_chunk[i] = if fired != 0 { 1 } else { 0 };
308                }
309            });
310    }
311
312    pub fn reset_all(&mut self) {
313        for n in &mut self.neurons {
314            n.reset();
315        }
316        self.spikes.fill(0);
317        self.currents.fill(0.0);
318    }
319
320    pub fn reset_currents(&mut self) {
321        self.currents.fill(0.0);
322    }
323
324    pub fn collect_voltages(&self) -> Vec<f64> {
325        self.neurons.iter().map(|n| n.soma_voltage()).collect()
326    }
327}
328
329// ── ProjectionRunner ────────────────────────────────────────────────
330
331/// CSR-stored synaptic projection with optional axonal delay.
332pub struct ProjectionRunner {
333    pub src_pop: usize,
334    pub tgt_pop: usize,
335    row_offsets: Vec<usize>,
336    col_indices: Vec<usize>,
337    values: Vec<f64>,
338    delay_steps: usize,
339    delay_buffer: Vec<Vec<u8>>,
340    buf_idx: usize,
341}
342
343impl ProjectionRunner {
344    pub fn new(
345        src_pop: usize,
346        tgt_pop: usize,
347        row_offsets: Vec<usize>,
348        col_indices: Vec<usize>,
349        values: Vec<f64>,
350        delay_steps: usize,
351    ) -> Self {
352        let n_delay = if delay_steps > 0 { delay_steps } else { 0 };
353        let n_src = if row_offsets.is_empty() {
354            0
355        } else {
356            row_offsets.len() - 1
357        };
358        let delay_buffer = if n_delay > 0 {
359            vec![vec![0u8; n_src]; n_delay]
360        } else {
361            Vec::new()
362        };
363        Self {
364            src_pop,
365            tgt_pop,
366            row_offsets,
367            col_indices,
368            values,
369            delay_steps: n_delay,
370            delay_buffer,
371            buf_idx: 0,
372        }
373    }
374
375    /// Scatter spikes through CSR connectivity into target current buffer.
376    pub fn propagate(&mut self, src_spikes: &[u8], tgt_currents: &mut [f64]) {
377        let spikes = if self.delay_steps > 0 {
378            let delayed = &self.delay_buffer[self.buf_idx];
379            let out: Vec<u8> = delayed.clone();
380            self.delay_buffer[self.buf_idx] = src_spikes.to_vec();
381            self.buf_idx = (self.buf_idx + 1) % self.delay_steps;
382            out
383        } else {
384            src_spikes.to_vec()
385        };
386
387        let n_src = self.row_offsets.len().saturating_sub(1);
388        for i in 0..n_src {
389            if spikes.get(i).copied().unwrap_or(0) == 0 {
390                continue;
391            }
392            let start = self.row_offsets[i];
393            let end = self.row_offsets[i + 1];
394            for k in start..end {
395                let j = self.col_indices[k];
396                if j < tgt_currents.len() {
397                    tgt_currents[j] += self.values[k];
398                }
399            }
400        }
401    }
402}
403
404// ── SimResults ──────────────────────────────────────────────────────
405
406pub struct SimResults {
407    pub spike_counts: Vec<usize>,
408    /// Per-population flat spike data: (neuron_id << 32) | timestep packed as u64.
409    /// Supports up to 2^32 neurons and 2^32 timesteps.
410    pub spike_data: Vec<Vec<u64>>,
411    pub voltages: Vec<Vec<f64>>,
412}
413
414// ── NetworkRunner ───────────────────────────────────────────────────
415
416pub struct NetworkRunner {
417    pub populations: Vec<PopulationRunner>,
418    pub projections: Vec<ProjectionRunner>,
419}
420
421impl NetworkRunner {
422    pub fn new() -> Self {
423        Self {
424            populations: Vec::new(),
425            projections: Vec::new(),
426        }
427    }
428
429    pub fn add_population(&mut self, pop: PopulationRunner) -> usize {
430        let idx = self.populations.len();
431        self.populations.push(pop);
432        idx
433    }
434
435    pub fn add_projection(&mut self, proj: ProjectionRunner) {
436        self.projections.push(proj);
437    }
438
439    pub fn run(&mut self, n_steps: usize) -> SimResults {
440        let n_pops = self.populations.len();
441        let mut spike_counts = vec![0usize; n_pops];
442        let mut spike_data: Vec<Vec<u64>> = vec![Vec::new(); n_pops];
443
444        for t in 0..n_steps {
445            // Reset currents
446            for pop in &mut self.populations {
447                pop.reset_currents();
448            }
449
450            // Propagate spikes through projections
451            // Must borrow populations mutably for target currents, immutably for source spikes.
452            // Use index-based access to avoid aliasing issues.
453            for proj_idx in 0..self.projections.len() {
454                let src = self.projections[proj_idx].src_pop;
455                let tgt = self.projections[proj_idx].tgt_pop;
456                if src == tgt {
457                    // Self-projection: copy spikes, then propagate
458                    let spikes_copy = self.populations[src].spikes.clone();
459                    let currents = &mut self.populations[tgt].currents;
460                    self.projections[proj_idx].propagate(&spikes_copy, currents);
461                } else {
462                    // Split borrow via raw pointers (safe because src != tgt)
463                    let pops_ptr = self.populations.as_mut_ptr();
464                    let src_spikes = unsafe { &(*pops_ptr.add(src)).spikes };
465                    let tgt_currents = unsafe { &mut (*pops_ptr.add(tgt)).currents };
466                    self.projections[proj_idx].propagate(src_spikes, tgt_currents);
467                }
468            }
469
470            // Step all populations
471            for (pop_idx, pop) in self.populations.iter_mut().enumerate() {
472                pop.step_all();
473                for (nid, &spike) in pop.spikes.iter().enumerate() {
474                    if spike != 0 {
475                        spike_counts[pop_idx] += 1;
476                        spike_data[pop_idx].push(((nid as u64) << 32) | (t as u64));
477                    }
478                }
479            }
480        }
481
482        let voltages: Vec<Vec<f64>> = self
483            .populations
484            .iter()
485            .map(|p| p.collect_voltages())
486            .collect();
487
488        SimResults {
489            spike_counts,
490            spike_data,
491            voltages,
492        }
493    }
494}
495
496impl Default for NetworkRunner {
497    fn default() -> Self {
498        Self::new()
499    }
500}
501
502// ── Factory ─────────────────────────────────────────────────────────
503
504/// Create a population of `n` identical neurons by model name.
505pub fn create_population(model_name: &str, n: usize) -> Result<PopulationRunner, String> {
506    let neurons: Vec<NeuronVariant> = (0..n)
507        .map(|_| create_neuron(model_name))
508        .collect::<Result<_, _>>()?;
509    Ok(PopulationRunner::new(neurons))
510}
511
512pub fn create_neuron(name: &str) -> Result<NeuronVariant, String> {
513    match name {
514        "Izhikevich" => Ok(NeuronVariant::Izhikevich(Izhikevich::regular_spiking())),
515        "AdEx" | "AdExNeuron" => Ok(NeuronVariant::AdEx(AdExNeuron::new())),
516        "ExpIF" | "ExpIfNeuron" => Ok(NeuronVariant::ExpIF(ExpIfNeuron::new())),
517        "Lapicque" | "LapicqueNeuron" => Ok(NeuronVariant::Lapicque(LapicqueNeuron::new(
518            20.0, 1.0, 1.0, 1.0,
519        ))),
520        "HomeostaticLif" => Ok(NeuronVariant::HomeostaticLif(
521            HomeostaticLif::with_defaults(),
522        )),
523        "HodgkinHuxley" | "HodgkinHuxleyNeuron" => {
524            Ok(NeuronVariant::HodgkinHuxley(HodgkinHuxleyNeuron::new()))
525        }
526        "TraubMiles" | "TraubMilesNeuron" => Ok(NeuronVariant::TraubMiles(TraubMilesNeuron::new())),
527        "WangBuzsaki" | "WangBuzsakiNeuron" => {
528            Ok(NeuronVariant::WangBuzsaki(WangBuzsakiNeuron::new()))
529        }
530        "ConnorStevens" | "ConnorStevensNeuron" => {
531            Ok(NeuronVariant::ConnorStevens(ConnorStevensNeuron::new()))
532        }
533        "DestexheThalamic" | "DestexheThalamicNeuron" => Ok(NeuronVariant::DestexheThalamic(
534            DestexheThalamicNeuron::new(),
535        )),
536        "HuberBraun" | "HuberBraunNeuron" => Ok(NeuronVariant::HuberBraun(HuberBraunNeuron::new())),
537        "GolombFS" | "GolombFSNeuron" => Ok(NeuronVariant::GolombFS(GolombFSNeuron::new())),
538        "Pospischil" | "PospischilNeuron" => Ok(NeuronVariant::Pospischil(PospischilNeuron::new())),
539        "MainenSejnowski" | "MainenSejnowskiNeuron" => {
540            Ok(NeuronVariant::MainenSejnowski(MainenSejnowskiNeuron::new()))
541        }
542        "DeSchutterPurkinje" | "DeSchutterPurkinjeNeuron" => Ok(NeuronVariant::DeSchutterPurkinje(
543            DeSchutterPurkinjeNeuron::new(),
544        )),
545        "PlantR15" | "PlantR15Neuron" => Ok(NeuronVariant::PlantR15(PlantR15Neuron::new())),
546        "Prescott" | "PrescottNeuron" => Ok(NeuronVariant::Prescott(PrescottNeuron::new())),
547        "MihalasNiebur" | "MihalasNieburNeuron" => {
548            Ok(NeuronVariant::MihalasNiebur(MihalasNieburNeuron::new()))
549        }
550        "GLIF" | "GLIFNeuron" => Ok(NeuronVariant::GLIF(GLIFNeuron::new())),
551        "GIFPopulation" | "GIFPopulationNeuron" => {
552            Ok(NeuronVariant::GIFPopulation(GIFPopulationNeuron::new(42)))
553        }
554        "AvRonCardiac" | "AvRonCardiacNeuron" => {
555            Ok(NeuronVariant::AvRonCardiac(AvRonCardiacNeuron::new()))
556        }
557        "DurstewitzDopamine" | "DurstewitzDopamineNeuron" => Ok(NeuronVariant::DurstewitzDopamine(
558            DurstewitzDopamineNeuron::new(),
559        )),
560        "HillTononi" | "HillTononiNeuron" => Ok(NeuronVariant::HillTononi(HillTononiNeuron::new())),
561        "BertramPhantom" | "BertramPhantomBurster" => {
562            Ok(NeuronVariant::BertramPhantom(BertramPhantomBurster::new()))
563        }
564        "Yamada" | "YamadaNeuron" => Ok(NeuronVariant::Yamada(YamadaNeuron::new())),
565        "FitzHughNagumo" | "FitzHughNagumoNeuron" => {
566            Ok(NeuronVariant::FitzHughNagumo(FitzHughNagumoNeuron::new()))
567        }
568        "MorrisLecar" | "MorrisLecarNeuron" => {
569            Ok(NeuronVariant::MorrisLecar(MorrisLecarNeuron::new()))
570        }
571        "HindmarshRose" | "HindmarshRoseNeuron" => {
572            Ok(NeuronVariant::HindmarshRose(HindmarshRoseNeuron::new()))
573        }
574        "ResonateAndFire" | "ResonateAndFireNeuron" => {
575            Ok(NeuronVariant::ResonateAndFire(ResonateAndFireNeuron::new()))
576        }
577        "FitzHughRinzel" | "FitzHughRinzelNeuron" => {
578            Ok(NeuronVariant::FitzHughRinzel(FitzHughRinzelNeuron::new()))
579        }
580        "McKean" | "McKeanNeuron" => Ok(NeuronVariant::McKean(McKeanNeuron::new())),
581        "TermanWang" | "TermanWangOscillator" => {
582            Ok(NeuronVariant::TermanWang(TermanWangOscillator::new()))
583        }
584        "GutkinErmentrout" | "GutkinErmentroutNeuron" => Ok(NeuronVariant::GutkinErmentrout(
585            GutkinErmentroutNeuron::new(),
586        )),
587        "WilsonHR" | "WilsonHRNeuron" => Ok(NeuronVariant::WilsonHR(WilsonHRNeuron::new())),
588        "Chay" | "ChayNeuron" => Ok(NeuronVariant::Chay(ChayNeuron::new())),
589        "ChayKeizer" | "ChayKeizerNeuron" => Ok(NeuronVariant::ChayKeizer(ChayKeizerNeuron::new())),
590        "ShermanRinzelKeizer" | "ShermanRinzelKeizerNeuron" => Ok(
591            NeuronVariant::ShermanRinzelKeizer(ShermanRinzelKeizerNeuron::new()),
592        ),
593        "ButeraRespiratory" | "ButeraRespiratoryNeuron" => Ok(NeuronVariant::ButeraRespiratory(
594            ButeraRespiratoryNeuron::new(),
595        )),
596        "EPropALIF" | "EPropALIFNeuron" => Ok(NeuronVariant::EPropALIF(EPropALIFNeuron::default())),
597        "SuperSpike" | "SuperSpikeNeuron" => {
598            Ok(NeuronVariant::SuperSpike(SuperSpikeNeuron::default()))
599        }
600        "LearnableNeuron" | "LearnableNeuronModel" => {
601            Ok(NeuronVariant::LearnableNeuron(LearnableNeuronModel::new()))
602        }
603        "Pernarowski" | "PernarowskiNeuron" => {
604            Ok(NeuronVariant::Pernarowski(PernarowskiNeuron::new()))
605        }
606        "QuadraticIF" | "QuadraticIFNeuron" => {
607            Ok(NeuronVariant::QuadraticIF(QuadraticIFNeuron::default()))
608        }
609        "Theta" | "ThetaNeuron" => Ok(NeuronVariant::Theta(ThetaNeuron::default())),
610        "PerfectIntegrator" | "PerfectIntegratorNeuron" => Ok(NeuronVariant::PerfectIntegrator(
611            PerfectIntegratorNeuron::default(),
612        )),
613        "GatedLIF" | "GatedLIFNeuron" => Ok(NeuronVariant::GatedLIF(GatedLIFNeuron::default())),
614        "NonlinearLIF" | "NonlinearLIFNeuron" => {
615            Ok(NeuronVariant::NonlinearLIF(NonlinearLIFNeuron::new()))
616        }
617        "SFA" | "SFANeuron" => Ok(NeuronVariant::SFA(SFANeuron::new())),
618        "MAT" | "MATNeuron" => Ok(NeuronVariant::MAT(MATNeuron::new())),
619        "KLIF" | "KLIFNeuron" => Ok(NeuronVariant::KLIF(KLIFNeuron::default())),
620        "InhibitoryLIF" | "InhibitoryLIFNeuron" => {
621            Ok(NeuronVariant::InhibitoryLIF(InhibitoryLIFNeuron::default()))
622        }
623        "ComplementaryLIF" | "ComplementaryLIFNeuron" => Ok(NeuronVariant::ComplementaryLIF(
624            ComplementaryLIFNeuron::default(),
625        )),
626        "ParametricLIF" | "ParametricLIFNeuron" => {
627            Ok(NeuronVariant::ParametricLIF(ParametricLIFNeuron::default()))
628        }
629        "NonResettingLIF" | "NonResettingLIFNeuron" => {
630            Ok(NeuronVariant::NonResettingLIF(NonResettingLIFNeuron::new()))
631        }
632        "AdaptiveThresholdIF" | "AdaptiveThresholdIFNeuron" => Ok(
633            NeuronVariant::AdaptiveThresholdIF(AdaptiveThresholdIFNeuron::new()),
634        ),
635        "SigmaDelta" | "SigmaDeltaNeuron" => {
636            Ok(NeuronVariant::SigmaDelta(SigmaDeltaNeuron::default()))
637        }
638        "EnergyLIF" | "EnergyLIFNeuron" => Ok(NeuronVariant::EnergyLIF(EnergyLIFNeuron::new())),
639        "ClosedFormContinuous" | "ClosedFormContinuousNeuron" => Ok(
640            NeuronVariant::ClosedFormContinuous(ClosedFormContinuousNeuron::new()),
641        ),
642        "ChialvoMap" | "ChialvoMapNeuron" => Ok(NeuronVariant::ChialvoMap(ChialvoMapNeuron::new())),
643        "RulkovMap" | "RulkovMapNeuron" => Ok(NeuronVariant::RulkovMap(RulkovMapNeuron::new())),
644        "IbarzTanakaMap" | "IbarzTanakaMapNeuron" => {
645            Ok(NeuronVariant::IbarzTanakaMap(IbarzTanakaMapNeuron::new()))
646        }
647        "MedvedevMap" | "MedvedevMapNeuron" => {
648            Ok(NeuronVariant::MedvedevMap(MedvedevMapNeuron::default()))
649        }
650        "CazellesMap" | "CazellesMapNeuron" => {
651            Ok(NeuronVariant::CazellesMap(CazellesMapNeuron::new()))
652        }
653        "CourageNekorkinMap" | "CourageNekorkinMapNeuron" => Ok(NeuronVariant::CourageNekorkinMap(
654            CourageNekorkinMapNeuron::new(),
655        )),
656        "BrainScaleSAdEx" | "BrainScaleSAdExNeuron" => {
657            Ok(NeuronVariant::BrainScaleSAdEx(BrainScaleSAdExNeuron::new()))
658        }
659        "SpiNNakerLIF" | "SpiNNakerLIFNeuron" => {
660            Ok(NeuronVariant::SpiNNakerLIF(SpiNNakerLIFNeuron::new()))
661        }
662        "NeuroGrid" | "NeuroGridNeuron" => Ok(NeuronVariant::NeuroGrid(NeuroGridNeuron::new())),
663        "DPI" | "DPINeuron" => Ok(NeuronVariant::DPI(DPINeuron::new())),
664        "MarderSTG" | "MarderSTGNeuron" => Ok(NeuronVariant::MarderSTG(MarderSTGNeuron::new())),
665        "RallCable" | "RallCableNeuron" => Ok(NeuronVariant::RallCable(RallCableNeuron::new(5))),
666        "BoothRinzel" | "BoothRinzelNeuron" => {
667            Ok(NeuronVariant::BoothRinzel(BoothRinzelNeuron::new()))
668        }
669        "Dendrify" | "DendrifyNeuron" => Ok(NeuronVariant::Dendrify(DendrifyNeuron::new())),
670        "LiquidTimeConstant" | "LiquidTimeConstantNeuron" => Ok(NeuronVariant::LiquidTimeConstant(
671            LiquidTimeConstantNeuron::new(),
672        )),
673        "ParallelSpiking" | "ParallelSpikingNeuron" => Ok(NeuronVariant::ParallelSpiking(
674            ParallelSpikingNeuron::new(4, 0.5),
675        )),
676        "FractionalLIF" | "FractionalLIFNeuron" => Ok(NeuronVariant::FractionalLIF(
677            FractionalLIFNeuron::new(0.8, 50),
678        )),
679        "StochasticIF" | "StochasticIFNeuron" => {
680            Ok(NeuronVariant::StochasticIF(StochasticIFNeuron::new(42)))
681        }
682        "GalvesLocherbach" | "GalvesLocherbachNeuron" => Ok(NeuronVariant::GalvesLocherbach(
683            GalvesLocherbachNeuron::new(42),
684        )),
685        "SpikeResponse" | "SpikeResponseNeuron" => {
686            Ok(NeuronVariant::SpikeResponse(SpikeResponseNeuron::new()))
687        }
688        "GLM" | "GLMNeuron" => Ok(NeuronVariant::GLM(GLMNeuron::new(5, 10, 42))),
689        "Arcane" | "ArcaneNeuron" => Ok(NeuronVariant::Arcane(ArcaneNeuron::new())),
690        _ => Err(format!("Unsupported model: '{name}'")),
691    }
692}
693
694/// List all supported model names.
695pub fn supported_models() -> Vec<&'static str> {
696    vec![
697        "Izhikevich",
698        "AdEx",
699        "ExpIF",
700        "Lapicque",
701        "HomeostaticLif",
702        "HodgkinHuxley",
703        "TraubMiles",
704        "WangBuzsaki",
705        "ConnorStevens",
706        "DestexheThalamic",
707        "FitzHughNagumo",
708        "MorrisLecar",
709        "HindmarshRose",
710        "ResonateAndFire",
711        "FitzHughRinzel",
712        "McKean",
713        "TermanWang",
714        "GutkinErmentrout",
715        "WilsonHR",
716        "Chay",
717        "ChayKeizer",
718        "ShermanRinzelKeizer",
719        "ButeraRespiratory",
720        "EPropALIF",
721        "SuperSpike",
722        "LearnableNeuron",
723        "Pernarowski",
724        "QuadraticIF",
725        "Theta",
726        "PerfectIntegrator",
727        "GatedLIF",
728        "NonlinearLIF",
729        "SFA",
730        "MAT",
731        "KLIF",
732        "InhibitoryLIF",
733        "ComplementaryLIF",
734        "ParametricLIF",
735        "NonResettingLIF",
736        "AdaptiveThresholdIF",
737        "SigmaDelta",
738        "EnergyLIF",
739        "ClosedFormContinuous",
740        "ChialvoMap",
741        "RulkovMap",
742        "IbarzTanakaMap",
743        "MedvedevMap",
744        "CazellesMap",
745        "CourageNekorkinMap",
746        "BrainScaleSAdEx",
747        "SpiNNakerLIF",
748        "NeuroGrid",
749        "DPI",
750        "MarderSTG",
751        "RallCable",
752        "BoothRinzel",
753        "Dendrify",
754        "LiquidTimeConstant",
755        "ParallelSpiking",
756        "FractionalLIF",
757        "StochasticIF",
758        "GalvesLocherbach",
759        "SpikeResponse",
760        "GLM",
761        "ArcaneNeuron",
762    ]
763}
764
765// ── Tests ───────────────────────────────────────────────────────────
766
767#[cfg(test)]
768mod tests {
769    use super::*;
770
771    #[test]
772    fn izhikevich_population_spikes() {
773        let mut pop = create_population("Izhikevich", 10).unwrap();
774        let mut total_spikes = 0usize;
775        for _ in 0..100 {
776            pop.currents.fill(10.0);
777            pop.step_all();
778            total_spikes += pop.spikes.iter().filter(|&&s| s != 0).count();
779        }
780        assert!(
781            total_spikes > 0,
782            "10 Izhikevich neurons must spike with I=10"
783        );
784    }
785
786    #[test]
787    fn projection_propagates_spikes() {
788        let row_offsets = vec![0, 2, 4];
789        let col_indices = vec![0, 1, 0, 1];
790        let values = vec![5.0, 3.0, 2.0, 4.0];
791
792        let mut proj = ProjectionRunner::new(0, 1, row_offsets, col_indices, values, 0);
793        let src_spikes = vec![1u8, 0];
794        let mut tgt_currents = vec![0.0; 2];
795        proj.propagate(&src_spikes, &mut tgt_currents);
796
797        assert!((tgt_currents[0] - 5.0).abs() < 1e-10);
798        assert!((tgt_currents[1] - 3.0).abs() < 1e-10);
799    }
800
801    #[test]
802    fn all_to_all_network_100_steps() {
803        let mut runner = NetworkRunner::new();
804        let pop = create_population("Izhikevich", 4).unwrap();
805        runner.add_population(pop);
806
807        // All-to-all CSR: 4 src -> 4 tgt
808        let mut row_offsets = Vec::new();
809        let mut col_indices = Vec::new();
810        let mut values = Vec::new();
811        let mut offset = 0;
812        for _i in 0..4 {
813            row_offsets.push(offset);
814            for j in 0..4 {
815                col_indices.push(j);
816                values.push(2.0);
817                offset += 1;
818            }
819        }
820        row_offsets.push(offset);
821
822        let proj = ProjectionRunner::new(0, 0, row_offsets, col_indices, values, 0);
823        runner.add_projection(proj);
824
825        // Inject external current by pre-filling
826        for n in &mut runner.populations[0].neurons {
827            if let NeuronVariant::Izhikevich(iz) = n {
828                iz.v = -50.0;
829            }
830        }
831
832        let results = runner.run(100);
833        assert_eq!(results.spike_counts.len(), 1);
834        assert_eq!(results.voltages.len(), 1);
835        assert_eq!(results.voltages[0].len(), 4);
836    }
837
838    #[test]
839    fn mixed_hh_adex_network() {
840        let mut runner = NetworkRunner::new();
841
842        let hh_pop = create_population("HodgkinHuxley", 3).unwrap();
843        let adex_pop = create_population("AdEx", 3).unwrap();
844        let hh_idx = runner.add_population(hh_pop);
845        let adex_idx = runner.add_population(adex_pop);
846
847        // HH -> AdEx projection
848        let row_offsets = vec![0, 3, 6, 9];
849        let col_indices = vec![0, 1, 2, 0, 1, 2, 0, 1, 2];
850        let values = vec![100.0; 9];
851        let proj = ProjectionRunner::new(hh_idx, adex_idx, row_offsets, col_indices, values, 0);
852        runner.add_projection(proj);
853
854        // Drive HH with external current
855        runner.populations[0].currents.fill(15.0);
856
857        let results = runner.run(50);
858        assert_eq!(results.spike_counts.len(), 2);
859        assert_eq!(results.voltages.len(), 2);
860    }
861
862    #[test]
863    fn large_network_performance() {
864        let n = 1000;
865        let mut pop = create_population("Izhikevich", n).unwrap();
866        // Run 1000 steps with constant drive — should complete quickly
867        for _ in 0..1000 {
868            pop.currents.fill(10.0);
869            pop.step_all();
870        }
871        let total: usize = pop.spikes.iter().map(|&s| s as usize).sum();
872        // Sanity: spike count should be deterministic and nonzero after 1000 driven steps
873        let _ = total;
874        // Check voltages are finite
875        let voltages = pop.collect_voltages();
876        assert_eq!(voltages.len(), n);
877        for v in &voltages {
878            assert!(v.is_finite(), "voltage must be finite");
879        }
880    }
881
882    #[test]
883    fn batch_simulate_single_neuron() {
884        let mut neuron = create_neuron("AdEx").unwrap();
885        let n_steps = 1000;
886        let current = 500.0;
887        let mut voltages = Vec::with_capacity(n_steps);
888        let mut spikes = Vec::new();
889        for t in 0..n_steps {
890            let fired = neuron.step(current);
891            voltages.push(neuron.soma_voltage());
892            if fired != 0 {
893                spikes.push(t);
894            }
895        }
896        assert_eq!(voltages.len(), n_steps);
897        assert!(voltages.iter().all(|v| v.is_finite()));
898        assert!(!spikes.is_empty(), "AdEx with I=10 should spike");
899    }
900
901    #[test]
902    fn create_neuron_all_supported() {
903        for name in supported_models() {
904            let result = create_neuron(name);
905            assert!(result.is_ok(), "create_neuron({name}) failed: {:?}", result.err());
906        }
907    }
908}