Architecture

spice-rs is a faithful port of ngspice's core simulation engine in Rust. Every algorithm -- Newton-Raphson convergence, device model evaluation, timestep control, sparse factorization -- is translated directly from the ngspice C source code, preserving the same logic, the same control flow, and in many cases the same variable names.

The architecture mirrors ngspice's structure:

The solver layer is decoupled: sparse-rs provides two independent backends (KLU and Markowitz), both ported from their respective C reference implementations. spice-rs uses the Markowitz backend by default, matching ngspice's default solver.

Key design decisions

One load() method per device. ngspice uses a single DEVload function pointer per device type that handles all modes (DC, transient, AC initialization). spice-rs follows this pattern exactly. Devices check mode flags internally to determine behavior.

Persistent MNA matrix. The matrix is created once, elements are allocated during setup, and values are cleared/restamped each NR iteration. The matrix is never moved, copied, or rebuilt. This matches ngspice's architecture and is critical for Markowitz solver correctness.

State vectors with history. Device state (charges, fluxes, junction voltages) is stored in flat arrays with 8 history levels, matching ngspice's CKTstates[0..7]. The arrays are rotated between timesteps using O(1) pointer swaps.

Crate Structure

The simulation system is split across five crates in the sim/ directory.

spice-rs (sim/spice-rs/)

The core SPICE engine. Contains:

sparse-rs (sim/sparse-rs/)

Pure Rust sparse direct solver. Two independent backends:

ngspice-ffi (sim/ngspice-ffi/)

Safe Rust wrapper around libngspice (shared library). Used exclusively by the eval harness -- never by the simulation engine itself.

spice-eval (sim/spice-eval/)

Validation harness that runs test circuits through both spice-rs and ngspice (via ngspice-ffi), then compares results.

spice-rs-wasm (sim/spice-rs-wasm/)

wasm-bindgen bindings exposing spice-rs to the browser.

Dependency graph

spice-eval
|-- spice-rs
|   +-- sparse-rs
+-- ngspice-ffi
      +-- libngspice (C, linked at build time)

spice-rs-wasm
|-- spice-rs
|   +-- sparse-rs
|-- ferrite-schematic-render
+-- ferrite-data-model

Data Flow

The simulation pipeline from netlist text to results follows ngspice's structure: parse, setup, temperature, analyze, extract.

                    +-------------------+
  netlist text ---> |  parse_netlist()  | ---> ParseResult
                    +-------------------+
                              |
                    +-------------------+
                    |  circuit.setup()  |  allocate state vectors
                    +-------------------+
                              |
                    +--------------------------------+
                    |  resolve_coupled_inductors()   |  link K <-> L
                    +--------------------------------+
                              |
                    +------------------------+
                    |  circuit.temperature()  |  temp-dependent params
                    +------------------------+
                              |
                    +------------------------------+
                    |  analysis engine              |
                    |  (dc_operating_point /        |
                    |   transient / ac_analysis /   |
                    |   dc_sweep / ...)             |
                    +------------------------------+
                              |
                    +------------------------------+
                    |  extract_node_values()       |  solution -> HashMap
                    +------------------------------+

Step-by-step

1. Parse: The parser reads each line and builds the circuit topology. Device lines create device instances and allocate nodes. .MODEL lines populate model parameter structs. Analysis directives set the Analysis enum.

2. Setup: Calls device.setup(&mut states) on every device. Each device allocates contiguous state vector slots. After all devices have allocated, states.finalize() resizes all 8 history arrays.

3. Coupled inductors: Links K elements to their referenced L elements by name lookup.

4. Temperature: Calls device.temperature(temp, tnom) on every device. Devices compute temperature-dependent parameters.

5. Analysis engine: Creates a SimState, calls device.setup_matrix(&mut mna) to pre-allocate matrix elements, then runs the appropriate analysis.

6. Extract: Maps equation numbers back to node names.

The NR iteration loop

