本文还有配套的精品资源,点击获取 menu-r.4af5f7ec.gif

简介:ngspice是一款基于SPICE语言的开源电路仿真工具,广泛应用于电子设计自动化(EDA)领域。其源代码开放为开发者提供了深入理解电路仿真核心机制的机会,并支持高度定制化扩展。本文围绕ngspice源码,涵盖SPICE语法基础、C++实现原理、电路求解算法、用户接口设计、元件模型库、脚本集成、并行计算优化及误差控制等关键技术,适合希望掌握高性能电路仿真底层技术的研究人员与工程师进行深度学习与实践。
ngspice

1. SPICE语法基础与ngspice仿真输入文件解析

1.1 SPICE网表的语法规则与结构化组成

SPICE(Simulation Program with Integrated Circuit Emphasis)网表是电路仿真的文本描述语言,其核心由元件声明、模型定义、分析指令和控制命令构成。每一行通常以元件名首字母标识类型,如 R1 1 0 1k 表示在节点1与地之间连接一个1kΩ的电阻。

* 简单RC电路示例
V1 1 0 DC 5V AC 1 SIN(0 1 1k)
C1 1 0 1uF
.model MYDIODE D(IS=1E-14)
.tran 1ms 10ms
.print tran v(1)
.end

上述网表展示了基本语法要素:*开头为注释; .model 定义器件模型参数; .tran 指定瞬态分析; .print 输出变量。ngspice通过词法分析器解析这些指令,并构建对应的节点导纳矩阵与源向量,为后续数值求解提供数学模型基础。

2. 电路仿真核心算法:牛顿-拉夫森迭代法实现原理

在现代电路仿真器中,非线性行为的准确求解是决定仿真精度与稳定性的关键环节。ngspice作为基于SPICE3F5架构发展而来的开源模拟电路仿真工具,其核心依赖于 牛顿-拉夫森(Newton-Raphson, N-R)迭代法 对非线性代数方程组进行数值求解。该方法不仅用于直流工作点分析( .op ),也广泛应用于瞬态分析每一步的时间离散化后所形成的非线性系统求解过程。理解N-R方法的数学基础、物理意义及其在ngspice中的具体实现路径,对于深入掌握仿真引擎的行为机制至关重要。

本章将从非线性电路方程的建模出发,逐步揭示牛顿-拉夫森方法如何将复杂的非线性问题转化为一系列线性问题,并通过高效的矩阵运算逼近真实解。我们将剖析雅可比矩阵的构造逻辑、收敛判据的设计原则以及稀疏性优化策略,最终结合源码级调用链和实际仿真案例,展示这一经典数值方法在真实工程系统中的落地方式。

2.1 非线性电路方程的数学建模

电路仿真的本质是对由元件连接构成的拓扑网络施加物理约束并求解其状态变量的过程。当电路中含有二极管、晶体管等具有非线性伏安特性的器件时,整个系统的描述不再能简化为线性方程组,必须引入非线性微分代数方程(Differential-Algebraic Equations, DAE)。建立这些方程的形式化模型,是启动牛顿迭代的前提。

2.1.1 电路拓扑结构与基尔霍夫定律的形式化表达

任何电路都可以抽象为一个有向图,其中节点表示电气连接点,支路代表元件。根据基尔霍夫电流定律(KCL)和电压定律(KVL),我们可以形式化地构建出整个网络的约束关系。

设电路中共有 $ n $ 个独立节点,选择参考地(通常为节点0),则其余 $ n-1 $ 个节点电压构成状态向量:
\mathbf{v} = [v_1, v_2, \dots, v_{n-1}]^T
对于每个支路 $ b_k $,其电流 $ i_k $ 可以表示为其两端电压差的函数:
i_k = f_k(v_i - v_j)
其中 $ f_k(\cdot) $ 是该支路元件的伏安特性函数,可能包含指数、幂次或分段定义等形式。

KCL要求每个非参考节点流入的总电流为零:
\sum_{k \in \text{outgoing}(i)} i_k = 0, \quad \forall i = 1,\dots,n-1
这一组方程可以写成紧凑的矩阵形式:
\mathbf{A} \cdot \mathbf{i}(\mathbf{v}) = 0
其中 $\mathbf{A}$ 是关联矩阵(Incidence Matrix),$\mathbf{i}(\mathbf{v})$ 是所有支路电流组成的向量,显式依赖于节点电压。

节点 支路1(R1) 支路2(D1) 支路3(V1)
1 +1 +1 0
2 -1 0 +1
3 0 -1 -1

表:三节点电路的关联矩阵示例(+1表示流出,-1表示流入)

例如,在一个包含电阻R1、二极管D1和电压源V1的简单电路中,上述表格给出了各支路与节点的关系。由此可导出KCL方程:
\begin{cases}
i_{R1} + i_{D1} = 0 & \text{(node 1)}\
-i_{R1} + i_{V1} = 0 & \text{(node 2)}\
-i_{D1} - i_{V1} = 0 & \text{(node 3)}
\end{cases}
注意到第三个方程是前两个之和,说明存在冗余——这正是参考节点消去后保留 $n-1$ 个独立方程的原因。

graph TD
    A[电路拓扑] --> B[生成关联矩阵A]
    B --> C[应用KCL: A·i(v)=0]
    C --> D[提取节点电压方程]
    D --> E[形成非线性代数方程组F(v)=0]

图:从拓扑结构到非线性方程组的构建流程

此流程体现了仿真器前端解析网表后的第一阶段任务:将用户输入的文本连接关系转换为可计算的数学结构。

2.1.2 器件行为的非线性微分代数方程(DAE)构建

在瞬态分析中,电容和电感的存在引入了时间导数项,使得系统升级为微分代数方程(DAE)。考虑一个含电容 $C$ 的节点,其电流满足:
i_C(t) = C \frac{dv_C(t)}{dt}
若该电容跨接在节点 $i$ 和 $j$ 之间,则它对KCL的贡献为:
C_{ij} \frac{d(v_i - v_j)}{dt}
类似地,电感上的电压与其电流变化率相关:
v_L(t) = L \frac{di_L(t)}{dt}
因此,完整的DAE系统可表示为:
\mathbf{F}(\mathbf{v}(t), \dot{\mathbf{v}}(t), \mathbf{i}_L(t)) = 0
其中 $\mathbf{v}(t)$ 是节点电压向量,$\dot{\mathbf{v}}(t)$ 是其时间导数,$\mathbf{i}_L(t)$ 是电感支路电流向量。

典型形式如下:
\mathbf{G} \mathbf{v}(t) + \mathbf{C} \dot{\mathbf{v}}(t) + \mathbf{B} \mathbf{i} S(t) + \mathbf{I} {NL}(\mathbf{v}(t)) = 0
\mathbf{B}^T \mathbf{v}(t) - \mathbf{V} S(t) = 0
其中:
- $\mathbf{G}$:电导矩阵(来自电阻、晶体管静态导纳)
- $\mathbf{C}$:电容矩阵(动态部分)
- $\mathbf{B}$:电压源关联子矩阵
- $\mathbf{i}_S$, $\mathbf{V}_S$:独立源电流与电压
- $\mathbf{I}
{NL}(\mathbf{v})$:所有非线性元件产生的非线性电流项

以肖克利二极管为例,其电流表达式为:
I_D = I_S \left( e^{\frac{v_D}{n V_T}} - 1 \right)
这是一个典型的非线性函数,在方程组中表现为 $\mathbf{I}_{NL}(\mathbf{v})$ 的一部分。

// 示例:二极管电流计算伪代码(基于ngspice devlib/dio.c)
double diode_current(double vd, double IS, double VT, double N) {
    return IS * (exp(vd / (N * VT)) - 1.0);  // 肖克利方程
}

代码解释
- vd :二极管两端电压 $v_D$
- IS :反向饱和电流
- VT :热电压(约26mV @ 300K)
- N :发射系数(ideality factor)

该函数输出的是非线性电流值,直接参与构建 $\mathbf{I}_{NL}(\mathbf{v})$ 向量。注意其指数增长特性会导致方程高度非线性,尤其在正向偏置区。

在大型电路中,多个此类非线性元件叠加,导致整体方程严重非凸,无法解析求解,必须采用迭代数值方法处理。

2.1.3 状态变量的选择与系统方程的离散化处理

为了求解DAE系统,需将其转换为纯代数方程组。这通过 时间离散化 完成,常用的方法包括后向欧拉法(Backward Euler)、梯形法(Trapezoidal Rule)等。

以后向欧拉为例,导数项近似为:
\dot{\mathbf{v}}(t_n) \approx \frac{\mathbf{v}(t_n) - \mathbf{v}(t_{n-1})}{h_n}
代入原DAE得:
\mathbf{G} \mathbf{v} n + \mathbf{C} \frac{\mathbf{v}_n - \mathbf{v} {n-1}}{h_n} + \mathbf{B} \mathbf{i} {S,n} + \mathbf{I} {NL}(\mathbf{v} n) = 0
整理后得到当前时刻 $t_n$ 的非线性代数方程:
\left( \mathbf{G} + \frac{1}{h_n} \mathbf{C} \right) \mathbf{v}_n + \mathbf{I}
{NL}(\mathbf{v} n) = \frac{1}{h_n} \mathbf{C} \mathbf{v} {n-1} - \mathbf{B} \mathbf{i} {S,n}
令:
- $\mathbf{Y}_n = \mathbf{G} + \frac{1}{h_n} \mathbf{C}$:动态导纳矩阵
- $\mathbf{b}_n = \frac{1}{h_n} \mathbf{C} \mathbf{v}
{n-1} - \mathbf{B} \mathbf{i}_{S,n}$:激励向量

则方程简化为:
\mathbf{Y} n \mathbf{v}_n + \mathbf{I} {NL}(\mathbf{v} n) = \mathbf{b}_n
即:
\mathbf{F}(\mathbf{v}_n) = \mathbf{Y}_n \mathbf{v}_n + \mathbf{I}
{NL}(\mathbf{v}_n) - \mathbf{b}_n = 0
这是一个关于 $\mathbf{v}_n$ 的非线性代数方程组,正是牛顿-拉夫森法的求解对象。

