PaperlayerView exact revision on arXiv

Source-linked AI summary

Fast simulation of nonlinear deep resistive networks for energy-based computation

Filip Osana, Julie Grollier, Damien Querlioz

arXiv:2609.07356v1cs.ET

TL;DR

Repeated steady-state circuit solves make nonlinear deep resistive networks costly to simulate and train, while prior coordinate-descent speedups covered only ideal-diode models. This work extends coordinate descent to realistic monotone device characteristics and matches SPICE voltages with substantially lower runtime across the reported Digits networks.

  • Problem

    Repeated steady-state solves make nonlinear DRNs costly to study, and prior fast coordinate-descent simulation was limited to ideal-diode models rather than realistic device characteristics.

    Source: Introduction · p. 2

  • Method

    The solver formulates monotone nonlinear DRNs as convex content minimization and uses local scalar node updates for realistic device characteristics.

    Source: Introduction · pp. 2–3, Discussion · p. 17

  • Results

    For 90% of validation samples, CD reproduces SPICE steady-state voltages with relative L1 errors below 10^-4 while reducing runtime by up to 1750×.

    Source: Discussion · p. 17

  • Takeaways & Limitations

    The results demonstrate training and SPICE validation of a realistic nonlinear DRN at a scale where direct SPICE training is prohibitive.

    Source: Discussion · p. 17

    The formulation assumes linear conductances, ideal sources, and one-port monotone nonlinearities.

    Source: Discussion · p. 17

Abstract

from arXiv · show

Deep resistive networks are electronic energy-based systems in which computation is performed by the steady-state voltages of nonlinear circuits. Nonlinear devices enable expressive input-output transformations, but make the circuit equilibria costly to compute during simulation and training. Recent coordinate-descent solvers have achieved large speedups over SPICE-class circuit simulators, but only for nonlinearities modeled as ideal-diode models. This mathematical simplification is not sufficient for practical analog circuits. Here we extend coordinate-descent simulation to realistic monotone nonlinearities, including Shockley diodes, antiparallel diode pairs, and piecewise-linear current-voltage characteristics. With neighboring voltages fixed, each node update remains a scalar Kirchhoff-law solve. Single-exponential characteristics admit closed-form Lambert-(W) updates, while more general monotone characteristics can be handled with scalar root-finding methods. Across networks with one to three hidden layers and hidden widths from 64 to 1024, the solver reproduces matched SPICE steady-state voltages with relative $L_1$ errors below $1.1\times10^{-4}$ for 90% of validation samples, while achieving speedups up to $(1.7\times10^3)$. We further train a $1568\times100\times20$ double-Shockley network on MNIST, reaching about 3.0% test error and reducing the per-epoch training time by roughly $4.4\times10^2$. Together, these results establish a practical route to the circuit-level design and training of large-scale analog energy-based systems incorporating realistic nonlinear devices.

Introduction

Deep resistive networks offer a trainable electronic route to energy-based computation, but repeated steady-state solves make them costly to simulate and train. This work extends coordinate descent beyond ideal-diode models by exploiting monotone realistic device characteristics to preserve local node updates.

Introduction

  • Physical energy-based machines compute through states evolving under physical laws toward equilibria or energy-like minima.In electronic implementations, Kirchhoff’s laws, device nonlinearities, and dissipation perform energy minimization intrinsically.Source: Introduction · p. 1
  • DRNs encode trainable parameters as interlayer conductances, and their steady-state voltages minimize an electrical content function compatible with equilibrium propagation.Learning uses local differences between free and nudged equilibria.Source: Introduction · p. 2
  • SPICE’s repeated generic nonlinear solves make DRNs costly to study, while training can require millions of steady-state solves.This cost has limited SPICE-based nonlinear DRN studies to modest sizes, and scalable hardware demonstrations have often focused on linear networks.Source: Introduction · p. 2
  • Earlier coordinate-descent speedups applied to ideal-diode models, whose quadratic-program formulation simplifies updates but does not capture practical nonlinear device behavior.The paper targets realistic characteristics such as Shockley curves, antiparallel diode pairs, and transistor-derived responses.Source: Introduction · p. 2
  • Monotone realistic device characteristics preserve convexity and reduce each coordinate update to a scalar Kirchhoff-law solve, retaining layer-wise coordinate descent.With neighboring voltages fixed, the update enforces Kirchhoff’s current law at one node.Source: Introduction · p. 2
  • The solver formulates monotone nonlinear DRNs as convex content minimization and derives Lambert-W or scalar updates for analytical and tabulated device curves.It is validated against matched SPICE simulations across Digits networks and used to train a double-Shockley network on MNIST.Source: Introduction · p. 3

Results

The results evaluate coordinate descent for monotone nonlinear networks, reporting voltage agreement with SPICE across Digits configurations and matched configurations where the solver is faster. They also examine its use for MNIST-scale training and its agreement near the MNIST data distribution.

