简介：这是一份用纯Python实现的轻量级SPICE风格电路仿真工具，面向电路原理学习者、教学演示者以及需要在Python项目中嵌入简单仿真能力的开发者。它无需外部EDA软件或C编译环境，完整覆盖网表解析、电路拓扑建模、节点电压法矩阵构建，并支持工作点、直流扫描、交流小信号与时域瞬态四类分析，同时能处理二极管等非线性器件。资源包共20个文件，约370KB，以6个Python源码模块为主，分别承担网表解析、器件模型定义、分析算法封装、求解流程协调与图形界面入口；另含10个PNG图标、2个TXT说明文档及配置文件，便于直接运行与二次开发。目前已有87人学习下载。读者可借此理解SPICE核心机制，将其作为教学演示、电路验证或嵌入自有项目的参考实现，代码结构清晰，适合初学者按模块逐步研读。
1. 用 Python 写一个能跑直流/交流/瞬态分析的轻量级电路仿真工具，到底值不值得做
很多做硬件或嵌入式的朋友都有过这种经历：板子焊好了，上电前想先看看某个节点的电压电流大概是多少，打开商业仿真软件，光是建工程、选模型、配库就要折腾半小时，最后为了算一个电阻分压，杀鸡用了牛刀。更麻烦的是，商业工具是黑匣子，你想在仿真循环里加一段自己的控制逻辑、跑个参数扫描、或者把结果直接喂给 Python 的数据分析栈，几乎无从下手。这就是「Python 写的轻量级电路仿真工具，支持直流/交流/瞬态分析和非线性器件」这个方向真正要解决的问题：不是替代工业级 EDA，而是把电路求解这件事拉回到几十行到几百行代码就能掌控的尺度。
它适合三类人：一是想搞懂 SPICE 类工具内部到底怎么解方程的工程师，二是需要在自动化流程里嵌入电路求解的开发者，三是教学场景下想让学生看到每一步矩阵运算的导师。核心能力就三块——直流工作点、交流小信号、瞬态时域，再加上二极管、BJT、MOS 这类非线性器件的迭代求解。下面我按「先立住理论、再动手复现、最后讲坑」的顺序，把这条路走一遍。
2. 从节点方程到矩阵求解：轻量级仿真器的骨架怎么搭
2.1 为什么选改进节点分析（MNA）而不是纯节点法
纯节点法（Nodal Analysis）只对节点列 KCL，遇到电压源就傻眼——理想电压源的电流是未知量，没法直接写进节点方程。改进节点分析（Modified Nodal Analysis, MNA）的做法是把电压源支路电流也当成未知量加进未知向量，同时补一条电压约束方程。这样未知量从「节点电压」扩展成「节点电压 + 部分支路电流」，矩阵维度变大，但电压源、受控源、电感（在特定积分方式下）都能统一处理。
我一般会把未知向量排成[V1, V2, ..., Vn, I_v1, I_v2, ...]，前 n 个是除参考地以外的节点电压，后面是电压源支路电流。MNA 矩阵的标准形式是：
[ G B ] [ v ] [ i ] [ C D ] [ j ] = [ e ]其中 G 是节点导纳子矩阵，B 和 C 描述电压源与节点的连接关系，D 通常为零（除非有受控源之间的耦合）。对线性电阻、电容、电感，G 的填充规则很固定；对非线性器件，G 里的元素会随工作点变化，需要迭代。
2.2 用 Python 把 MNA 矩阵组装出来
先定义一个最小可用的线性求解器。下面这段代码只处理电阻、独立电压源、独立电流源，目的是把矩阵组装逻辑讲清楚，后面再往上叠非线性。
import numpy as np class Circuit: def __init__(self): self.nodes = {} # 节点名 -> 索引 self.elements = [] # 元件列表 self.vsources = [] # 电压源索引 def add_node(self, name): if name not in self.nodes: self.nodes[name] = len(self.nodes) return self.nodes[name] def add_resistor(self, n1, n2, R): self.elements.append(('R', n1, n2, R)) def add_vsource(self, n_pos, n_neg, V): idx = len(self.vsources) self.vsources.append((n_pos, n_neg, V)) self.elements.append(('V', n_pos, n_neg, idx)) def add_isource(self, n1, n2, I): self.elements.append(('I', n1, n2, I)) def build_mna(self): n = len(self.nodes) m = len(self.vsources) size = n + m A = np.zeros((size, size)) z = np.zeros(size) for el in self.elements: if el[0] == 'R': _, a, b, R = el ia = self.nodes.get(a, None) ib = self.nodes.get(b, None) g = 1.0 / R if ia is not None: A[ia, ia] += g if ib is not None: A[ib, ib] += g if ia is not None and ib is not None: A[ia, ib] -= g A[ib, ia] -= g elif el[0] == 'I': _, a, b, I = el ia = self.nodes.get(a, None) ib = self.nodes.get(b, None) if ia is not None: z[ia] -= I if ib is not None: z[ib] += I elif el[0] == 'V': _, a, b, k = el ia = self.nodes.get(a, None) ib = self.nodes.get(b, None) row = n + k if ia is not None: A[ia, row] += 1 A[row, ia] += 1 if ib is not None: A[ib, row] -= 1 A[row, ib] -= 1 z[row] = el[2] if len(el) > 2 else 0 return A, z这段代码里几个关键点：self.nodes用字典把节点名映射到矩阵索引，地节点不加入字典，等价于索引 -1，所有涉及地的项直接跳过。电阻的 stamp 是标准的四元素填充，电压源在矩阵右上和左下各填 ±1，右端向量填电压值。build_mna返回的 A 和 z 就是待解的线性方程组A x = z。
参数说明：R单位欧姆，V单位伏特，I单位安培，全部用国际单位制，避免中途换算引入错误。节点名用字符串，方便和网表对应。
2.3 直流工作点求解与结果校验
有了矩阵，直流工作点就是一次线性求解：
def solve_dc(circ): A, z = circ.build_mna() x = np.linalg.solve(A, z) n = len(circ.nodes) voltages = {name: x[idx] for name, idx in circ.nodes.items()} currents = {f"I_V{i}": x[n + i] for i in range(len(circ.vsources))} return voltages, currents拿一个分压电路验证：12V 电源、两个 1k 电阻串联。理论值中间节点 6V，电流 6mA。跑出来如果对不上，先检查地节点有没有被误加进self.nodes，再检查电压源 stamp 的符号。我见过最常见的翻车是电压源正负端填反，导致解出来是 -6V，这种错误在复杂网表里很难肉眼发现，所以每加一个新元件类型，都要用一个已知答案的最小电路回归一次。
3. 交流小信号与瞬态分析：从复数矩阵到时间步进
3.1 交流分析：把实矩阵换成复数矩阵
交流小信号分析的本质是在直流工作点附近做线性化，然后把所有元件用复阻抗表示，扫频求解。电阻还是 R，电容变成1/(jωC)，电感变成jωL。矩阵组装逻辑和直流几乎一样，只是数据类型从 float 变成 complex，右端向量里独立源变成相量。
def build_ac(self, freq): omega = 2 * np.pi * freq n = len(self.nodes) m = len(self.vsources) size = n + m A = np.zeros((size, size), dtype=complex) z = np.zeros(size, dtype=complex) for el in self.elements: if el[0] == 'R': _, a, b, R = el y = 1.0 / R elif el[0] == 'C': _, a, b, C = el y = 1j * omega * C elif el[0] == 'L': _, a, b, L = el y = 1.0 / (1j * omega * L) else: continue ia = self.nodes.get(a, None) ib = self.nodes.get(b, None) if ia is not None: A[ia, ia] += y if ib is not None: A[ib, ib] += y if ia is not None and ib is not None: A[ia, ib] -= y A[ib, ia] -= y # 电压源在交流下幅值设为 1，相位 0，做归一化 for k, (a, b, _) in enumerate(self.vsources): ia = self.nodes.get(a, None) ib = self.nodes.get(b, None) row = n + k if ia is not None: A[ia, row] += 1 A[row, ia] += 1 if ib is not None: A[ib, row] -= 1 A[row, ib] -= 1 z[row] = 1.0 return A, z扫频时对每个频点调用一次，解出节点电压相量，取模得到幅频，取角度得到相频。这里有个容易忽略的点：电感的导纳在直流时是无穷大（短路），在代码里如果 freq=0 会除零，所以交流扫描的起点一般从 1Hz 或更低但非零的频率开始。我一般会在扫频函数里加一句保护，freq 小于 1e-3 时直接跳过或钳位。
3.2 瞬态分析：后向欧拉与时间步长的选择
瞬态分析要把电容、电感离散化。最省事的是后向欧拉（Backward Euler），数值稳定但精度一般；要求高一点用梯形法（Trapezoidal），精度好但可能引入数值振荡。轻量级工具我建议先用后向欧拉把流程跑通。
电容的后向欧拉伴随模型是一个电阻并联一个电流源：等效电导G = C/dt，等效电流源I = C/dt * v_prev。电感类似，等效电导G = dt/L，等效电流源I = i_prev。每推进一步，把当前时刻的等效元件 stamp 进矩阵，解一次线性方程，更新状态变量。
def transient(circ, tstop, dt): t = 0.0 # 初始化电容电压、电感电流为 0 v_cap = {id(el): 0.0 for el in circ.elements if el[0] == 'C'} i_ind = {id(el): 0.0 for el in circ.elements if el[0] == 'L'} results = [] while t < tstop: A, z = circ.build_mna() n = len(circ.nodes) # 把动态元件等效 stamp 叠加进去 for el in circ.elements: if el[0] == 'C': _, a, b, C = el g = C / dt ia = circ.nodes.get(a, None) ib = circ.nodes.get(b, None) if ia is not None: A[ia, ia] += g if ib is not None: A[ib, ib] += g if ia is not None and ib is not None: A[ia, ib] -= g A[ib, ia] -= g # 等效电流源 i_eq = g * v_cap[id(el)] if ia is not None: z[ia] += i_eq if ib is not None: z[ib] -= i_eq elif el[0] == 'L': _, a, b, L = el g = dt / L ia = circ.nodes.get(a, None) ib = circ.nodes.get(b, None) if ia is not None: A[ia, ia] += g if ib is not None: A[ib, ib] += g if ia is not None and ib is not None: A[ia, ib] -= g A[ib, ia] -= g i_eq = i_ind[id(el)] if ia is not None: z[ia] -= i_eq if ib is not None: z[ib] += i_eq x = np.linalg.solve(A, z) # 更新状态 for el in circ.elements: if el[0] == 'C': _, a, b, C = el va = x[circ.nodes[a]] if a in circ.nodes else 0.0 vb = x[circ.nodes[b]] if b in circ.nodes else 0.0 v_cap[id(el)] = va - vb elif el[0] == 'L': _, a, b, L = el # 电感电流由等效模型反推 va = x[circ.nodes[a]] if a in circ.nodes else 0.0 vb = x[circ.nodes[b]] if b in circ.nodes else 0.0 i_ind[id(el)] += (dt / L) * (va - vb) results.append((t, x.copy())) t += dt return results时间步长 dt 的选择直接决定精度和速度。经验做法是取最小时间常数的 1/10 到 1/20。比如 RC 电路 R=1k、C=1uF，时间常数 1ms，dt 取 50us 比较稳。dt 太大，后向欧拉会把波形抹平；dt 太小，Python 循环会慢到怀疑人生。我一般会先跑一个粗步长看趋势，再对关键区间细化。
3.3 非线性器件的迭代：牛顿-拉夫逊怎么接进去
二极管、BJT、MOS 的 I-V 特性是非线性的，不能一次求解。标准做法是牛顿-拉夫逊迭代：在当前猜测的工作点附近把非线性器件线性化，得到等效电导和等效电流源，stamp 进矩阵，解一次，比较新旧解，收敛就退出，不收敛就更新工作点再来。
以二极管为例，肖克利方程I = Is*(exp(V/(n*Vt))-1)，在猜测电压 V0 处求导得到等效电导gd = Is/(n*Vt)*exp(V0/(n*Vt))，等效电流源Ieq = I(V0) - gd*V0。每次迭代把 gd 和 Ieq 更新进矩阵。
def solve_dc_nonlinear(circ, max_iter=100, tol=1e-6): # 初始猜测全部节点电压为 0 x = np.zeros(len(circ.nodes) + len(circ.vsources)) for it in range(max_iter): A, z = circ.build_mna() for el in circ.elements: if el[0] == 'D': _, a, b, Is, n, Vt = el ia = circ.nodes.get(a, None) ib = circ.nodes.get(b, None) va = x[ia] if ia is not None else 0.0 vb = x[ib] if ib is not None else 0.0 Vd = va - vb # 限制指数项防止溢出 Vd = max(min(Vd, 0.8), -5.0) exp_term = np.exp(Vd / (n * Vt)) Id = Is * (exp_term - 1) gd = Is / (n * Vt) * exp_term Ieq = Id - gd * Vd if ia is not None: A[ia, ia] += gd z[ia] -= Ieq if ib is not None: A[ib, ib] += gd z[ib] += Ieq if ia is not None and ib is not None: A[ia, ib] -= gd A[ib, ia] -= gd x_new = np.linalg.solve(A, z) if np.max(np.abs(x_new - x)) < tol: return x_new, it x = x_new raise RuntimeError("牛顿迭代未收敛")这里有两个血泪经验：一是指数项必须钳位，否则迭代初期电压猜测离谱时exp直接溢出成 inf；二是初始猜测不能全零就完事，对多二极管电路，全零猜测可能让所有二极管都处于关断区，迭代在原地打转。常见做法是先用线性电阻网络解一次作为初值，或者给每个 PN 结一个 0.6V 左右的正向偏置猜测。
4. 避坑与排查：轻量级仿真器最容易翻车的五个地方
4.1 矩阵奇异，解出来全是 nan
现象：np.linalg.solve抛LinAlgError: Singular matrix，或者解出来是 nan/inf。
原因：最常见的是某个节点悬空，没有任何支路连接到它，导致该节点对应的矩阵行全零。其次是电压源形成回路，比如两个理想电压源并联且电压不等，MNA 矩阵会出现线性相关行。
解决：组装完矩阵后加一个检查，如果某行对角线元素为 0 且该行其他元素也为 0，说明该节点悬空，报错时把节点名打出来。对电压源回路，要么在网表层面禁止，要么给每个电压源串一个小电阻（比如 1mΩ）做数值上的正则化。
4.2 瞬态波形在换步长后对不上
现象：dt 从 50us 改成 10us，瞬态波形幅值差了一大截。
原因：后向欧拉是一阶方法，本身有数值阻尼，dt 越大阻尼越明显。如果电路里有高频振荡，粗步长会把振荡完全抹掉。
解决：对精度要求高的场景换梯形法，或者用自适应步长——根据相邻两步解的变化量动态调整 dt。轻量级工具里我一般先用固定小步长跑通，确认物理趋势对了，再考虑优化步长。
4.3 非线性迭代不收敛，卡在某个迭代次数
现象：牛顿迭代跑到 max_iter 还没收敛，或者解在两个值之间来回跳。
原因：一是指数项没钳位，数值溢出；二是初始猜测太差；三是电路本身没有直流解（比如正反馈结构）。
解决：先钳位指数项，再给合理的初值。如果还不收敛，把迭代过程的前几步电压打印出来，看是哪个节点在震荡。有时候加一个小的源极/发射极电阻（比如 0.1Ω）能显著改善收敛性，这是工程上常用的技巧。
4.4 交流扫描在低频段结果异常
现象：扫频从 0Hz 开始时，电感支路导纳无穷大，矩阵出现 inf。
原因：1/(jωL)在 ω=0 时除零。
解决：扫频起点设为 1e-3 Hz 或更高，或者在代码里对 ω 加一个极小值保护。另外，如果电路里有纯电容隔直路径，低频段节点电压可能浮空，也需要检查。
4.5 节点命名冲突导致 stamp 错位
现象：电路拓扑明明对，但解出来的电压和手算对不上，差一个符号或差一个系数。
原因：节点名字符串在字典里映射时，大小写或空格不一致，导致同一个物理节点被当成两个不同节点。
解决：在add_node里统一做strip()和lower()处理，或者在网表解析阶段就规范化。我一般会在组装矩阵前打印一次节点映射表，肉眼扫一遍，这个习惯救过我好几次。
5. 把仿真器用起来：参数扫描、结果导出与一个收敛技巧
轻量级工具真正的价值不在单次仿真，而在你能用 Python 把仿真嵌进更大的流程。比如做参数扫描：固定电路拓扑，让某个电阻从 100Ω 扫到 10kΩ，每个值跑一次直流，把结果收集成数组，直接喂给 matplotlib 或 pandas。这段代码不到二十行，但商业工具里往往要写脚本或手动点很多次。
import numpy as np def sweep_resistor(base_circ, res_index, values): results = [] for R in values: # 深拷贝电路，修改指定电阻 import copy c = copy.deepcopy(base_circ) el = list(c.elements[res_index]) el[3] = R c.elements[res_index] = tuple(el) v, i = solve_dc(c) results.append((R, v.get('out', 0.0))) return results # 用法示例：扫描分压电路的下臂电阻 vals = np.logspace(1, 4, 20) data = sweep_resistor(circ, res_index=1, values=vals)参数说明：res_index是元件在elements列表里的位置，改之前先打印一遍列表确认。copy.deepcopy保证每次扫描不污染原始电路。扫出来的(R, Vout)可以直接画成半对数曲线，看分压比随负载变化的趋势。
再讲一个收敛技巧：对含多个二极管的电路，如果牛顿迭代反复不收敛，可以先把所有二极管替换成一个小电阻（比如 1Ω）跑一次线性解，把得到的节点电压作为非线性迭代的初值。这个「先线性后非线性」的两段式做法，比全零初值的收敛率高很多。我一般会把它写成solve_dc_nonlinear的一个可选参数warm_start=True，默认开启。
最后说结果导出。瞬态跑完的results是一个(t, x)列表，直接转成 numpy 数组存成.npz或.csv，后面用 pandas 做后处理、用 scipy 做 FFT 看谐波，都比在仿真器内部硬编码要灵活。我自己的习惯是仿真器只负责解方程和吐数据，所有分析、绘图、报告生成都放到外面，这样仿真核心可以保持足够小，小到你能在一杯咖啡的时间里通读一遍。
这套东西我从最早只能跑电阻网络，到后来能处理二极管整流、BJT 共射放大器的直流偏置，前后迭代了七八个版本。最大的教训是：不要一上来就追求器件模型多全、精度多高，先把 MNA 组装、线性求解、非线性迭代这三步用最小电路跑通，每加一个器件类型就用已知答案回归一次。仿真器这种东西，错一个符号，后面全错，而且错得很隐蔽。希望帮到你。
本文还有配套的精品资源，点击获取