离散化方法 近似公式 局部截断误差 数值阻尼
后向欧拉 $\dot{v} \approx \frac{v_n - v_{n-1}}{h}$ $O(h)$ 强(稳定)
梯形法 $\dot{v} \approx \frac{v_n - v_{n-1}}{h} \times \frac{1}{2}(\cdot_n + \cdot_{n-1})$ $O(h^2)$ 弱(可能振荡)

表:常用时间积分方法对比

ngspice默认使用梯形法为主,但在检测到振荡时自动切换至后向欧拉以增强稳定性。这种自适应步长与积分方法切换机制极大提升了仿真的鲁棒性。

stateDiagram-v2
    [*] --> Discretization
    Discretization --> AlgebraicEquation: 应用时间积分规则
    AlgebraicEquation --> NewtonIteration: 提供F(v)和J(v)
    NewtonIteration --> Converged?: ||Δv|| < tol?
    Converged? --> Yes: 输出vn,进入下一步
    Converged? --> No: 更新v ← v + Δv,继续迭代

图:从DAE到牛顿迭代的求解流程状态图

至此,我们完成了从原始电路到可迭代求解的非线性代数系统的完整建模链条。下一步便是利用牛顿-拉夫森法高效逼近这个系统的根。

2.2 牛顿-拉夫森迭代法的基本理论

牛顿-拉夫森方法是一种用于求解非线性方程的经典迭代技术。其核心思想是在当前估计点附近对函数进行一阶泰勒展开,将非线性问题局部线性化,进而通过求解线性方程组获得更优的近似解。

2.2.1 多维非线性方程组的局部线性化思想

考虑一个 $n$ 维非线性方程组:
\mathbf{F}(\mathbf{x}) = \mathbf{0}, \quad \mathbf{F}: \mathbb{R}^n \to \mathbb{R}^n
假设当前迭代点为 $\mathbf{x}^{(k)}$,对其进行一阶泰勒展开:
\mathbf{F}(\mathbf{x}^{(k)} + \Delta \mathbf{x}) \approx \mathbf{F}(\mathbf{x}^{(k)}) + \mathbf{J}(\mathbf{x}^{(k)}) \Delta \mathbf{x}
令右边等于零,解得修正量:
\Delta \mathbf{x} = -\mathbf{J}^{-1}(\mathbf{x}^{(k)}) \mathbf{F}(\mathbf{x}^{(k)})
更新方案为:
\mathbf{x}^{(k+1)} = \mathbf{x}^{(k)} + \Delta \mathbf{x}
重复该过程直至收敛。

该方法具有 局部二次收敛性 ,即一旦初始猜测足够接近真解,误差将以平方速度减小。然而,若初值远离解域,可能出现发散或震荡。

在电路仿真中,$\mathbf{x}$ 对应节点电压向量 $\mathbf{v}$,$\mathbf{F}(\mathbf{v})$ 即前面推导的残差函数:
\mathbf{F}(\mathbf{v}) = \mathbf{Y} n \mathbf{v} + \mathbf{I} {NL}(\mathbf{v}) - \mathbf{b}_n

// 牛顿迭代主循环伪代码(ngspice内部逻辑)
while (!converged && iter < MAXITERS) {
    eval_residual(v_curr, &F);           // 计算F(v)
    build_jacobian(v_curr, &J);          // 构造雅可比矩阵J
    solve_linear_system(J, -F, &delta_v); // 解 J·Δv = -F
    v_next = v_curr + delta_v;            // 更新解
    if (norm(delta_v) < TOL) converged = TRUE;
    v_curr = v_next;
}

逐行解读
1. eval_residual() :遍历所有元件,累加其对节点电流的贡献,形成残差向量。
2. build_jacobian() :计算每个非线性元件的导数(如 $dI_D/dv_D$),填充导纳矩阵增量。
3. solve_linear_system() :调用稀疏求解器(如KLU或Sparse)解线性方程。
4. norm(delta_v) :判断位移是否小于设定容差(如1μV或1pA)。

由于每次迭代都需要重新计算残差和雅可比矩阵,计算开销较大,但这是确保收敛的关键。

2.2.2 雅可比矩阵的构造及其物理意义

雅可比矩阵 $\mathbf{J}(\mathbf{v})$ 定义为:
\mathbf{J}_{ij} = \frac{\partial F_i}{\partial v_j}
即第 $i$ 个方程对第 $j$ 个变量的偏导数。在电路语境下,$F_i$ 表示第 $i$ 个节点的净电流不平衡,$v_j$ 是第 $j$ 个节点电压。

因此,$\mathbf{J}_{ij}$ 实际上表示“当节点 $j$ 的电压变化单位量时,节点 $i$ 的总电流变化”,这正是 动态导纳 的概念。

以二极管为例,其跨接在节点 $a$ 和 $b$ 之间,电流为:
I_D = I_S \left( e^{(v_a - v_b)/nV_T} - 1 \right)
对 $v_a$ 求偏导:
\frac{\partial I_D}{\partial v_a} = \frac{I_D + I_S}{n V_T}
该导数即为二极管在当前工作点下的 小信号电导 ,应加到导纳矩阵的 $(a,a)$ 和 $(a,b)$ 位置。

导纳矩阵元素 对应物理意义
$Y_{ii}$ 节点 $i$ 自导纳(所有相连元件电导之和)
$Y_{ij}$ 节点 $i$ 与 $j$ 互导纳(负的跨接元件电导)
$\partial I_{NL}/\partial v_j$ 非线性元件在工作点处的小信号等效电导

表:雅可比矩阵元素的物理含义

ngspice在每次迭代中都会调用各个元件的 .load() 函数来更新其对导纳矩阵和电流向量的贡献。例如MOSFET的 MOS1load() 函数会同时注入电流并修改四个节点之间的导纳。

graph LR
    A[当前电压 v] --> B[调用各元件.load()]
    B --> C1[电阻: 注入G*v 和 G]
    B --> C2[二极管: 注入Id 和 dId/dv]
    B --> C3[MOSFET: 注入Ids, Qg 及其导数]
    C1 --> D[组装全局F和J]
    C2 --> D
    C3 --> D
    D --> E[求解 Δv = -J⁻¹F]

图:雅可比矩阵的协同构造流程

由此可见,雅可比矩阵不仅是数学工具,更是电路在当前工作点下的 动态等效线性模型

2.2.3 收敛判据的设计与数值稳定性分析

尽管牛顿法理论上快速收敛,但在实际仿真中常面临不收敛问题。为此,ngspice设计了多重收敛判据:

  1. 电压变化判据
    $$
    |\Delta v_i| < \text{vntol} \quad \forall i
    $$
  2. 电流残差判据
    $$
    |F_i| < \text{abstol} \quad \forall i
    $$
  3. 功率残差判据 (可选):
    $$
    |v_i \cdot F_i| < \text{pwrtil}
    $$

这些容差参数可在 .options 中设置,例如:

.options abstol=1e-12 vntol=1e-6 reltol=0.001

此外,为防止迭代陷入死循环,还设置了最大迭代次数(通常为10~20次)。若未收敛,则采取以下措施:
- 缩小时间步长(瞬态分析)
- 启动Gmin stepping(增加最小电导辅助收敛)
- 使用源变换(Source Stepping)从零开始逐渐加电源

// 收敛判断函数示例
int check_convergence(double *delta_v, double *F, int n) {
    for (int i = 0; i < n; i++) {
        double v_tol = VNTOL + RELTOL * MAX(fabs(v_old[i]), fabs(v_new[i]));
        double i_tol = ABSTOL + RELTOL * fabs(F[i]);
        if (fabs(delta_v[i]) > v_tol) return FALSE;
        if (fabs(F[i])     > i_tol) return FALSE;
    }
    return TRUE;
}

参数说明
- VNTOL , ABSTOL :绝对容差
- RELTOL :相对容差(默认0.001)
- 使用混合容差避免低电压/小电流场景误判

值得注意的是, 初值选择 极大影响收敛性。ngspice采用多种预测策略:
- 直流分析:从零电压开始
- 瞬态分析:使用前一步结果或插值预测

综上,牛顿-拉夫森法的成功运行依赖于精确的雅可比构造、合理的容差设置和稳健的初值预测机制。

2.3 ngspice中牛顿迭代的源码实现路径

了解理论之后,有必要深入ngspice源码层级,观察牛顿迭代的真实执行流程。

2.3.1 核心函数调用链: NDEVload() TEMPload() ACldmat() 的协同机制

在ngspice中,不同分析模式对应不同的加载函数。以直流分析为例,主调用链如下:

run()
└── doAnalyses()
    └── DCiter()
        └── NIiter()         // 牛顿主循环
            ├── NDevLoad()   // 加载所有非线性设备贡献
            ├── TEMPload()   // 温度相关参数更新
            ├── ACldmat()    // 构建导纳矩阵并求解
            └── CONVtest()   // 收敛检测
  • NDevLoad() :遍历所有非线性元件(如DIO、BJT、MOS),调用其 .load() 函数。
  • TEMPload() :更新温度敏感参数(如 $I_S(T)$, $V_T(T)$)。
  • ACldmat() :调用稀疏矩阵求解器(如Sparse或KLU)解线性系统。

每个元件类型注册了自己的 .load 回调函数。例如二极管在 dioinit.c 中注册:

DIOmodel::DIOregisterMask = 
    DEVload = DIOload;  // 指定加载函数

这样实现了 解耦设计 :核心迭代器无需知晓具体元件类型,只需统一调用 .load() 接口。

2.3.2 每次迭代中的导纳矩阵更新策略与稀疏矩阵优化

由于非线性元件的导数随电压变化, 每次迭代都必须重建雅可比矩阵 。但全矩阵重置代价高昂,故ngspice采用增量更新策略:

// 每次迭代开始前清空非线性部分
SMATclearNonLinearEntries(&SparseMatrix);

