用Python造一个轻量级电路仿真器：直流/交流/瞬态与非线性器件全解析

简介：这是一套用纯Python开发的轻量级SPICE风格电路仿真工具，面向电路课程教学、电路原理验证，以及希望在Python项目中嵌入仿真能力的开发者。工具不需要借助外部EDA软件，也无需C语言编译环境，完整实现了网表解析、电路拓扑建模和节点电压法矩阵构建，可处理电阻、电容、二极管等线性与非线性元件；分析类型涵盖工作点、直流扫描、交流小信号和时域瞬态四种基础场景，能支撑教学演示到中级仿真练习。压缩包共20个文件、约370KB，代码部分包含6个Python源码文件，按职责拆分为网表解析、器件模型、分析算法、流程控制、图形界面和启动入口，另有10个图标用于界面按钮、2个文本说明和配套配置文件，结构清晰便于阅读和二次开发。目前已有83人学习下载；通过源码可直观理解SPICE仿真器中解析、建模与求解的完整链路，既适合初学者从零掌握核心机制，也可作为自定义轻量仿真组件的参考起点。 想把一个电路仿真器塞进Python脚本里，还要同时搞定直流、交流、瞬态分析和非线性器件，听起来像是个自找麻烦的事。但这个麻烦我踏踏实实踩了两个月，最后写出来的小工具只有一千多行，能跑出和Ngspice对得上的波形。工具没有GUI，交互方式就是把网表文件丢进去，然后拿到一个字典：里面有各节点工作点、频响点列、瞬态波形序列。对我这种经常要批量验证电源网络、偶尔手算核对结果的人来说，已经够用了。

做这个东西之前我翻过不少开源代码，最大感受是：SPICE类工具不缺，缺的是“电路规模小但能塞进自动化流程”的工具。所以我想聊的也不是要写一个完整商业软件，而是分享如何用Python搭出一个轻量电路仿真核心，真正覆盖直流/交流/瞬态分析和非线性器件建模。适合两类人：一是想理解仿真器内部原理的Python开发者，二是想在测试脚本里塞一个可控仿真器的硬件工程师。

1. 为什么我要写一个Python电路仿真器而不是直接用Ngspice

1.1 现有工具的重与轻

先说实话，如果只是偶尔仿真一个三极管放大电路，直接用Ngspice就好，没必要再造轮子。我写这个工具的场景是自动化测试：同一份电路要扫温度、扫容差、换参数，跑几百次，每次还要把输出和理想公式做比对。Ngspice命令行虽然能批量执行，但输出是文本，解析起来绕；PySpice虽然封装得不错，底层依然是SPICE，出了问题很难判断是调用方式错还是模型本身行为。自己写一个轻量核心，所有参数都在Python字典里流转，出了问题可以直接打断点，这是最吸引我的地方。

还有一点是我个人需求：我想在仿真过程中随时拿到中间矩阵，看看哪些节点数值开始异常，哪些元件贡献了小到离谱的导纳。这在现成工具里做起来麻烦，但自研就是在自己的代码里插几行print的事。

1.2 我定义的轻量级边界

既然轻量，功能边界必须明确。我的工具只支持R、C、L、独立电压源/电流源、二极管、以及一个简化到“能用就行”的NMOS/PMOS开关模型。分析类型固定四种：直流工作点、直流扫描、交流扫描、瞬态。不支持传输线、噪声、失真、工艺角统计。网表语法参考SPICE里最常用的那两三成子集。

这个范围听着窄，但覆盖了教科书里大部分模拟电路实验，以及电源设计里最常见的RC滤波、钳位、整流、稳压结构。架构上分三层：最底层是从网表解析出的组件列表；中间是按分析类型组织的求解器；最上层是统一的输出结构。新加一种器件，只需要在“线性化器件模型”里补一个函数，返回它在当前电压下的等效电导和电流源，其他代码完全不用动。

2. 直流工作点：从网表到矩阵，再让非线性收敛

2.1 MNA矩阵里每个元件怎么填

直流分析是所有分析的地基，交流小信号要在直流工作点上线性化，瞬态也需要工作点做初始条件。我用的是修正节点电压法，也就是MNA。未知量不仅是节点电压，还包括电压源支路电流。为什么必须这样？因为电压源在节点方程里是“恒压未知电流”，如果不把电流设成未知量，矩阵就少一行，解出来不唯一。MNA把方程补全，代价是变量数增加，但对轻量工具来说完全可接受。

元件填充矩阵的规则很机械。电阻最简单，导纳加到对应节点的自导纳上，互导纳取负。我常用LIL稀疏矩阵装配，最后转CSR求解：

