Case Study 1 — Circuit Simulation: One Factorization, Ten Thousand Timesteps

Every time you open a circuit simulator — SPICE, the engine behind essentially all integrated-circuit design — and watch a waveform unfold, you are watching linear algebra run in a loop. Underneath the animated voltages is a single matrix being solved against a fresh right-hand side at every instant of simulated time. If that matrix were re-eliminated from scratch at each timestep, simulating a complex chip for a microsecond could take hours. The reason it takes seconds instead is precisely the idea of this chapter: factor the matrix once, reuse the factors for every timestep. This case study shows exactly how, with a small circuit worked by hand and then scaled up in code.

A circuit becomes a linear system

When you analyze a resistor network by nodal analysis, you write one equation per node expressing Kirchhoff's current law — the currents flowing into a node must sum to the currents flowing out. Collecting these equations produces a linear system

$$G\mathbf{v} = \mathbf{i},$$

where $\mathbf{v}$ is the vector of unknown node voltages, $\mathbf{i}$ is the vector of currents injected into the nodes by sources, and $G$ is the conductance matrix built from the resistor values (conductance is $1/\text{resistance}$). The matrix $G$ depends only on how the circuit is wired and what its components are — not on the inputs. The vector $\mathbf{i}$ is the input: change the source, change $\mathbf{i}$, but $G$ stays put. That separation — fixed structure $G$, varying input $\mathbf{i}$ — is the whole reason LU is the perfect tool here.

Take a small three-node circuit whose conductance matrix is

$$G = \begin{bmatrix} 2 & -1 & 0 \\ -1 & 3 & -1 \\ 0 & -1 & 2 \end{bmatrix}.$$

(The diagonal entries are the total conductance attached to each node; the off-diagonals are the negated conductances of the resistors linking pairs of nodes — node 1 connects to node 2, node 2 to node 3, with grounding resistors at the ends. The exact wiring does not matter for our purposes; what matters is that $G$ is a fixed matrix we will solve repeatedly.) Suppose first that we inject $3$ amperes of current into node 1 and nothing elsewhere, so $\mathbf{i}_1 = (3, 0, 0)$.

Factor the conductance matrix once

We factor $G = LU$ by the elimination of §10.3 — no row swaps are needed because $G$ is symmetric and diagonally dominant, so every pivot is safely nonzero (this is typical of conductance matrices, which is why circuit solvers can often skip pivoting or use a cheap symmetric variant).

Column 1. Pivot $2$. The only entry to clear below it is the $-1$ in row 2; the multiplier is $\ell_{21} = -1/2 = -0.5$. Subtracting $-0.5$ times row 1 from row 2 (i.e. adding half of row 1) turns row 2 into $(-1,3,-1) + 0.5(2,-1,0) = (0,\,2.5,\,-1)$. Row 3 already has a $0$ in column 1, so $\ell_{31}=0$ and it is untouched:

$$\begin{bmatrix} 2 & -1 & 0 \\ 0 & 2.5 & -1 \\ 0 & -1 & 2 \end{bmatrix}.$$

Column 2. Pivot $2.5$. Clear the $-1$ in row 3: the multiplier is $\ell_{32} = -1/2.5 = -0.4$. Row 3 becomes $(0,-1,2) + 0.4(0,2.5,-1) = (0,0,1.6)$:

$$U = \begin{bmatrix} 2 & -1 & 0 \\ 0 & 2.5 & -1 \\ 0 & 0 & 1.6 \end{bmatrix}, \qquad L = \begin{bmatrix} 1 & 0 & 0 \\ -0.5 & 1 & 0 \\ 0 & -0.4 & 1 \end{bmatrix}.$$

That single factorization — the expensive step — is now done, and we will never repeat it for this circuit.

Each source is one cheap pair of substitutions

To find the node voltages under $\mathbf{i}_1 = (3,0,0)$, we do the two triangular solves of §10.6.

Forward-substitute $L\mathbf{y} = \mathbf{i}_1$. Reading $L$ row by row: $y_1 = 3$; then $-0.5y_1 + y_2 = 0 \Rightarrow y_2 = 1.5$; then $-0.4y_2 + y_3 = 0 \Rightarrow y_3 = 0.6$. So $\mathbf{y} = (3, 1.5, 0.6)$.

Back-substitute $U\mathbf{x} = \mathbf{y}$. From the bottom: $1.6\,v_3 = 0.6 \Rightarrow v_3 = 0.375$; then $2.5\,v_2 - v_3 = 1.5 \Rightarrow v_2 = (1.5+0.375)/2.5 = 0.75$; then $2v_1 - v_2 = 3 \Rightarrow v_1 = (3+0.75)/2 = 1.875$. The node voltages are $\mathbf{v} = (1.875, 0.75, 0.375)$ volts.