Resistive networks with nonlinear one-port elements as energy-based machines

  • Monotone nonlinear current laws make their energy terms convex, stationarity on free nodes is equivalent to KCL, and strict convexity gives a unique equilibrium.The full network content is convex on the feasible node potentials under the stated conductance and monotonicity conditions.Source: Resistive networks with nonlinear… · p. 4

Coordinate descent reduces each update to a local KCL solve

  • Each exact coordinate-descent update solves the selected node’s KCL equation with neighboring voltages fixed, although other free-node residuals generally remain nonzero until convergence.At convergence, all free-node KCL equations are satisfied.Source: Coordinate descent reduces each update… · p. 6

Device-specific scalar updates

  • Shockley coordinate updates have a closed-form Lambert-W solution, while the antiparallel double-Shockley update uses an active-branch approximation because the full law has two exponential terms.The approximation retains only the forward-biased branch; the active branch is selected using one-sided residual signs.Source: Device-specific scalar updates · pp. 7–8
  • For tabulated device laws, the solver represents the characteristic as a monotone piecewise-linear interpolant and solves each scalar update with Newton–Raphson using the active segment’s slope.The active segment is updated after each step, and iterations continue until convergence.Source: Device-specific scalar updates · p. 9

Application to deep resistive networks

  • Deep resistive networks use the same scalar coordinate updates, with independent nodes in a layer updated in parallel; bidirectional amplification allows forward voltage gain A and backward current gain B to be chosen independently.Odd- and even-indexed layers can be updated in alternating half-iterations when there are no within-layer conductive couplings.Source: Application to deep resistive networks · p. 10
  • Final P90 relative voltage error was of order 10^-5 across reported depths and nonlinearities; Shockley networks reached their error floors within 16 sweeps and PWL networks within 32.Deeper networks required more sweeps to reach the same accuracy; the reported sweep study used width-128 Digits networks.Source: Application to deep resistive networks · p. 12

Runtime scaling

  • Runtime experiments use adaptive stopping based on relative voltage change because the sweep count required for convergence depends on network conductances, nonlinearity, width, and depth.The supplied passage describes the stopping-rule motivation but does not report runtime results.Source: Runtime scaling · p. 12
  • At vtol = 10−5, P90 error is already on the order of 10−5, while tighter tolerances require more sweeps for marginal further improvement.The study uses vtol = 10−5 for its main runtime experiments as an accuracy–cost trade-off.Source: Runtime scaling · p. 12, Figure 4 · p. 14
  • CD runtime scales more moderately than SPICE as width and depth increase, with its cost set by outer sweeps and local nonlinear updates.For exponential nonlinearities, local updates use a closed-form Lambert-W expression, with clipping for extreme arguments and Newton-Raphson refinement when needed.Source: Runtime scaling · p. 15

MNIST-scale training

  • MNIST training uses four CD sweeps per equilibrium solve; increasing the budget to six or eight sweeps did not meaningfully improve final test error.A fixed sweep budget is more convenient for batched training than adaptive stopping.Source: MNIST-scale training · p. 15
  • Around 3% final test error was reached in roughly 3.5 hours; reported about 3.4% after 10 epochs and roughly one week for the same network class.The CD-trained model's learned weights were exported to SPICE, which recovered the same accuracy.Source: MNIST-scale training · p. 15

Runtime scaling

  • Speedups range from 1.2× to 1750×, while P90 relative L1 error stays below 1.0 × 10−4 and is typically in the 10−5 range.CD is faster than SPICE for every configuration with a reliable matched SPICE timing.Source: Runtime scaling · p. 15, Figure 5 · p. 16

MNIST-scale training

  • On a 30 × 30 PCA-plane input sweep, CD–SPICE discrepancies are smallest near the center and increase mainly at the edges and corners, indicating larger deviations mostly out of distribution.The comparison measures relative L1 voltage error separately for hidden and output layers; inputs near the center remain close to the MNIST data distribution.Source: MNIST-scale training · pp. 15, 17

Discussion

The method extends coordinate-descent simulation to realistic monotone nonlinear devices while retaining local node updates and the layer-parallel odd-even structure. It reproduces matched SPICE voltages at Digits scale and supports MNIST training, but remains bounded by idealizations in the circuit model.

Discussion

  • For grounded voltage-controlled devices with monotone current-voltage characteristics, convex electrical content preserves local scalar node updates and the efficient layer-parallel odd-even coordinate-descent structure.With neighboring voltages fixed, each update enforces Kirchhoff’s current law at one node.Source: Discussion · p. 17
  • Digits simulations reproduce SPICE steady-state voltages with relative L1 errors below 10^-4 for 90% of validation samples, while reducing runtime by up to 1750×.The MNIST demonstration trains a 1568 × 100 × 20 nonlinear network for 100 epochs and exports it to SPICE.Source: Discussion · p. 17
  • The formulation assumes linear conductances, ideal sources, and one-port monotone nonlinearities, leaving nonlinear conductances, non-ideal sources, and non-grounded two-terminal devices outside its stated scope.The authors also identify deeper architectures and harder inference tasks as future tests.Source: Discussion · p. 17