import numpy as np import scipy.sparse as sp import scipy.sparse.linalg as spla # 节点编号：0为地，真实节点从1开始 def stamp_resistance(A, n1, n2, value): g = 1.0 / value if n1 >= 0: A[n1, n1] += g if n2 >= 0: A[n2, n2] += g if n1 >= 0 and n2 >= 0: A[n1, n2] -= g A[n2, n1] -= g def stamp_voltage_source(A, b, n1, n2, vs_idx, voltage): # vs_idx 是电压源对应的未知电流变量索引 if n1 >= 0: A[n1, vs_idx] += 1 A[vs_idx, n1] += 1 if n2 >= 0: A[n2, vs_idx] -= 1 A[vs_idx, n2] -= 1 b[vs_idx] = voltage

这里用-1表示接地端，处理起来省事。电压源那几行和几列全是±1，最后求解得到的就是各节点电压和电压源支路电流。这一步是后面所有分析的基础。

2.2 二极管迭代的收敛技巧

非线性器件让线性方程组变成固定点问题。我的迭代流程是：先假设所有二极管开路，解出初始节点电压；然后把每个二极管当前位置的电流和电导代入，作为线性化器件；重新解方程，如此反复。每一步本质上是更新Jacobian矩阵，也就是更新电导和电流源。

二极管用Shockley方程：

i_d = I_s * (exp(v_d / (n * V_t)) - 1)

对v_d求导得到小信号电导：

g_d = I_s / (n * V_t) * exp(v_d / (n * V_t))

这里最大的坑是指数溢出。v_d超过0.7V后，指数项轻轻松松超过1e12，double直接爆炸。我用的办法是电压限幅：每次迭代前把二极管两端电压限制在-0.5V到0.9V区间内，如果修正量想超过这个区间，就按比例缩小步长。

另一个必做操作是给每个PN结并联一个很大电阻，比如1e9欧姆。目的不是模拟漏电，而是让Jacobian矩阵在二极管未导通时依然有非零导纳路径，否则矩阵会奇异。

收敛判定也要注意。不要只看电压增量小于1e-6就认为收敛，因为二极管电流随电压指数变化，电压差很小，电流差可能很大。我会同时检查每个非线性器件的电流残差：当|Δi|小于1e-9且|Δv|小于1e-6才算收敛。这个双判据帮我少走了很多弯路。

3. 交流小信号：把电容电阻变成复数矩阵里的一个元素

3.1 在工作点旁边做线性化

交流小信号分析其实是在直流工作点上做一次线性化。二极管用它在工作点的电导代替，直流源全部置零，电压源短路、电流源开路，电容和电感根据频率变成阻抗或导纳。于是每个频率点只需要求解一个复线性方程组：

A(ω) x(ω) = b(ω)

矩阵元素是复数。比如电阻导纳是1/R，电容是jωC，电感是1/(jωL)。注意MNA中的节点电压是相量，b向量通常只保留一个被激励的源，其他独立源置零。所以我在做AC之前，会先跑一遍DC工作点，把非线性器件的小信号参数算出来，然后锁定不动。

代码层面，只需要把之前的实数矩阵换成复数矩阵：

def stamp_capacitor_ac(A, n1, n2, value, freq): y = 1j * 2 * np.pi * freq * value # 电容器件的贡献和电阻类似，只是用复数导纳 y if n1 >= 0: A[n1, n1] += y if n2 >= 0: A[n2, n2] += y if n1 >= 0 and n2 >= 0: A[n1, n2] -= y A[n2, n1] -= y

求解后得到复数电压相量，幅值用np.abs，相位用np.angle。这个结果和手算的复阻抗分压公式非常容易对拍。

3.2 频率扫描怎么跑才快

AC扫描一次会跑几十到几百个频率点，如果每个点都重新做非线性工作点线性化，纯属浪费时间。我的做法是，先用DC求解器把所有非线性的线性化参数存起来，AC求解时只更新与频率相关的电容电感项。因为电阻和受控源的值不变，可以预先填充到稀疏矩阵，再在循环里复制，只修改动态元件对应的位置，能省一半时间。

频率点用logspace生成，从10Hz到10MHz扫50个点很常见。我通常只关心某些探针节点的幅频响应，所以求解完直接存成字典。

测试时用一个RC低通和理论公式比对，cutoff=1/(2πRC)。R=1k，C=100nF，算出来截止频率约1.592kHz，AC扫描结果和理论值误差在0.1%以内，说明求解器没有问题。

4. 瞬态分析：梯形积分背后的等价电流源

4.1 动态元件在时间轴上长什么样

瞬态分析把时间离散成步长h。电容的电流-电压关系是微分关系，直接解微分方程很麻烦，但在每个时间步上，它可以变成一个线性电阻并联电流源。最简单的后向欧拉：

i_C(t) ≈ C * (v(t) - v(t-h)) / h