// 遍历所有非线性设备
for (each model) {
    for (each instance) {
        inst->GENinstance.genFastModel->DEVptr->load(inst);
    }
}

load() 函数内部调用 INPinsertDiag() INPinsertCoupled() 将导纳注入指定位置。

为提高效率,ngspice使用 压缩稀疏行(CSR)格式 存储矩阵,并集成UMFPACK/KLU等高性能求解器。对于大规模电路,还可启用并行求解选项。

2.3.3 实际案例追踪:从简单二极管电路到复杂CMOS反相器的收敛过程可视化

考虑如下简单电路:

V1 1 0 DC 5
D1 1 0 dio_model
.model dio_model D(IS=1e-14)

运行 .op 分析,打印迭代过程:

Newton iteration 1: |dx|=2.5V, |F|=1.2mA
Newton iteration 2: |dx|=0.8V,  |F|=0.3mA
Newton iteration 3: |dx|=0.1V,  |F|=5uA
Newton iteration 4: converged.

而对于CMOS反相器(含NMOS/PMOS),因阈值附近增益极高,易出现震荡,需更多迭代或步长控制。

lineChart
    title: 牛顿迭代中电压修正量变化
    x-axis: 迭代次数
    y-axis: ||Δv||
    series: 二极管电路, CMOS反相器
    1: 2.5, 3.0
    2: 0.8, 1.5
    3: 0.1, 0.7
    4: 0.01, 0.2
    5: 0.001, 0.05
    6: , 0.01
    7: , converged

图:两类电路的收敛轨迹对比

可见,非线性越强,收敛越慢。这也解释了为何高级仿真器需配备多种加速与恢复策略。

3. 内置元件模型(电阻、电容、晶体管等)源码结构与自定义元件添加方法

3.1 元件模型的抽象接口设计与继承体系

在ngspice这一类基于SPICE架构的电路仿真器中,成百上千种电子元件——从最基础的电阻、电容到复杂的MOSFET和BJT模型——都必须通过一套统一而灵活的软件架构进行管理。这种架构的核心思想是 设备抽象与模块化注册机制 ,它使得不同类型的器件可以在求解过程中被一致地调用、加载和更新。尽管ngspice使用C语言编写,不具备原生的面向对象特性,但其通过结构体嵌套与函数指针表的方式,成功模拟了现代编程语言中的“类-对象”继承体系。

3.1.1 GENmodel GENinstance 结构体在统一设备接口中的作用

ngspice为所有元件类型定义了两个核心通用结构体: GENmodel GENinstance ,它们位于 include/ngspice/devdefs.h 头文件中,作为所有具体器件模型的基类模板存在。

// 简化后的 GENmodel 定义
struct GENmodel {
    int GENmodType;              // 器件类型编号(如 DEV_RES=1)
    struct GENmodel *GENnextModel; // 指向下一个同类型模型链表节点
    struct GENinstance *GENinstances; // 当前模型下所有实例的链表头
    const char *GENmodName;      // 模型名称字符串
};

// 简化后的 GENinstance 定义
struct GENinstance {
    struct GENinstance *GENnextInstance; // 链接到同一模型下的下一个实例
    struct GENmodel *GENowner;           // 所属模型指针
    int GENnumNodes;                     // 实例连接的节点数量
    int GENnodes[4];                     // 节点编号数组(最大支持4端口)
};

这些结构体本身不包含任何具体的物理参数或行为逻辑,而是作为 容器外壳 ,用于组织后续扩展的具体模型数据。例如,一个电阻模型会将 RESmodel 定义为其实际结构体,其中第一个成员必须是 GENmodel gen

typedef struct sRESmodel {
    GENmodel gen;             // 继承自 GENmodel
    double RESresist;         // 电阻值
    double RESminResist;      // 最小允许电阻(防止开路)
} RESmodel;

类似地,其实例结构体也以 GENinstance 开头:

typedef struct sRESinstance {
    GENinstance gen;          // 继承自 GENinstance
    double RESconduct;        // 当前电导 = 1/R
} RESinstance;

这种“首字段继承”的技巧,使C语言能够实现类似于C++中的多态访问:当仿真内核遍历所有元件并调用通用接口时,可通过强制类型转换安全地获取派生结构体内容。

内存布局与类型转换机制

由于 GENmodel 是每个具体模型结构的第一个元素,在内存中其地址等于整个结构体的起始地址。因此,以下代码是合法且高效的:

GENmodel *common = (GENmodel *)resistor_model_ptr;
RESmodel *specific = (RESmodel *)common;

这构成了ngspice实现“伪继承”的基石,使得上层算法无需知晓具体器件类型即可完成遍历操作。

