"""Stage C: SPICE-in-the-loop training of a toy 2-layer net (8 -> 8 tanh -> 4 softmax). Every MVM — forward AND transpose (the error transport J^T read) — is solved by ngspice on a differential resistive column bank with: per-bit ladder mismatch (fixed per device), per-column input offsets, 8-bit code quantization of a digital fp32 master (the T64 word-streaming story), and the v1.1-measured read noise (139 uV) added at each read. Three arms: SPICE loop | ideal numpy | behavioral numpy (same non-idealities, algebraic). Curve overlap SPICE ~= behavioral certifies the fidelity ladder used at 72M+. """ import os, time os.environ.setdefault('NGSPICE_LIBRARY_PATH', '/home/yurenh2/miniconda3/lib/libngspice.so') import numpy as np import matplotlib matplotlib.use('Agg') import matplotlib.pyplot as plt from PySpice.Spice.Netlist import Circuit rng = np.random.default_rng(3) D_IN, D_H, D_OUT = 8, 8, 4 R0 = 11e3 RF = 1e3 VOS_SIG = 500e-6 READ_NOISE = 139e-6 SIG_BIT = 0.004 STEPS = 300 BATCH = 8 LR = 0.15 BETA_SCALE = 1.0 # top nudge folded into d2 analytically (single-layer-linear blocks -> K=1 exact) # ---- the "device": fixed mismatch/offset samples ---- def make_device(rows, cols): return { 'eps_p': rng.normal(0, SIG_BIT, (rows, cols, 8)), 'eps_n': rng.normal(0, SIG_BIT, (rows, cols, 8)), 'vos_col': rng.normal(0, VOS_SIG, cols), 'vos_row': rng.normal(0, VOS_SIG, rows), } DEV1 = make_device(D_IN, D_H) DEV2 = make_device(D_H, D_OUT) def code_eff(codes, eps): bits = ((codes[..., None].astype(int) >> np.arange(8)) & 1) w = (2 ** np.arange(8)) * (1 + eps) return (bits * w).sum(-1) / 256.0 def quant_codes(W): Wp = np.clip(W, 0, None); Wn = np.clip(-W, 0, None) scale = max(np.abs(W).max(), 1e-6) return (np.round(Wp / scale * 255).astype(int), np.round(Wn / scale * 255).astype(int), scale) def spice_mvm(W, x, dev, transpose=False): """y = W^T x (forward, x over rows) or y = W d (transpose, d over cols) via ngspice DC.""" cp, cn, scale = quant_codes(W) gp = code_eff(cp, dev['eps_p']) / R0 # rows x cols conductances gn = code_eff(cn, dev['eps_n']) / R0 rows, cols = W.shape c = Circuit('mvm') if not transpose: for i in range(rows): c.V(f'x{i}', f'r{i}', c.gnd, float(x[i])) outs = [] for j in range(cols): for i in range(rows): if gp[i, j] > 1e-9: c.R(f'p{i}_{j}', f'r{i}', f'sp{j}', 1.0 / gp[i, j]) if gn[i, j] > 1e-9: c.R(f'n{i}_{j}', f'r{i}', f'sn{j}', 1.0 / gn[i, j]) for tag in ('sp', 'sn'): c.B(f'amp{tag}{j}', f'o{tag}{j}', c.gnd, v=f'-1e6*(v({tag}{j})-{dev["vos_col"][j]})') c.R(f'f{tag}{j}', f'o{tag}{j}', f'{tag}{j}', RF) outs.append(j) sim = c.simulator(temperature=27, nominal_temperature=27) an = sim.operating_point() y = np.array([float(an[f'osp{j}']) - float(an[f'osn{j}']) for j in outs]) else: d = x for j in range(cols): c.V(f'd{j}', f'c{j}', c.gnd, float(d[j])) for i in range(rows): for j in range(cols): if gp[i, j] > 1e-9: c.R(f'p{i}_{j}', f'c{j}', f'sp{i}', 1.0 / gp[i, j]) if gn[i, j] > 1e-9: c.R(f'n{i}_{j}', f'c{j}', f'sn{i}', 1.0 / gn[i, j]) for i in range(rows): for tag in ('sp', 'sn'): c.B(f'amp{tag}{i}', f'o{tag}{i}', c.gnd, v=f'-1e6*(v({tag}{i})-{dev["vos_row"][i]})') c.R(f'f{tag}{i}', f'o{tag}{i}', f'{tag}{i}', RF) sim = c.simulator(temperature=27, nominal_temperature=27) an = sim.operating_point() y = np.array([float(an[f'osp{i}']) - float(an[f'osn{i}']) for i in range(rows)]) y = -y / RF * R0 * scale # undo -RF*G scaling back to W-units return y + rng.normal(0, READ_NOISE / RF * R0 * scale, y.shape) def behav_mvm(W, x, dev, transpose=False): cp, cn, scale = quant_codes(W) gp = code_eff(cp, dev['eps_p']); gn = code_eff(cn, dev['eps_n']) Weff = (gp - gn) * scale y = (Weff.T @ x) if not transpose else (Weff @ x) return y + rng.normal(0, READ_NOISE / RF * R0 * scale, y.shape) def ideal_mvm(W, x, dev, transpose=False): return (W.T @ x) if not transpose else (W @ x) # ---- data: 4 Gaussian blobs in 8-D ---- NTR = 256 centers = rng.normal(0, 1.0, (D_OUT, D_IN)) Xtr = np.concatenate([c + 0.35 * rng.normal(0, 1, (NTR // D_OUT, D_IN)) for c in centers]) Ytr = np.concatenate([np.full(NTR // D_OUT, k) for k in range(D_OUT)]) perm = rng.permutation(NTR); Xtr, Ytr = Xtr[perm], Ytr[perm] def softmax(z): e = np.exp(z - z.max()); return e / e.sum() def train(mvm, tag): r = np.random.default_rng(42) W1 = r.normal(0, 0.3, (D_IN, D_H)) W2 = r.normal(0, 0.3, (D_H, D_OUT)) losses = [] t0 = time.time() for step in range(STEPS): idx = r.integers(0, NTR, BATCH) gW1 = np.zeros_like(W1); gW2 = np.zeros_like(W2); L = 0.0 for s in idx: x, yl = Xtr[s], Ytr[s] a1 = mvm(W1, x, DEV1) # analog MVM 1 h = np.tanh(a1) # digital activation (residency theorem) a2 = mvm(W2, h, DEV2) # analog MVM 2 p = softmax(a2) L -= np.log(max(p[yl], 1e-9)) d2 = p.copy(); d2[yl] -= 1.0 # top nudge (digital loss head) d1 = mvm(W2, d2, DEV2, transpose=True) * (1 - h ** 2) # ANALOG transpose read gW2 += np.outer(h, d2); gW1 += np.outer(x, d1) W1 -= LR * gW1 / BATCH; W2 -= LR * gW2 / BATCH losses.append(L / BATCH) if step % 50 == 0: print(f'[{tag}] step {step} loss {losses[-1]:.4f} ({time.time()-t0:.0f}s)', flush=True) # final train accuracy acc = 0 for s in range(NTR): h = np.tanh(mvm(W1, Xtr[s], DEV1)); p = softmax(mvm(W2, h, DEV2)) acc += int(p.argmax() == Ytr[s]) print(f'[{tag}] final loss {np.mean(losses[-20:]):.4f} | train acc {acc/NTR*100:.1f}% ' f'| wall {time.time()-t0:.0f}s', flush=True) return np.array(losses), acc / NTR l_ideal, a_ideal = train(ideal_mvm, 'ideal') l_behav, a_behav = train(behav_mvm, 'behavioral') l_spice, a_spice = train(spice_mvm, 'SPICE') fig, ax = plt.subplots(figsize=(8.6, 4.8)) k = np.ones(10) / 10 for l, name, col in ((l_ideal, f'ideal numpy (acc {a_ideal*100:.0f}%)', '#888888'), (l_behav, f'behavioral non-idealities (acc {a_behav*100:.0f}%)', '#2c6fbb'), (l_spice, f'SPICE-in-the-loop (acc {a_spice*100:.0f}%)', '#d95f02')): ax.plot(np.convolve(l, k, 'valid'), label=name, lw=1.6, color=col, alpha=0.9) ax.set_xlabel('training step'); ax.set_ylabel('train CE (smoothed)') ax.set_title('Stage C: every MVM (forward + transpose) through ngspice — 8-bit codes, mismatch,\n' 'offsets, measured read noise. Ladder certified if SPICE ≈ behavioral.') ax.legend(fontsize=9) fig.tight_layout() fig.savefig('/home/yurenh2/ept/assets/figs/fig_spice_loop.png', dpi=150) print('DONE_SPICE_LOOP')