ni_iter() is the innermost loop and the performance-critical path. Each iteration:

  1. mna.clear() -- zero all matrix elements and RHS.
  2. device.pre_load() -- inductor flux computation (first pass).
  3. device.load() -- stamp conductances and currents.
  4. Add diag_gmin to matrix diagonal (if nonzero).
  5. mna.solve() -- Markowitz LU factorization + forward/backward substitution.
  6. Convergence check.
  7. If converged and NEWCONV: run device.conv_test() for per-device checks.
  8. Swap rhs and rhs_old for the next iteration.

Device Trait

All circuit components implement the Device trait, matching ngspice's SPICEdev function pointer table.

pub trait Device: std::fmt::Debug + Any {
    fn name(&self) -> &str;

    // --- Setup phase ---
    fn setup(&mut self, states: &mut StateVectors) -> usize { 0 }
    fn setup_matrix(&mut self, mna: &mut MnaSystem) {}
    fn setic(&mut self, rhs: &[f64]) {}
    fn temperature(&mut self, temp: f64, tnom: f64) {}

    // --- NR iteration ---
    fn pre_load(&mut self, mna: &mut MnaSystem,
                states: &mut StateVectors, mode: Mode) {}
    fn load(&mut self, mna: &mut MnaSystem,
            states: &mut StateVectors, mode: Mode,
            src_fact: f64, gmin: f64,
            noncon: &mut bool) -> Result<(), SimError>;
    fn conv_test(&self, mna: &MnaSystem,
                 states: &StateVectors,
                 reltol: f64, abstol: f64) -> bool { true }

    // --- Transient ---
    fn truncate(&self, states: &StateVectors) -> f64 { f64::INFINITY }
    fn accept(&mut self, states: &StateVectors) {}

    // --- AC ---
    fn ac_load(&mut self, mna: &mut MnaSystem,
               states: &StateVectors,
               omega: f64) -> Result<(), SimError> { Ok(()) }

    // --- Pole-Zero ---
    fn pz_load(&mut self, mna: &mut MnaSystem,
               s_re: f64, s_im: f64) -> Result<(), SimError> { Ok(()) }
}

Device inventory

Prefix Device Source States
R Resistor resistor.rs 0
C Capacitor capacitor.rs 2
L Inductor inductor.rs 2
K Mutual Inductor mutual_inductor.rs 0
V Voltage Source vsource.rs 0
I Current Source isource.rs 0
D Diode diode.rs 5
M (Level 1) MOSFET Level 1 mosfet1.rs 17
M (Level 2) MOSFET Level 2 mosfet2.rs 17
M (Level 3) MOSFET Level 3 mosfet3.rs 17
M (BSIM3) BSIM3v3 bsim3.rs 17
M (BSIM4) BSIM4 bsim4.rs 17
Q BJT bjt.rs 13
J JFET jfet.rs 7
E VCVS vcvs.rs 0
G VCCS vccs.rs 0
F CCCS cccs.rs 0
H CCVS ccvs.rs 0
T Transmission Line tline.rs 6

State Management

Device state -- capacitor charges, inductor fluxes, junction voltages, terminal currents -- is stored in a flat arena with multiple history levels. This matches ngspice's CKTstates[0..7] arrays.

StateVectors

pub struct StateVectors {
    states: [Vec<f64>; 8],
    num_states: usize,
}

Eight parallel arrays, each of length num_states. Indexed by (level, offset):

Allocation

During circuit.setup(), each device calls states.allocate(count) to claim a contiguous block of slots. The allocator is a simple bump allocator.

History rotation

Between accepted transient timesteps, the state vectors are rotated:

states.rotate(max_order);

This is O(1) -- it swaps Vec ownership without copying data.

Numerical integration

The integration system converts stored charges to currents using the companion model. ni_integrate() (port of ngspice NIintegrate) takes the ag coefficients, a charge state offset, and the capacitance value, and returns the companion model :

For trapezoidal:

Charges and currents are stored at adjacent offsets: qcap = charge, qcap + 1 = current (derivative of charge).

Truncation error

ckt_terr() (port of ngspice CKTterr) estimates the local truncation error for a charge state and returns the maximum safe timestep. The transient engine takes the minimum across all devices.