Methods

The methods formulate coordinate descent as energy minimization with local node updates, and incorporate bidirectional amplification into the energy and solver framework. Simulation details specify the hardware and SPICE setup used for reported experiments.

Simulation details

  • The fixed-sweep and fixed-voltage-tolerance experiments, along with MNIST experiments, used an NVIDIA RTX 3090 GPU; Digits runtime comparisons used a dual-socket Intel Xeon Gold 5220 server.CD used batch size 1 for the SPICE comparison, and reference simulations used ELDO 17.1 interactive DC mode with default convergence parameters.Source: Simulation details · pp. 17, 19

CD equations for deep resistive networks

  • The DRN energy is expressed layerwise, with nonlinear-device dissipation potentials and nudging-source terms incorporated so that stationarity yields the coordinate-update conditions.The update condition is specified for hidden layers.Source: CD equations for deep resistive networks · pp. 19–20

Bidirectional amplification

  • Bidirectional amplification combines a VCVS forward-voltage gain A and a CCCS backward-current gain B; the method allows A and B to be chosen independently and incorporates the rescaling into the energy and coordinate-descent framework.The earlier choice in ref. used B = 1/A.Source: Bidirectional amplification · p. 20, Figure 7 · p. 21

S1 Supplementary Material - Derivations

The supplementary derivations reduce single- and double-Shockley coordinate updates to single-exponential equations solvable with the principal Lambert-W branch. For the antiparallel pair, branch selection follows from the derivatives and yields a unique global minimizer.

S1.1.1 Single Shockley diode

  • For a single Shockley diode, the stationarity condition is rearranged into a canonical single-exponential form and solved for the node voltage using the principal Lambert-W branch.The derivation identifies the principal branch W0 through W0(x)e^W0(x) = x.Source: S1.1.1 Single Shockley diode · p. 23

S1.2 Piecewise Lambert-W update for the double Shockley (antiparallel) diode

  • For the double Shockley diode, assuming only one diode conducts at a time, each open branch has a strictly convex subproblem and its update reduces to a single-exponential equation solvable by Lambert-W.The upward derivative jump at zero makes the global minimizer unique, with its location determined by the one-sided derivative signs.Source: S1.2 Piecewise Lambert-W update for… · p. 24

S2 Supplementary material - Solver Details

The supplementary material explains how overrelaxation, clipped Lambert-W updates, and Newton refinement affect solver runtime, and analyzes local convergence of odd–even sweeps. It also reports timing and accuracy data, including fewer sweeps with overrelaxation to reach specified SPICE-relative error thresholds.

S2.1 Overrelaxed diode solver

  • At a P90 error threshold of 10^-4, the first qualifying sweep budget falls from 64 at ω = 0.8 to 12 at ω = 1.4.At a 2 × 10^-5 threshold, the corresponding budgets are 64, 32, 32, and 16 for ω = 0.8, 1.0, 1.2, and 1.4; each overrelaxed sweep costs nearly the same as a standard sweep.Source: S2.1 Overrelaxed diode solver · p. 25

S2.2 CD runtime decomposition

  • Total coordinate-descent runtime combines outer odd–even sweeps needed to meet the stopping criterion with inner nonlinear-update work within each sweep.The outer count reflects global convergence; inner work is reported as Newton iterations after a clipped Lambert initialization.Source: S2.2 CD runtime decomposition · p. 26

S2.2.1 Clipped Lambert solve and Newton polish

  • Clipping the Lambert exponent at 165 avoids overflow and instability, but perturbs the exact update, so a short Newton polish refines the initialization.The reported inner iterations are the Newton steps used.Source: S2.2.1 Clipped Lambert solve and… · p. 27

S2.2.2 Odd–even sweep map and local convergence rate

  • When the odd–even sweep map has spectral radius below 1, its fixed point is locally linearly stable, and a smaller radius implies faster asymptotic error decay.The average-Hessian proxy 1/(−log ρ(T)) is used to predict the average number of outer sweeps.Source: S2.2.2 Odd–even sweep map and local… · pp. 29–30
  • The average-Hessian sweep analysis compares its contraction proxy with average outer-sweep counts across one-, two-, and three-hidden-layer DRNs.The Hessian combines linear-conductance and grounded-nonlinearity contributions; the linear current-source term does not contribute.Source: S2.2.2 Odd–even sweep map and local… · pp. 29–30

S2.3 Timing summary tables

  • Table S1 reports runtime, outer-sweep, and inner-iteration data, while Table S2 reports CD–SPICE errors; configurations without generated SPICE netlists are marked with dashes.For PWL networks with three hidden layers, the reported p90 relative errors are 1.47e-05 at width 128 and 1.19e-05 at width 256.Source: S2.3 Timing summary tables · pp. 32, 33, 36, Table S1 · p. 33, Table S2 · p. 35
Loading 2609.07356v1…