让v(t)作为未知数，这一项就等效于一个电导G_eq = C/h，并联一个已知电流源C*v(t-h)/h。电感做类似处理，但用诺顿等效更容易。

后向欧拉精度一般，我用得比较多的是梯形法。梯形法对电容的伴随模型是：

i_C(t) ≈ (2C/h) * (v(t) - v(t-h)) - i_C(t-h)

也就是说G_eq = 2C/h，等效电流源I_eq = -(2C/h) * v(t-h) - i_C(t-h)。这里的i_C(t-h)是上一个步长流过电容的电流，需要记录。电感类似，只是公式里的电压电流互换。

用这种思路，每个瞬态步长最终还是回到解同一个MNA矩阵，只是右侧电流源向量每步都在变。对非线性器件，每个时间点还要再跑一遍Newton-Raphson迭代，和直流分析的过程几乎一样。

4.2 步长自适应与振荡抑制

固定步长在RC电路上没什么问题，但遇到LC谐振或开关电路就麻烦了。我一开始用固定1us步长仿真RLC串联电路，结果波形出现明显的高频振铃，一开始以为是算法写错了，后来才知道是梯形法的数值寄生振荡。判断办法很简单：看振荡频率是否和时间步长相关，而非电路本身。

我的处理方式是自适应步长，必要时切换到后向欧拉。局部截断误差可以这样估计：先用步长h算一步，再用两个h/2步算到同一时刻，两个结果的差值就是误差信号。如果误差超过容差，就把h减半重算；如果误差很小，则适当加大h。这个流程会引入约1.5倍计算量，但对几千个时间点的轻量场景无所谓。

当检测到连续三个步长内节点电压交替跳跃、且幅度大于阈值，我就把该时间窗口的积分公式暂时切成后向欧拉。虽然精度掉到一阶，但能有效抑制寄生振荡，等波形平稳再切回梯形法。这个混合策略让我的半波整流仿真稳定了很多。

5. 非线性器件建模：从PN结到简单MOS开关

5.1 二极管的温度与过压处理

非线性器件是标题里的重头戏。我在仿真器中把二极管模型做成一个类，构造函数带Is、n、Vt、Cjo、Rs。Cjo是结电容，瞬态时需要参与电荷存储；如果完全忽略结电容，开关波形会很假。对轻量工具，至少加一个固定结电容并联在理想PN结上，仿真结果就靠谱不少。

温度处理我直接用了SPICE的通用思路：Vt = k*T/q，Is随温度指数上升。实际操作中最麻烦的是，模型参数和温度对收敛的影响比想象中大。比如在25℃时Is=1e-14，温度升到85℃，Is可能变成1e-12，指数项更敏感，原来限压范围可能要调整。我的做法是把电压限幅上限和温度关联，温度越高上限越低。

MOS管部分，我一开始尝试完整Level 1模型，结果收敛性相当差。后来改成一段式平滑模型：用Vgst = log1p(exp(K*(VGS - Vth))) / K这类平滑函数近似VGS-Vth的过零，避免导数跳变。漏电流公式简化为I_D = Kp * W/L * Vgst^2 * (1 + lambda*VDS)，再对VDS做限幅处理。这个模型在8V以内帮我把NMOS开关仿得比较像，代价是牺牲了细节，但换来稳定收敛。

5.2 实测：半波整流、RC低通和稳压管钳位

我拿这个工具做了几个经典电路的复现，方便大家验证。

半波整流，网表很简单：

V1 in 0 SINE(0 5 1k) D1 in out diode_1n4148 R1 out 0 1k .TRAN 0.5u 5m

用1k电阻负载，二极管用1N4148近似模型，仿真得到的输出峰值约4.52V，和Ngspice结果差不到2%。主要误差来源是二极管模型参数和结电容，尤其在高频时差异会拉大。

RC低通前面提过，R=1k，C=100nF，AC扫描从10Hz到10kHz，-3dB点算出来约1.592kHz，和理论值一致。

稳压管钳位用两个背靠背二极管模拟，输入峰值12V，输出被钳在正负6.8V左右。这个实验检验的是瞬态下非线性器件的鲁棒性。我一开始在零点附近反复不收敛，后来加了1e9并联电阻和电压限幅才稳定。

三个实验的结果对比：

所有实验都是直接跑脚本，输出到Python字典，再交给matplotlib绘图。对我个人来说，这比打开一个EDA软件点半天选波形更顺手。

最后说一个我自己总结的小习惯：不管接下来要做交流还是瞬态，先把DC工作点跑通。DC不收敛的时候，其他分析全都不用看；DC收敛了但结果和手算对不上，也先别急着调瞬态参数，多半是线性化模型的问题。仿真器虽小，但这套排查顺序和大工具完全一致。用这个思路，多复杂的非线性电路也能一点点啃下来。

本文还有配套的精品资源，点击获取