# Verify the hand factorization and solve, then reuse the factors for new sources.
import numpy as np
from scipy.linalg import lu_factor, lu_solve
G = np.array([[ 2., -1.,  0.],
              [-1.,  3., -1.],
              [ 0., -1.,  2.]])
L = np.array([[1., 0., 0.], [-0.5, 1., 0.], [0., -0.4, 1.]])
U = np.array([[2., -1., 0.], [0., 2.5, -1.], [0., 0., 1.6]])
print(np.allclose(L @ U, G))               # True -> hand factorization exact
lu_piv = lu_factor(G)                       # factor ONCE
for current in [np.array([3., 0., 0.]),     # inject at node 1
                np.array([0., 4., 0.]),     # inject at node 2
                np.array([1., -2., 1.])]:   # a mixed source
    v = lu_solve(lu_piv, current)           # reuse: only two substitutions each
    print(np.round(v, 4), np.allclose(G @ v, current))
# [1.875 0.75  0.375] True
# [1.    2.    1.   ] True
# [ 0.25 -0.5   0.25] True

The first solve reproduces our hand answer $(1.875, 0.75, 0.375)$ exactly, and the second and third sources are solved by reusing the same factorizationlu_factor(G) is called once, before the loop, and each lu_solve inside the loop is just the cheap forward-and-back substitution. We never re-factor $G$. Three sources here; in a real simulation it would be ten thousand.

Why this scales: the simulation loop

A transient circuit simulation discretizes time into thousands of tiny steps. At each step, the simulator solves a system whose matrix is built from the circuit's conductances and capacitances — and crucially, for a linear circuit with a fixed timestep, that matrix is the same at every step. Only the right-hand side changes, because it carries the previous step's state and the current source values. So the simulator's inner loop is exactly our pattern:

Factor once, before the loop: compute the LU (or, exploiting symmetry, the Cholesky) factorization of the system matrix. Then, every timestep: assemble the new right-hand side from the present state, and solve by two triangular substitutions using the stored factors.

The payoff is dramatic. Factoring an $n\times n$ matrix costs $\sim\tfrac23 n^3$; each triangular solve costs $\sim n^2$. If a simulation runs $T$ timesteps, re-factoring each step would cost $\sim T\cdot\tfrac23 n^3$, while factor-once-reuse costs $\sim \tfrac23 n^3 + T\cdot 2n^2$. For $T = 10{,}000$ timesteps and a circuit with $n = 1000$ nodes, that is the difference between ten thousand cubic-cost jobs and one cubic-cost job followed by ten thousand cheap quadratic ones — easily a thousandfold speedup. This is not a micro-optimization; it is what makes transient simulation feasible at all.

# A toy "simulation": one factorization, many timesteps with evolving right-hand sides.
import numpy as np
from scipy.linalg import lu_factor, lu_solve
rng = np.random.default_rng(1)
n = 200
M = rng.standard_normal((n, n))
G = M @ M.T + n * np.eye(n)         # a fixed SPD "conductance-like" matrix
lu_piv = lu_factor(G)               # FACTOR ONCE (the expensive O(n^3) step)
v = np.zeros(n)
for step in range(10_000):          # 10,000 timesteps
    i = rng.standard_normal(n) + 0.1 * v   # rhs depends on the evolving state v
    v = lu_solve(lu_piv, i)         # REUSE: only O(n^2) per step
print("final-state norm:", round(float(np.linalg.norm(v)), 4))   # finishes fast

This loop runs ten thousand solves in a fraction of a second, because the only $O(n^3)$ work — the lu_factor — happens once, outside the loop. Move that line inside the loop (re-factoring every step) and the same computation would crawl. That single placement decision is the entire lesson of this case study, and it is the Common Pitfall from §10.6 in its natural habitat.

What you should take away

A circuit is a fixed structure stimulated by changing inputs, and that is the exact shape of problem LU was built for. The conductance matrix $G$ is the asset; its LU factorization is how you make the asset reusable; and every timestep of a simulation is one cheap query against that reusable factorization. The hand-worked three-node example and the ten-thousand-step toy loop are the same idea at two scales — and the real SPICE simulators that design the chips in your phone are the same idea at a third. Factor once, reuse forever.

The deeper current here connects to the chapter's framing: factoring a matrix exposes structure you can reuse. In a circuit, the structure $L$ and $U$ encode is literally the circuit's elimination order, frozen so it never has to be recomputed. When you meet the symmetric positive-definite matrices of Chapter 28, you will learn that conductance matrices have an even cheaper factorization (Cholesky, $G = R^{\mathsf{T}}R$, roughly half the work of LU), which is what production circuit solvers actually use — but the principle is identical to the one you just practiced. The matrix is solved once into factors; the factors serve every input that ever comes.

Forward references: symmetric positive-definite matrices and the Cholesky factorization (Chapter 28); the condition number and when even a good factorization strains (Chapter 38). Backward: nodal analysis as a linear system was previewed with Kirchhoff's laws in Chapter 4's circuit application.