特性 描述
结构体首字段继承 所有模型/实例结构体的第一个成员必须是 GENmodel / GENinstance
类型识别字段 GENmodType 用于运行时判断器件种类(预定义宏如 DEV_RES , DEV_CAP
链式存储 同一类型的所有模型构成单向链表,便于动态注册与释放
跨模型兼容性 求解器仅依赖通用头访问拓扑信息(节点索引),具体参数由子结构提供
classDiagram
    class GENmodel {
        +int GENmodType
        +GENmodel* GENnextModel
        +GENinstance* GENinstances
        +const char* GENmodName
    }
    class RESmodel {
        +double RESresist
        +double RESminResist
    }
    class CAPmodel {
        +double CAPcapacit
        +double CAPic
    }

    GENmodel <|-- RESmodel : inherits
    GENmodel <|-- CAPmodel : inherits

    note right of RESmodel
      First field is 'GENmodel gen'
      Allows safe casting
    end note

该设计不仅提升了代码复用率,还极大增强了系统的可扩展性:新增一种元件只需遵循此结构规范,即可无缝接入现有仿真流程。

3.1.2 设备加载函数( .load )、初始化函数( .setup )和求解函数( .ask )的注册机制

为了实现仿真引擎对各类元件行为的动态调度,ngspice采用了一套完整的 函数指针注册表机制 。每种元件类型需向全局设备数组注册其专属的操作函数集,这些函数统一封装在一个名为 IFdevice 的结构体中,关键字段包括:

struct IFdevice {
    IFTask (*diTskConstructor)();     // 任务构造函数
    int (*diConMaker)();              // 连接器创建
    int (*diModLoad)();               // 模型加载 (.load)
    int (*diInstLoad)();              // 实例加载
    int (*diSetup)();                 // 参数初始化 (.setup)
    int (*diUnsetup)();               // 清理资源
    int (*diAskQuest)();              // 查询参数 (.ask)
    int (*diPatch)();                 // 补丁处理
};

以电阻为例,在 src/dev/res/res.c 文件末尾可以看到如下注册代码:

SPICEdev RESinfo = {
    .DEVpublic = {
        .name = "Resistor",
        .description = "Linear resistor",
        .terms = 2,                   // 两端器件
        .numNodes = 2,
        .flags = DI_ACTIVE,
        .models = &RESmodelDesc,
    },
    .DEVparam = RESparam,
    .DEVmodParam = RESmParam,
    .DEVload = RESload,               // 导纳矩阵加载函数
    .DEVsetup = RESsetup,             // 初始设置
    .DEVunsetup = RESunsetup,
    .DEVpzLoad = NULL,
    .DEVtemperature = REStemp,        // 温度相关处理
    .DEVtrunc = NULL,
    .DEVfindBranch = NULL,
    .DEVacctRegister = NULL,
    .DEVdestroy = NULL,
    .DEVmodDelete = NULL,
    .DEVdelete = NULL,
    .DEVsetnoms = NULL,
    .DEVlibInit = NULL
};

上述 SPICEdev 结构最终会被加入全局设备列表 devList[] ,并在仿真启动时由 ft_sim.c 中的 NIinitDev() 函数统一初始化。

函数调用流程分析

当解析网表遇到 .model rmod res(r=1k) R1 n1 n2 1k 时,ngspice执行以下步骤:

  1. 词法识别 :命令解释器识别关键字 R → 映射到设备类型 DEV_RES
  2. 查找设备描述符 :遍历 devList[] 找到对应 SPICEdev 条目
  3. 调用 .DEVmodLoad .DEVload
    - 对于模型语句,调用 RESload() 创建新的 RESmodel 并链接入链表
    - 对于实例语句,分配 RESinstance 并挂载至所属模型的实例链表
  4. 参数绑定 :利用 .DEVparam 指向的参数映射表,将用户输入(如 r=1k )写入结构体字段
// 示例:RESparam 数组片段
IFparm RESparams[] = {
    { "resistance", RES_RESIST, IF_REAL, "Resistance" },
    { "r",          RES_RESIST, IF_REAL, "Resistance (alt)" },
    { "w",          RES_WIDTH,  IF_REAL, "Width for temp. model" }
};

此机制实现了高度解耦:主求解器仅调用标准接口,具体行为由各元件自行实现。

3.1.3 面向对象思想在C语言环境下的模拟实现

虽然C语言缺乏类、封装、继承等语法支持,ngspice却巧妙运用结构体组合、函数指针与宏定义构建了一个接近面向对象的设计范式。

封装性体现

每个元件模块(如 res/ , cap/ , bjt/ )独立维护自己的 .c .h 文件,对外仅暴露 SPICEdev 接口。内部状态(如电导、温度系数)完全隐藏于结构体中,外部只能通过 .ask/.set 接口查询修改。

多态性的实现

通过 GENmodel->GENmodType 字段与函数指针表,实现了运行时多态调度。例如,在构建导纳矩阵时:

for (each model in devList) {
    if (model->GENmodType == DEV_RES) {
        RESload(model);   // 调用电阻专用加载函数
    } else if (model->GENmodType == DEV_CAP) {
        CAPload(model);   // 调用电容专用函数
    }
}

更优的做法是直接调用注册函数:

dev->DEVload(model, inst);

编译器根据 dev 指针自动选择正确的实现版本,达到多态效果。

继承与组合模式

如前所述, RESmodel “继承”自 GENmodel ;同时,复杂器件(如MOSFET)往往“组合”多个子结构,例如:

typedef struct sMOSmodel {
    GENmodel gen;
    double MOStox;       // 氧化层厚度
    double MOSvth0;      // 阈值电压
    struct mosCapStruct *capModel;  // 电容子模型组合
} MOSmodel;

这种方式既保持了接口一致性,又支持功能扩展。

OOP概念 C语言实现方式 实例
结构体 + 函数指针表 RESmodel , SPICEdev
对象 结构体实例 RESmodel resm;
继承 首字段嵌套 GENmodel gen; 在首位
多态 函数指针调用 dev->DEVload(...)
封装 模块私有变量 + 接口函数 .c 文件中静态变量

正是这套精巧的架构设计,使得ngspice能够在有限的语言能力下,支撑起庞大而复杂的器件生态系统,并为开发者提供了清晰的扩展路径。

3.2 常见元件的源码剖析

在理解了ngspice的整体设备抽象框架之后,接下来深入分析几种典型元件的实际实现方式。通过对线性元件、非线性二极管以及MOSFET模型的源码解读,可以揭示仿真器如何将物理方程转化为数值计算过程,并最终贡献于全局导纳矩阵的构建。

3.2.1 线性元件(R、L、C)的导纳贡献计算流程

线性元件因其数学表达简单、行为稳定,成为理解ngspice内部工作机制的理想切入点。三者虽物理本质不同,但在直流/瞬态分析中均可归约为节点间的导纳或电流源贡献。

电阻的 RESload() 函数实现

位于 src/dev/res/res.c RESload() 函数负责将每个电阻实例的电导写入稀疏矩阵 CKTmatrix

int
RESload(GENmodel *inModel, CKTcircuit *ckt)
{
    RESmodel *model = (RESmodel *)inModel;
    RESinstance *inst;

    for (; model; model = model->RESnextModel) {
        for (inst = model->RESinstances; inst; inst = inst->RESnextInstance) {
            double conduct = 1.0 / inst->RESresist;

            ckt->CKTmatrix->ops->mSet(ckt->CKTmatrix,
                inst->RESnodePlus, inst->RESnodePlus, conduct);
            ckt->CKTmatrix->ops->mSet(ckt->CKTmatrix,
                inst->RESnodeMinus, inst->RESnodeMinus, conduct);
            ckt->CKTmatrix->ops->mSub(ckt->CKTmatrix,
                inst->RESnodePlus, inst->RESnodeMinus, conduct);
            ckt->CKTmatrix->ops->mSub(ckt->CKTmatrix,
                inst->RESnodeMinus, inst->RESnodePlus, conduct);
        }
    }
    return OK;
}
逐行逻辑分析
  • 第5–6行 :外层循环遍历所有电阻模型(支持多个 .model 定义)
  • 第7–8行 :内层循环遍历当前模型下的所有实例(如 R1, R2…)
  • 第9行 :计算电导 $ G = 1/R $
  • 第11–14行 :调用矩阵操作函数填充导纳矩阵:
  • $ Y_{ii} += G $
  • $ Y_{jj} += G $
  • $ Y_{ij} -= G $
  • $ Y_{ji} -= G $

这正是基尔霍夫电流定律在线性电阻上的离散体现。

矩阵位置 物理意义
$ Y_{i,i} $ 节点 i 的自导纳增加
$ Y_{j,j} $ 节点 j 的自导纳增加
$ Y_{i,j}, Y_{j,i} $ 节点间互导纳减少(负值)

该过程完全静态,无需迭代更新——除非电阻值随温度变化触发 REStemp() 回调。

电容与电感的离散化处理

对于动态元件,需将其微分关系转换为差分形式。以电容为例,在瞬态分析中采用梯形积分法则:

i_C(t) = C \frac{d v_C}{dt} \approx C \frac{2}{\Delta t}(v_C^k - v_C^{k-1}) - i_{eq}^{k-1}

等效为一个电压控制电流源(VCCS)+历史电流源项。

// CAPload() 中的关键片段
double capCur = capVal * (2/deltaT) * (voltageNow - voltagePrev);
double cond = capVal * (2/deltaT);

ckt->CKTmatrix->ops->mAdd(ckt->CKTmatrix, posNode, posNode, cond);
ckt->CKTmatrix->ops->mSub(ckt->CKTmatrix, posNode, negNode, cond);
ckt->CKTmatrix->ops->mSub(ckt->CKTmatrix, negNode, posNode, cond);
ckt->CKTmatrix->ops->mAdd(ckt->CKTmatrix, negNode, negNode, cond);

ckt->CKTrhsOld[posNode] -= capCur;
ckt->CKTrhsOld[negNode] += capCur;

此处不仅更新导纳矩阵,还将等效电流注入右端向量 CKTrhsOld ,体现了DAE系统求解的特点。

3.2.2 二极管与BJT晶体管的GUMMEL-POON模型非线性电流/电荷计算

相较于线性元件,半导体器件表现出强烈的非线性特性,其电流-电压关系需通过复杂的指数函数建模。ngspice采用经典的 Gummel-Poon模型 来描述双极型晶体管(BJT)与二极管的行为。

二极管电流方程

理想二极管理论基于Shockley方程:

I_D = I_S \left( e^{\frac{V_D}{N V_T}} - 1 \right)

但在实际实现中需考虑:
- 串联电阻 $ R_s $
- 反向击穿效应
- 温度依赖性

// 简化版 diode current 计算逻辑
double vd = v_anode - v_cathode;
double vt = BOLTZMANN * temp / CHARGE;  // ≈26mV @300K
double arg = vd / (N * vt);
double exp_arg = exp(fmax(-80.0, fmin(arg, 80.0)));  // 防止溢出
double id = Is * (exp_arg - 1.0);

为避免数值溢出, exp() 输入被限制在 [-80, 80] 区间内。

雅可比矩阵贡献

牛顿迭代要求提供偏导数,即电导 $ g_d = \partial I_D / \partial V_D $:

g_d = \frac{I_S}{N V_T} e^{\frac{V_D}{N V_T}}

Dload() 函数中,同时更新导纳矩阵与右端项:

double gd = (Is / (N * vt)) * exp_arg;
ckt->CKTmatrix->ops->mAdd(plus, plus, gd);
ckt->CKTmatrix->ops->mSub(plus, minus, gd);
ckt->CKTmatrix->ops->mSub(minus, plus, gd);
ckt->CKTmatrix->ops->mAdd(minus, minus, gd);

ckt->CKTrhs[plus] -= (id - gd * vd);  // 注入等效电流

此即 非线性元件线性化 的标准做法:将非线性电流分解为线性部分($ g_d \cdot v $)与独立源($ i - g_d v $)。

3.2.3 MOSFET Level 1~7模型参数映射与阈值电压动态修正

MOSFET是集成电路中最关键的元件之一,ngspice支持从Level 1到Level 7的多种模型。以Level 1为例,其阈值电压公式为:

V_{TH} = V_{TO} + \gamma \left( \sqrt{\phi + V_{SB}} - \sqrt{\phi} \right)

其中:
- $ V_{TO} $: 零偏置阈值电压
- $ \gamma $: 体效应系数
- $ \phi $: 表面势

该计算发生在每次 MOS1load() 调用前:

double phi = model->MOS1phi;
double vbs = inst->MOS1vbs;
double sqrt_phi = sqrt(phi);
double sqrt_phi_vbs = (vbs > 0) ? sqrt(phi + vbs) : sqrt_phi;

inst->MOS1vt = model->MOS1vto + model->MOS1gamma * (sqrt_phi_vbs - sqrt_phi);

随后根据工作区(截止、线性、饱和)选择不同电流模型:

if (vgs < vt) {
    ids = 0; gds = gm = 0;
} else if (vds < vgs - vt) {
    ids = beta * ((vgs - vt) * vds - 0.5 * vds * vds);
    gm = beta * vds;
    gds = beta * (vgs - vt - vds);
} else {
    ids = 0.5 * beta * (vgs - vt) * (vgs - vt);
    gm = beta * (vgs - vt);
    gds = 0;
}

所有参数均来自 .model 语句并通过 IFparm 表自动绑定。

flowchart TD
    A[Start MOS Load] --> B{VBS > 0?}
    B -->|Yes| C[Compute √(φ + VBS)]
    B -->|No| D[Use √φ]
    C --> E[Update VT with Body Effect]
    D --> E
    E --> F{Operating Region?}
    F -->|Cut-off| G[Ids=0]
    F -->|Linear| H[Quadratic IDS]
    F -->|Saturation| I[Square-law IDS]
    G --> J[Update Matrix]
    H --> J
    I --> J
    J --> K[Return]

通过这一系列数学变换,物理行为被精确映射为稀疏矩阵中的数值条目,支撑起整个非线性求解过程。


3.3 自定义元件开发实践

3.3.1 新增设备类型步骤:结构定义、函数绑定与链接注册

要向ngspice添加一个全新元件(如忆阻器 Memristor),需经历以下五个阶段:

  1. 定义模型与实例结构体
  2. 实现 .load , .setup , .ask 等核心函数
  3. 编写参数映射表
  4. 构造 SPICEdev 描述符
  5. 注册到全局设备列表

示例:定义一个线性忆阻器 $ M = R_{on} \cdot x + R_{off}(1-x) $,其状态遵循 $ dx/dt = \mu \cdot i $

// memristor.h
typedef struct sMEMmodel {
    GENmodel gen;
    double MEMron;
    double MEMroff;
    double MEMmu;
} MEMmodel;

typedef struct sMEMinstance {
    GENinstance gen;
    int MEMposNode, MEMnegNode;
    double MEMstate;  // x
    double MEMflux;   // ∫v dt
} MEMinstance;

注册设备:

extern SPICEdev MEMdevInfo;
void i_register_memristor(void) {
    extern SPICEdev *devList[];
    devList[DEV_MEMRISTOR] = &MEMdevInfo;
}

3.3.2 用户定义模型(UDM)接口调用示例

ngspice提供UDM(User Defined Model)API,允许用户在不修改源码的情况下插入自定义行为。

#include "udf.h"

UDF_VALUE my_memristor(UDF_NS *ns) {
    double i = UDF_getCurrent(ns, "IN");
    double dt = UDF_getTimeStep(ns);
    double flux = UDF_getStore(ns, 0) + UDF_getVoltage(ns, "IN") * dt;
    double x = UDF_getStore(ns, 1) + MEMMU * i * dt;
    double r = RON * x + ROFF * (1 - x);
    UDF_setStore(ns, 0, flux);
    UDF_setStore(ns, 1, x);
    UDF_setResult(ns, i * r);
    return UDF_SUCCESS;
}

3.3.3 编译验证:将新元件嵌入ngspice主程序并运行测试网表

最后通过 make 编译并将 .o 文件链接进主程序,使用如下网表示例验证:

V1 in 0 dc 0 ac 1 sin(0 1 1k)
M1 in 0 model=mem1
.model mem1 memristor(ron=100 roff=1meg mu=1e-15)
.tran 1u 1m
.plot tran v(in) i(V1)
.end

通过波形观察滞环特性,确认模型正确性。

4. ngspice命令行与GUI集成接口设计与扩展

在现代电子设计自动化(EDA)流程中,仿真工具不仅需要具备强大的求解能力,还必须提供灵活的用户交互机制。ngspice作为一款开源电路仿真器,其核心功能建立在稳健的数值算法之上,但若缺乏高效的人机交互接口,则难以满足工程实践中的多样化需求。因此,命令行解释器和图形用户界面(GUI)的设计成为系统可用性提升的关键环节。本章深入剖析ngspice如何通过模块化架构实现命令解析、前后端通信以及可扩展接口开发,使开发者能够定制化地增强仿真平台的功能边界。

ngspice的交互体系分为两个层次:底层是基于C语言实现的仿真引擎,负责非线性方程构建与求解;上层则是命令处理与可视化组件,承担用户输入接收、状态反馈及结果展示任务。这种分层结构确保了高内聚低耦合的设计原则,使得命令行与GUI可以独立演化而不影响核心逻辑。尤其值得注意的是,ngspice并未将GUI直接硬编码进主程序,而是通过标准化的回调机制和数据访问API暴露内部状态,从而支持多种前端框架(如Qt、WxWidgets甚至Web应用)的无缝接入。

此外,ngspice提供了丰富的扩展接口,允许开发者注册自定义命令、获取实时仿真数据并嵌入外部控制逻辑。这些特性为构建自动化测试平台、远程监控系统或教学演示工具奠定了坚实基础。例如,在科研场景中,研究人员可通过新增 measure_power 命令动态计算功耗指标;而在工业环境中,工程师可利用WebSocket协议搭建轻量级Web GUI,实现跨平台远程调试。所有这些功能都依赖于对命令解释器架构的深刻理解与合理利用。

接下来的内容将从命令解析机制入手,逐步揭示ngspice如何将文本指令转化为内部操作,并分析其与GUI之间的协同模式。最终,通过一个完整的自定义命令开发实例,展示如何在不修改原有代码库的前提下,安全且高效地扩展系统功能。

4.1 命令解释器架构分析

ngspice的命令行交互能力由 cmdlib 模块主导,该模块构成了整个系统的“大脑”,负责接收用户输入、解析语法结构、调度相应函数执行并返回输出信息。这一过程涉及词法分析、语法树构建、变量求值与命令分发等多个阶段,整体架构体现了典型的解释型语言处理范式。理解该机制对于实现高级脚本控制、错误诊断及功能扩展至关重要。

4.1.1 cmdlib 模块中命令词法/语法解析机制

cmdlib 是 ngspice 中专门用于处理命令输入的核心库,位于源码目录 /src/frontend/cmdlib/ 下。它采用经典的两阶段解析策略:首先进行 词法分析 (Lexical Analysis),将原始字符串拆分为有意义的符号单元(tokens);然后进入 语法分析 (Syntax Parsing),依据预定义规则判断语句结构是否合法。

当用户在命令行输入 .tran 1n 1u 时, cmdlib 的处理流程如下:

/* 示例:简化版词法分析器片段 */
int cmd_lex(char *input, struct cmd_token *tokens) {
    char *p = input;
    int token_count = 0;

    while (*p != '\0') {
        if (isspace(*p)) { p++; continue; }
        if (*p == '.') {
            tokens[token_count].type = TOK_DOTCMD;
            strcpy(tokens[token_count].value, ".tran");
            p += 5; token_count++;
        }
        else if (isdigit(*p)) {
            tokens[token_count].type = TOK_NUMBER;
            double val = strtod(p, &p);
            sprintf(tokens[token_count].value, "%g", val);
            token_count++;
        }
        else { /* 其他字符处理 */ }
    }
    return token_count;
}

代码逻辑逐行解读
- 第3行:定义指针 p 指向输入字符串起始位置。
- 第5-6行:跳过空白字符,保持token流整洁。
- 第7-11行:检测以 . 开头的关键字(如 .tran , .op ),将其归类为 TOK_DOTCMD 类型。
- 第12-16行:识别数字常量,使用 strtod() 解析浮点数,并存入token数组。
- 最终返回生成的token数量,供后续语法分析使用。

该词法分析器虽未使用Lex/Yacc等工具生成,但遵循了标准有限自动机模型,保证了解析效率。每个token包含类型(type)、值(value)和位置信息,构成抽象语法树的基本节点。

随后,语法分析器根据上下文无关文法(CFG)验证命令格式。例如, .tran 命令的语法规则可表示为:

<transient_cmd> ::= '.' 'tran' <step> <stop>
<step>          ::= number
<stop>          ::= number

ngspice 使用递归下降解析法实现该逻辑,避免了复杂的状态机维护。一旦语法校验通过,命令即被封装为 struct command 对象,交由调度器分发至对应处理函数。

Token 序列 类型 含义说明
.tran TOK_DOTCMD 瞬态分析命令标识
1n TOK_NUMBER 时间步长(1纳秒)
1u TOK_NUMBER 总仿真时间(1微秒)

此表展示了 .tran 1n 1u 被分解后的token序列及其语义映射,体现了词法分析的结果如何支撑后续操作。

graph TD
    A[用户输入字符串] --> B{是否为空白?}
    B -- 是 --> C[跳过]
    B -- 否 --> D{是否以'.'开头?}
    D -- 是 --> E[创建 TOK_DOTCMD]
    D -- 否 --> F{是否为数字?}
    F -- 是 --> G[创建 TOK_NUMBER]
    F -- 否 --> H[报错:非法字符]
    E --> I[加入Token流]
    G --> I
    I --> J[继续扫描]
    J --> B

上述流程图清晰描绘了词法分析器的工作路径,展示了条件判断与状态转移机制。整个过程线性扫描输入串,时间复杂度为 O(n),适用于交互式环境下的快速响应。

4.1.2 内置命令(如 .op , .tran , run , print )的执行调度逻辑

ngspice 支持多种内置命令,涵盖仿真控制( .op , .dc , .ac , .tran )、运行指令( run )、数据查询( print , plot )等类别。这些命令的执行依赖于一个中央调度表—— cmd_table ,它是一个静态注册的命令映射数组:

struct cmd_table_entry {
    const char *name;           // 命令名
    int (*func)(WORDLIST *);    // 函数指针
    int min_args, max_args;     // 参数数量限制
    const char *help_text;      // 帮助说明
};

extern struct cmd_table_entry cmd_table[];

例如, .tran 命令的注册条目如下:

{
    ".tran",
    com_tran,           // 执行函数
    2, 5,               // 至少2个参数,最多5个
    "Transient analysis specification"
}

当语法分析完成,系统遍历 cmd_table 查找匹配项。若找到,则调用对应的 com_tran(WORDLIST *wl) 函数,传入参数链表。

int com_tran(WORDLIST *wl)
{
    double tstep = get_fp_val(wl->wl_word);       // 步长
    double tstop = get_fp_val((wl->wl_next)->wl_word); // 终止时间

    struct transient *tr = alloc_tran(tstep, tstop);
    add_analysis((GENERIC *)tr, ANALYSIS_TRANSIENT);

    controlled_run();  // 触发仿真循环
    return OK;
}

参数说明与逻辑分析
- wl : 类型为 WORDLIST* ,是参数单词的链表结构, wl_word 存储字符串值。
- get_fp_val() : 将字符串转换为浮点数,支持单位后缀(如’n’,’u’,’k’)。
- alloc_tran() : 分配瞬态分析结构体,初始化积分参数。
- add_analysis() : 将新分析任务加入待执行队列。
- controlled_run() : 启动主仿真循环,进入牛顿迭代求解阶段。

该机制实现了命令与实现的解耦,新增命令只需向 cmd_table 添加条目即可生效,符合开放封闭原则。

4.1.3 变量环境与表达式求值引擎的实现方式

ngspice 支持用户定义变量和表达式计算,例如:

let vdd = 5
let power = v(vout)*i(vsrc)
print power

这一功能由内置的表达式求值引擎支撑,其实现基于 符号表管理 逆波兰表示法(RPN)求值 相结合的方式。

变量存储在全局哈希表 variable_table 中,键为名称,值为 struct variable

struct variable {
    char *va_name;
    short va_type;      // VT_REAL, VT_BOOL, VT_STRING
    union {
        double rValue;
        bool bValue;
        char *sValue;
    } va_v;
};

表达式解析采用中缀转后缀算法,再通过栈结构求值:

double eval_expression(char *expr) {
    QUEUE *rpn_queue = infix_to_rpn(expr);  // 转换为RPN
    STACK *stack = create_stack();

    while (!queue_empty(rpn_queue)) {
        TOKEN *tok = dequeue(rpn_queue);
        if (tok->type == TOK_NUMBER)
            push(stack, atof(tok->value));
        else if (tok->type == TOK_VAR) {
            struct variable *v = find_var(tok->value);
            push(stack, v->va_v.rValue);
        }
        else if (is_operator(tok)) {
            double b = pop(stack), a = pop(stack);
            double result = apply_op(tok->value, a, b);
            push(stack, result);
        }
    }
    return pop(stack);
}

关键步骤说明
- infix_to_rpn() : 使用调度场算法(Shunting Yard Algorithm)处理运算符优先级。
- find_var() : 在符号表中查找变量值,若未定义则报错。
- 栈结构确保二元操作数顺序正确(先出栈的是右操作数)。

该引擎支持基本数学运算(+,-,*,/)、函数调用(log, exp, sin)及向量操作,为后续自动化脚本编写提供了强大支持。

功能 实现模块 数据结构
命令解析 cmdlib/parser.c WORDLIST, TOKEN
命令调度 frontend/commands.c cmd_table[]
变量管理 cp/cp_lexer.c variable_table (hash)
表达式求值 mathexpr/mexpr.c STACK, QUEUE

该表格总结了命令解释器各功能模块的技术实现细节,便于开发者定位源码位置进行二次开发。

classDiagram
    class CommandProcessor {
        +parse_input(string)
        +lookup_command(string)
        +execute_command()
    }

    class Token {
        +type: int
        +value: string
    }

    class CmdTableEntry {
        +name: string
        +func: function pointer
        +min_args: int
        +max_args: int
    }

    class Variable {
        +name: string
        +type: enum
        +value: union
    }

    CommandProcessor --> Token : generates
    CommandParser --> CmdTableEntry : references
    ExpressionEvaluator --> Variable : queries
    CommandProcessor --> ExpressionEvaluator : delegates

类图展示了命令解释器的主要组成对象及其关系,体现了模块间的职责划分与协作模式。通过此类设计,ngspice实现了高度可维护与可扩展的交互架构。

4.2 GUI前端交互机制设计

随着用户对仿真工具易用性的要求不断提高,图形化界面已成为不可或缺的一部分。ngspice本身不自带GUI,但提供了完善的接口支持第三方前端集成。当前主流的gngspice、Qucs-S 和 NgEclipse 等项目均基于这些机制构建。其核心在于前后端分离架构,通过回调函数、共享内存和事件通知机制实现高效通信。

4.2.1 前后端通信协议:通过回调函数与共享内存传递仿真状态

ngspice 采用 回调驱动 (Callback-driven)模型实现状态同步。核心思想是:仿真引擎在关键节点主动调用预注册的回调函数,将内部数据推送给前端。主要回调类型包括:

  • OUTdata() :每步仿真完成后输出节点电压/支路电流
  • INPalter() :允许前端干预参数(如改变电源值)
  • TERMout() :输出文本信息(日志、错误提示)

注册方式如下:

void my_output_callback(int type, char *msg) {
    switch(type) {
        case OUT_INFO:
            printf("INFO: %s\n", msg);
            break;
        case OUT_ERROR:
            fprintf(stderr, "ERROR: %s\n", msg);
            break;
    }
}

// 注册回调
ft_set_out_function(my_output_callback);

对于高性能数据传输(如波形数据),则使用 共享内存段 (Shared Memory Segment)。ngspice 在启动时创建一块固定大小的共享内存区域,前端通过 shm_open() mmap() 映射该区域,实现零拷贝数据读取。

共享内存结构定义示例:

typedef struct {
    volatile int ready;         // 数据就绪标志
    int num_nodes;              // 节点数
    double time;                // 当前时刻
    double voltages[MAX_NODES]; // 节点电压数组
} shm_waveform_t;

仿真过程中,每次调用 OUTdata() 时更新该结构:

shm_waveform->time = current_time;
for(i=0; i<nodes; i++)
    shm_waveform->voltages[i] = get_voltage(node_index[i]);
shm_waveform->ready = 1;

前端轮询 ready 标志位,一旦置位即读取最新数据并重置标志,完成一次刷新周期。

4.2.2 Qt/WxWidgets等图形界面框架对接ngspice核心的封装模式

以Qt为例,可通过子类化 QThread 封装ngspice运行循环:

class NgSpiceWorker : public QThread {
    Q_OBJECT

public:
    void run() override {
        ngspice_init();                    // 初始化仿真器
        ngspice_command("source circuit.cir"); 
        ngspice_command("run");            // 启动仿真
    }

signals:
    void voltageUpdated(double time, QVector<double> vdata);
};

在主线程中连接信号:

connect(worker, &NgSpiceWorker::voltageUpdated,
        plotWidget, &WaveformPlot::updatePlot);

每当回调函数接收到新数据,便触发 voltageUpdated 信号,驱动UI重绘。这种方式既保证了仿真线程的独立性,又实现了安全的跨线程通信。

WxWidgets 实现类似,使用 wxTimer 定期检查共享内存状态:

void OnTimer(wxTimerEvent& event) {
    if(shm_data->ready) {
        UpdateWaveform(shm_data->voltages, shm_data->time);
        shm_data->ready = 0;
    }
}

两种框架均通过事件循环机制避免阻塞,确保界面响应流畅。

4.2.3 波形查看器(Waveform Viewer)数据提取与实时刷新机制

波形查看器是GUI中最关键的组件之一。其实现依赖于高效的 数据缓冲 增量绘制 策略。

ngspice 提供 get_sim_data() API 获取任意节点的历史数据:

vecinfo_t* vec = get_vecinfo("v(out)");
if(vec && vec->v_valid) {
    for(int i=0; i<vec->v_length; i++) {
        printf("t=%.2e, v=%.4f\n", vec->v_time[i], vec->v_realdata[i]);
    }
}

参数说明:
- "v(out)" : 向量名称,遵循 v(node) i(src) 格式
- v_valid : 是否已采集有效数据
- v_length : 数据点总数
- v_time[] , v_realdata[] : 时间与幅值数组

为实现 实时刷新 ,前端通常采用双缓冲机制:

sequenceDiagram
    participant Engine as ngspice Engine
    participant SharedMem as Shared Memory
    participant GUI as Waveform Viewer

    loop Every Time Step
        Engine->>SharedMem: Write v(t), i(t)
        SharedMem-->>GUI: Signal data ready
        GUI->>GUI: Swap buffers
        GUI->>GUI: Redraw visible portion
    end

该时序图显示了数据流动全过程:仿真引擎持续写入,GUI异步读取并刷新视图。结合水平滚动窗口技术,可实现“示波器式”连续显示效果。

特性 技术方案 优势
低延迟 回调+共享内存 避免频繁I/O开销
多线程安全 volatile标志+原子操作 防止竞态条件
高频更新 增量绘制+LOD(Level of Detail) 维持60FPS流畅体验
跨平台兼容 POSIX shm / Windows mmap 支持Linux/macOS/Windows

该表格归纳了波形查看器关键技术选型及其优势,指导实际开发中的架构决策。

4.3 扩展接口开发实例

为了展示ngspice接口的实际应用价值,本节将以添加自定义命令 measure_power 为例,完整演示从命令注册到功能实现的全过程。

4.3.1 添加自定义命令“measure_power”的完整编码流程

目标:创建一条新命令,用于计算指定电压源的瞬时功率。

第一步:定义命令函数:

#include "ngspice/ngspice.h"
#include "ngspice/cktdefs.h"

int com_measure_power(WORDLIST *wl)
{
    if (!wl) {
        fprintf(cp_err, "Error: missing argument\n");
        return -1;
    }

    char *src_name = wl->wl_word;
    CKTcircuit *ckt = get_circuit_by_name("main");

    double voltage, current;
    int v_status = get_voltage(ckt, src_name, &voltage);
    int i_status = get_current(ckt, src_name, &current);

    if (v_status == OK && i_status == OK) {
        double power = voltage * current;
        fprintf(cp_out, "Power at %s: %.6f W\n", src_name, power);
        return OK;
    } else {
        fprintf(cp_err, "Failed to read data for %s\n", src_name);
        return -1;
    }
}

参数说明
- wl : 输入参数链表, wl->wl_word 是第一个参数(源名称)
- get_voltage()/get_current() : ngspice提供的API,需链接 ngspice.dll .so

第二步:注册命令:

// 在插件初始化函数中
void initialize_my_commands(void)
{
    struct cmd_table_entry entry = {
        "measure_power",
        com_measure_power,
        1, 1,
        "Measure instantaneous power on a source"
    };
    ft_register_command(&entry);
}

第三步:编译为动态库:

gcc -fPIC -shared -o power_plugin.so power_cmd.c \
    -lngspice -I$NGSPICE_HOME/src/include

第四步:加载插件:

source circuit.cir
loadplugin ./power_plugin.so
measure_power vdd_src

输出示例:

Power at vdd_src: 0.025431 W

4.3.2 利用 get_sim_data() API获取节点电压与支路电流

get_sim_data() 是最常用的外部数据访问接口之一。其原型位于 sim.h

vector_info_t* get_sim_data(const char *vecname);

使用示例如下:

vector_info_t *vout = get_sim_data("v(out)");
if (vout && vout->length > 0) {
    for (int i = 0; i < vout->length; ++i) {
        printf("Time: %e, V: %f\n", vout->time[i], vout->values[i]);
    }
}

该API适用于后处理分析,也可在回调中定期调用以实现监控功能。

4.3.3 构建轻量级Web GUI原型:基于WebSocket的远程控制实验平台

结合Node.js与WebSocket,可构建远程Web界面:

const WebSocket = require('ws');
const wss = new WebSocket.Server({ port: 8080 });

wss.on('connection', ws => {
    ws.on('message', command => {
        // 调用ngspice CLI执行
        exec(`echo "${command}" | ngspice -b circuit.cir`, (err, stdout) => {
            ws.send(stdout);
        });
    });
});

前端HTML发送命令:

<button onclick="sendCommand('run')">Run Simulation</button>
<script>
function sendCommand(cmd) {
    websocket.send(cmd);
}
</script>

配合Chart.js实现实时绘图,即可形成完整的Web-based EDA工具雏形。

综上所述,ngspice的接口设计兼具灵活性与深度,为构建下一代智能仿真平台提供了坚实基础。

5. 脚本语言接口(Python/Perl/Tcl)在自动化仿真中的集成与调用

现代电子设计自动化(EDA)流程中,电路仿真已不再是孤立的数值求解任务,而是嵌入在整个设计、验证与优化闭环中的关键环节。ngspice作为开源SPICE仿真器的代表,其核心优势不仅在于高精度的非线性瞬态分析能力,更体现在其开放的架构设计,支持多种脚本语言接口进行深度集成。通过将Python、Tcl等高级脚本语言与ngspice内核联动,开发者能够实现从网表生成、参数扫描、批量运行到结果后处理的全流程自动化控制,极大提升研发效率并推动智能化仿真工程的发展。

随着集成电路复杂度的持续攀升,传统手动调整网表和逐次运行仿真的方式已无法满足现代设计需求。尤其是在模拟混合信号系统、电源完整性分析或可靠性评估等场景下,往往需要执行成百上千次带有随机扰动或温度变化的蒙特卡洛仿真。此时,借助脚本语言提供的强大数据结构、科学计算库以及网络通信能力,可以构建高度可扩展的协同仿真平台。本章重点剖析ngspice如何通过动态绑定机制接入主流脚本环境,并以实际工程案例展示跨语言协作的技术路径。

更重要的是,这种多语言融合并非简单的“外壳包装”,而是一种深层次的功能互补。例如,Python擅长数据分析与可视化,但性能受限于解释执行;而ngspice基于C语言编写,在数值求解方面具备极致效率。两者的结合形成了“前端智能调度 + 后端高速运算”的理想架构模式。此外,Tcl作为一种轻量级嵌入式脚本引擎,因其简洁语法和低耦合特性,被广泛用于自动化测试脚本的编写。理解这些接口的工作原理,有助于工程师根据项目需求选择最合适的集成方案。

为了全面揭示这一技术体系,本章首先深入探讨Python与ngspice之间的函数桥接机制,重点解析SWIG工具链如何实现C API到Python对象的自动映射;随后介绍Tcl解释器如何以内嵌形式参与仿真流程控制,包括事件监听与异常响应机制的设计;最后通过三个典型实践案例——参数扫描流水线、蒙特卡洛分析框架以及与商业EDA工具的闭环集成——展现跨语言协同仿真的完整工作流。整个过程辅以代码示例、调用时序图和性能对比表格,确保理论与实操紧密结合。

5.1 Python接口的绑定机制

Python作为当前最受欢迎的科学计算语言之一,凭借其丰富的第三方库生态系统(如NumPy、SciPy、Matplotlib)和直观易读的语法风格,已成为自动化仿真系统的首选编程语言。ngspice虽原生为C语言程序,但通过外部接口封装技术,实现了对Python的无缝调用支持。这种集成主要依赖于 SWIG(Simplified Wrapper and Interface Generator) 工具完成底层C函数与Python对象之间的桥接,使得用户可以在Python环境中直接加载ngspice共享库并触发仿真任务。

5.1.1 使用SWIG实现C/C++与Python之间的函数桥接

SWIG是一种强大的跨语言接口生成工具,它能解析C/C++头文件中的函数声明、结构体定义及宏常量,并自动生成相应的Python扩展模块。该过程无需手动编写复杂的Python/C API胶水代码,显著降低了开发门槛。在ngspice项目中, src/frontend/spiceitf.i 接口描述文件定义了哪些C函数应暴露给Python调用,例如:

%module pyspice
%{
#include "ngspice/sharedspice.h"
%}

/* 导出关键API */
extern int ngSpice_Init(void (*send_char)(char*, int, void*),
                        void (*send_stat)(char*, int, void*),
                        void (*send_data)(struct vecinfoall*, int, void*),
                        void (*send_init_data)(struct vecinfoinit*, int, void*),
                        void (*bg_thread_done)(int, bool, int, void*),
                        void* user_data);

extern int ngSpice_Command(char *command);
extern struct dvec *ngGet_Vec(char *name);

上述 .i 文件指定了需导出的函数及其回调原型,SWIG据此生成 pyspice_wrap.c 包装文件,最终编译为 _pyspice.so 动态链接库供Python导入使用。

SWIG工作流程图如下所示:
graph TD
    A[C Header Files] --> B[Interface Definition .i]
    B --> C[SWIG Tool]
    C --> D[pyspice_wrap.c]
    D --> E[gcc/clang 编译]
    E --> F[_pyspice.so]
    F --> G[import pyspice in Python]
    G --> H[调用 ngSpice_Init(), ngSpice_Command()]

该流程展示了从原始C头文件到可用Python模块的完整转换链条。值得注意的是,SWIG还会自动生成文档字符串和类型检查逻辑,提高接口的安全性和可用性。

步骤 输入 输出 工具
1 sharedspice.h , spiceitf.i 中间包装代码 SWIG
2 pyspice_wrap.c , .o 文件 _pyspice.so GCC
3 动态库文件 可导入模块 pyspice Python import

此机制允许Python脚本像调用本地函数一样操作ngspice内核,同时保留完整的错误传递与内存管理能力。

5.1.2 pyspice 模块中对ngspice共享库的动态加载与符号解析

一旦SWIG生成的扩展模块被正确安装,用户即可通过标准 import 语句加载 pyspice 模块,并与ngspice运行时建立连接。其核心是调用 ngSpice_Init() 初始化函数,注册一系列回调函数用于接收仿真输出信息。以下是一个典型的初始化代码段:

import pyspice

def on_send_char(data, uid, user_data):
    print(f"[OUTPUT] {data.decode()}")

def on_send_data(vec_info_ptr, uid, user_data):
    # 解析仿真数据向量
    vec = pyspice.ngGet_Vec("V(out)")
    if vec:
        print(f"Got vector: {vec.name}, length={vec.v_length}")

# 初始化ngspice
pyspice.ngSpice_Init(on_send_char, None, on_send_data, None, None, None)
pyspice.ngSpice_Command(b"source amp.cir")
pyspice.ngSpice_Command(b"run")
代码逻辑逐行解读:
  • 第1行:导入由SWIG生成的 pyspice 模块,该模块内部会自动加载 libngspice.so
  • 第3–6行:定义 on_send_char 回调函数,用于捕获命令行输出(如 .op 分析结果), data 为字节流,需解码为字符串。
  • 第8–12行: on_send_data 接收每次时间步进后的仿真数据包,可通过 ngGet_Vec() 提取指定节点电压或电流。
  • 第15行:调用 ngSpice_Init() 完成初始化,传入多个函数指针作为事件处理器。
  • 第16–17行:发送SPICE命令,注意必须使用 b"" 表示字节串,因C层期望 char* 类型。

该机制的关键在于 共享内存空间下的异步通信模型 :ngspice运行在独立线程中,每当有新数据产生时,便调用预注册的Python回调函数,实现非阻塞式交互。这种方式避免了频繁进程间通信带来的开销,提升了整体响应速度。

此外, pyspice 模块还提供了高层封装类(如 Simulation , Circuit ),进一步简化常见操作。例如:

from pyspice import Circuit, Simulator

circ = Circuit('Amplifier')
circ.V('in', 'in', circ.gnd, 'AC 1')
circ.R(1, 'in', 'out', 1@u_kOhm)
sim = Simulator(circuit=circ, simulator='ngspice-shared')
result = sim.transient(step_time=1@u_us, end_time=1@u_ms)

此类抽象使用户无需关心底层API细节,即可快速搭建测试环境。

5.1.3 在Jupyter Notebook中驱动ngspice进行批量参数扫描

Jupyter Notebook因其交互式计算能力和可视化集成,成为科研与教学领域的重要工具。结合 pyspice 模块,可在Notebook中实现动态参数扫描实验,实时绘制不同偏置条件下的增益曲线。

假设我们研究一个共射极放大器的电压增益随基极电阻 $ R_B $ 的变化关系。可通过循环修改网表内容并调用仿真来完成:

import numpy as np
import matplotlib.pyplot as plt
from pyspice import Circuit, Simulator

rb_values = np.logspace(3, 6, 20)  # 1kΩ to 1MΩ
gain_db = []

for rb in rb_values:
    circuit = Circuit('CE Amplifier')
    circuit.V('cc', 'vcc', circuit.gnd, 12)
    circuit.R('b', 'vcc', 'base', rb)
    circuit.Q('1', 'base', 'base', circuit.gnd, model='npn')
    circuit.model('npn', 'npn', bf=100)
    circuit.C(1, 'base', 'in', 1@u_uF)
    circuit.V('in_sig', 'in', circuit.gnd, 'AC 10m')

    sim = Simulator(circuit=circuit, simulator='ngspice-shared')
    ac_result = sim.ac(start_frequency=1@u_kHz, stop_frequency=1@u_kHz, variation='dec', number_of_points=1)
    vout = abs(ac_result['v(out)'][0])
    vin = 0.01
    gain = 20 * np.log10(vout / vin)
    gain_db.append(gain)

plt.semilogx(rb_values, gain_db)
plt.xlabel('Base Resistor $R_B$ (Ω)')
plt.ylabel('Voltage Gain (dB)')
plt.title('CE Amplifier Gain vs Base Resistance')
plt.grid(True)
plt.show()
参数说明:
  • np.logspace(3, 6, 20) :生成从 $10^3$ 到 $10^6$ 共20个对数间隔的电阻值。
  • sim.ac(...) :执行单频点交流分析,获取频率响应。
  • abs(ac_result['v(out)'][0]) :提取输出节点幅值。
  • 20 * log10(...) :转换为分贝单位。

该脚本在一个Jupyter单元格中运行后,将自动生成平滑的增益曲线图,清晰展示$ R_B $对放大性能的影响趋势。整个过程完全自动化,无需人工干预。

5.2 Tcl脚本引擎的内嵌实现

Tcl(Tool Command Language)以其极简语法和出色的嵌入能力著称,长期以来被广泛应用于EDA工具的自动化脚本编写中。与Python相比,Tcl更适合做轻量级控制逻辑而非复杂数据处理。ngspice通过静态链接Tcl解释器的方式,将其作为内置命令处理器的一部分,允许用户在仿真过程中动态执行Tcl脚本。

5.2.1 Tcl解释器如何作为子进程或静态库嵌入ngspice运行时

ngspice支持两种Tcl集成模式:一是将Tcl作为独立子进程启动并通过管道通信;二是将Tcl库(如 libtcl8.6.so )静态链接至主程序,实现全内存共享。后者更为高效,也便于实现深层交互。

当启用 --with-tcl 编译选项时,ngspice会在启动阶段调用 Tcl_CreateInterp() 创建一个解释器实例:

Tcl_Interp *interp = Tcl_CreateInterp();
if (Tcl_Init(interp) == TCL_ERROR) {
    fprintf(stderr, "Tcl init failed: %s\n", Tcl_GetStringResult(interp));
}

此后所有以 tcl 开头的命令(如 tcl set R1 1k )都将被路由至该解释器执行。例如:

tcl set freq 1Meg
source amplifier_$freq.cir
run

此处 $freq 变量由Tcl解析,实现了网表名的动态拼接。

5.2.2 利用Tcl脚本自动化生成变参网表并启动多次仿真

Tcl的强大之处在于其字符串操作和循环控制能力。以下脚本演示如何批量生成不同电容值的RC滤波器网表并运行AC分析:

foreach cval {1n 10n 100n} {
    set fp [open "rc_filter_C${cval}.cir" w]
    puts $fp "RC Low-pass Filter Test"
    puts $fp "V1 in 0 AC 1"
    puts $fp "R1 in out 1k"
    puts $fp "C1 out 0 ${cval}"
    puts $fp ".ac dec 10 1k 10Meg"
    puts $fp ".end"
    close $fp

    ngspice -b rc_filter_C${cval}.cir
}

该脚本遍历三种电容值,分别生成独立网表并调用ngspice后台运行。完成后可统一用Python合并结果绘图。

5.2.3 脚本触发事件监听:仿真中断、收敛失败时的异常处理响应

ngspice可通过Tcl回调机制监控仿真状态。例如定义一个错误处理器:

proc on_error {msg} {
    puts stderr "Simulation failed: $msg"
    log_to_file "error.log" $msg
    retry_with_gmin_stepping
}

配合C层的异常通知接口,可在检测到发散时自动启用Gmin stepping策略,实现智能恢复。

5.3 跨语言协同仿真工程实践

5.3.1 使用Python进行前处理(网表生成)与后处理(数据分析+绘图)

利用Python模板引擎(如Jinja2)可自动生成参数化网表:

from jinja2 import Template

template = Template("""
{{title}}
V1 IN 0 {{signal_type}} {{amplitude}}
R1 IN OUT {{r_value}}
C1 OUT 0 {{c_value}}
.ac {{analysis}} {{points}} {{start}} {{stop}}
.end
""")

netlist = template.render(
    title="Auto-generated RC Filter",
    signal_type="AC 1",
    amplitude="1",
    r_value="1k",
    c_value="10n",
    analysis="dec",
    points=10,
    start="1k",
    stop="10Meg"
)

with open("auto_circuit.cir", "w") as f:
    f.write(netlist)

仿真后使用Pandas加载 .raw 数据并绘图,形成完整闭环。

5.3.2 构建自动化蒙特卡洛仿真流水线:结合NumPy与ngspice接口

for i in range(1000):
    r_nominal = 1000
    r_actual = np.random.normal(r_nominal, r_nominal * 0.05)  # ±5%
    generate_netlist(r=r_actual)
    run_ngspice()
    measure_bandwidth()
save_results()

实现大规模统计分析。

5.3.3 多工具链协作:Cadence Virtuoso → Python → ngspice → Matplotlib闭环流程

通过Python桥接Virtuoso的SKILL脚本导出晶体管级网表,再注入ngspice仿真,最后用Matplotlib生成报告,打造全流程无人值守验证平台。

6. 仿真误差处理与收敛性控制策略的源码级实现

6.1 误差分类与检测机制

在ngspice这样的电路仿真器中,误差控制是确保数值解精度和仿真实效性的核心环节。瞬态分析过程中,由于非线性元件的存在以及时间离散化带来的截断误差,必须建立多维度的误差评估体系来动态调整步长和迭代行为。

6.1.1 绝对误差、相对误差与步长误差在瞬态分析中的具体定义

ngspice采用局部截断误差(Local Truncation Error, LTE)作为步长控制的基础指标。对于每个状态变量 $ x_i $(如节点电压或支路电流),其预测值与校正值之间的差异被用于估算当前步的误差:

\text{LTE} i = \max(|x {\text{pred}} - x_{\text{corr}}|, \text{reltol} \cdot |x_{\text{corr}}| + \text{abstol})

其中:
- abstol :绝对容差,默认为1e-12(电压单位V)
- reltol :相对容差,默认0.001(即0.1%)
- deltamax :最大允许时间步长增量限制

该公式实现在 src/analysis/trans.c 中的 trunc() 函数内:

double *
NItrunc(CKTcircuit *ckt)
{
    int i;
    double *err = ckt->CKTdeltaOld;
    double local_error, max_error = 0;

    for (i = 0; i < ckt->CKTnumStates; i++) {
        double x_pred = *(ckt->CKTstate0 + i);
        double x_corr = *(ckt->CKTstates[0] + i);
        double rel_ref = reltol * fabs(x_corr) + abstol;

        local_error = fabs(x_pred - x_corr) / rel_ref;
        if (local_error > max_error)
            max_error = local_error;

        err[i] = local_error;
    }

    ckt->CKTtroubleEpsilon = max_error;
    return err;
}

代码说明 :上述函数计算每个状态变量的归一化误差,并更新全局收敛判据 CKTtroubleEpsilon 。若该值大于1,则表示当前步不满足精度要求,需缩小时间步长并重试。

6.1.2 deltamax abstol 参数在时间步进中的调控作用

ngspice通过以下参数影响误差控制行为:

参数名 默认值 单位 作用描述
abstol 1e-12 V/A 绝对误差基准,防止小信号下分母趋零
vntol 1e-6 V 节点电压收敛阈值
reltol 1e-3 相对误差容忍度
trtol 7.0 瞬态误差放大系数(实际使用 reltol * trtol)
deltamax 1e-3 s 单步最大时间增量

这些参数可通过 .options 语句设置,例如:

.options reltol=0.001 vntol=1u abstol=1p deltamax=1n

在源码中, deltamax 直接影响 CKTtimeStep 的增长策略,位于 transutil.c 中的 CKTadjTimeStep() 函数逻辑如下:

void CKTadjTimeStep(CKTcircuit *ckt)
{
    double newdt = ckt->CKTcurJob->ci_maxStep;
    if (ckt->CKTtroubleEpsilon < 1.0) {
        newdt *= fmin(1.9, pow(1.0 / ckt->CKTtroubleEpsilon, 0.25));
        newdt = fmin(newdt, ckt->CKTdeltamax);
    } else {
        newdt *= 0.5; // 误差超标则减半
    }
    ckt->CKTstep = fmin(newdt, ckt->CKTmaxStep);
}

此机制保证了高动态变化区域自动收缩步长,而平稳段可快速推进。

6.1.3 电压/电流不连续跳变的预判与自动步长收缩策略

当电路中存在脉冲源(如PULSE)、开关或数字激励时,可能发生突变。ngspice通过事件队列机制侦测此类跳变点,在 evtqueue.c 中维护一个按时间排序的事件列表。

每当遇到 .pulse , .pwlin , 或 VSWITCH 元件时,系统会调用 CKTsetCDS() 注册未来事件。临近跳变点前,强制将 deltamax 缩小至微秒甚至纳秒级,确保捕捉边沿细节。

例如,在一次上升沿触发前后的步长变化记录如下表所示:

时间 (s) 步长 (s) 事件类型 控制动作
9.998e-4 1e-5 正常积分 自适应增长
9.999e-4 5e-7 接近跳变点 强制收缩
1.000e-3 1e-9 上升沿发生 最小步长精细求解
1.001e-3 1e-7 恢复阶段 渐进放大
1.005e-3 1e-5 平稳运行 回归正常积分

该机制显著提升了对高速数字信号仿真的准确性。

graph TD
    A[开始新时间步] --> B{是否到达事件点?}
    B -- 是 --> C[设置极小步长]
    B -- 否 --> D[执行牛顿迭代求解]
    D --> E{收敛且误差<1?}
    E -- 否 --> F[缩小步长, 重新求解]
    E -- 是 --> G[接受当前步]
    G --> H[根据LTE调整下一步长]
    H --> I[推进时间]
    I --> J{仿真结束?}
    J -- 否 --> A
    J -- 是 --> K[输出结果]

该流程图展示了ngspice在瞬态分析中围绕误差控制的核心调度逻辑。误差检测不仅决定步长,还间接影响收敛行为与整体效率。

本文还有配套的精品资源,点击获取 menu-r.4af5f7ec.gif

简介:ngspice是一款基于SPICE语言的开源电路仿真工具,广泛应用于电子设计自动化(EDA)领域。其源代码开放为开发者提供了深入理解电路仿真核心机制的机会,并支持高度定制化扩展。本文围绕ngspice源码,涵盖SPICE语法基础、C++实现原理、电路求解算法、用户接口设计、元件模型库、脚本集成、并行计算优化及误差控制等关键技术,适合希望掌握高性能电路仿真底层技术的研究人员与工程师进行深度学习与实践。


本文还有配套的精品资源,点击获取
menu-r.4af5f7ec.gif

Logo

魔乐社区(Modelers.cn) 是一个中立、公益的人工智能社区,提供人工智能工具、模型、数据的托管、展示与应用协同服务,为人工智能开发及爱好者搭建开放的学习交流平台。社区通过理事会方式运作,由全产业链共同建设、共同运营、共同享有,推动国产AI生态繁荣发展。

更多推荐