CppCon 2014 学习:Rolling Your Own Circuit Simulator

这段话讲述了一个背景和动机，目的是阐明为什么开源C++库变得越来越复杂且在科学和工程领域有很大的应用潜力。

关键点：

1. 开源库的成熟：
   • 近年来，开源C++库在许多科学和工程领域变得越来越成熟和强大。
   • 这些库不再仅仅是简单的工具，而是具有足够复杂度和功能，能够处理一些复杂的应用需求。
2. 曾经是专有技术的技术：
   • 在过去，一些复杂的算法和技术（例如电路模拟、数值优化等）通常只存在于学术代码或商业专有软件中。
   • 这些技术现在通过开源库被暴露出来，任何有需求的人都可以使用这些工具，而不再需要依赖昂贵的商业软件或只能通过难以获取的学术代码来实现。
3. 对独立开发者的好处：
   • 这使得独立的开发者、研究人员、工程师等可以在没有过多依赖商业软件或进行大量重复工作时，直接利用这些库作为构建模块来开发他们自己的应用。
   • 这种开放性降低了进入门槛，让更多的人可以利用已有的工具来进行创新和研究。

小结：

这段话强调了通过开源库，复杂的技术和方法变得更容易访问和使用，促使了更多创新的可能性。尤其是在需要进行电路模拟等专业任务时，开发者可以使用现有的成熟库来加速开发，避免从头开始编写所有功能。这种趋势对于科学、工程和软件开发等领域都具有积极意义。

通过介绍如何利用两种开源库，希望展示如何从自己的领域——电气工程出发，通过这两个库解决问题，从而实现强大的效果。

关键点：

1. 目标展示：
   • 演讲者计划展示如何利用这两个开源库来解决实际问题。这里暗示了演讲者将展示一些具体的、实用的应用，可能涉及到复杂的技术或工程问题。
2. 领域背景：
   • 演讲者的背景是电气工程，这可能意味着演讲者将通过这些库来解决电气工程领域中的某些复杂问题，可能涉及电路模拟、信号处理等方面。
3. 强大效果：
   • 通过这两个开源库，演讲者期望展示如何达成强大的技术结果，展示这些库的实际应用价值。

小结：

这段话概述了演讲的目的：展示如何利用现代开源库来处理复杂的电气工程问题，解决实际的挑战，并展示这种方法能够带来强大的效果。演讲者的目标是通过实例或实际应用，向听众证明这些库的强大功能。

如何将电路视为一种图（graph）结构，并描述了图的组成元素。

关键点：

1. 电路作为图：
   • 电路可以被视为一种图结构。
   • 图中的边（Edges）被称为支路（branches），每个支路代表一个电子组件（如电阻、电容等），并且支路中的电流取决于该组件的类型、数值以及它所连接的节点的电压。
2. 节点（Nodes）：
   • 节点是电路中的连接点，电流在这些节点之间流动，形成连接各个组件的通路。
3. 欧姆定律：
   • 公式 $I=\frac{(V_{\text{Node1}}-V_{\text{Node2}})}{R}$ 描述了电流（I）的计算方式，其中：
     • $V_{\text{Node1}}$ 和 $V_{\text{Node2}}$ 是两个节点的电压。
     • $R$ 是节点间组件的电阻。
     • $I$ 是流过电阻 $R$ 的电流。

小结：

• 电路 = 图：将电路看作一个图，其中节点表示电路的连接点，支路表示电路中的各个元件（如电阻、电容等）。
• 电流计算：通过欧姆定律计算电流，电流取决于连接节点的电压差以及组件的阻抗（电阻）。
  这个图结构化的方式有助于程序员以更系统的方式理解电路，并可以在模拟和计算电路时应用图算法。

这是一个简单的低通滤波器电路，包含一个100Ω电阻、一个20nF电容和一个20µH电感。以下是分析：

• 输入电压 (Vin) 施加在电阻和电容并联的两端。
• 电容 (20nF) 连接到地。
• 电感 (20µH) 与输出电压 (Vout) 串联，输出从电感两端取出。
  这种配置允许低频信号通过到输出，同时衰减高频信号。截止频率可以近似使用RC低通滤波器的公式计算：( f_c = \frac{1}{2\pi R C} )，其中 ( R = 100Ω )，( C = 20nF )。计算得出截止频率约为79.6 kHz。电感可能会在更高频率时提供额外的滤波效果。
  如果需要更精确的分析（包括电感的影响），需要进一步计算或仿真。
  

让我们逐步分析这个电路和给定的问题。这是一个包含100Ω电阻、20nF电容和20µH电感的低通滤波器电路，输入 $V_{in}(t)$ 为阶跃函数。我们需要理解给定的方程和输出电压 $V_{out}(t)$ 的表达式。

电路分析

1. 电路结构：
   • 输入电压 $V_{in}$ 施加在电阻 $R=100\Omega$ 和电容 $C=20nF$ 并联的两端。
   • 电容连接到地。
   • 电感 $L=20\mu H$ 与输出 $V_{out}$ 串联，输出从电感两端取出。
2. 节点方程：
   在节点1处，根据基尔霍夫电流定律（KCL），流入和流出节点的电流之和为0：
   $$\frac{V_{out}-V_{in}}{R}+C\frac{d{V}_{out}}{dt}+I_L=0$$
   • 第一项是流过电阻的电流。
   • 第二项是流过电容的电流（电容电流为 $I_C=C\frac{dV}{dt}$）。
   • 第三项是流过电感的电流 $I_L$。
3. 电感方程：
   电感上的电压为：
   $$V_{out}=L\frac{d{I}_L}{dt}$$
   这表示输出电压 $V_{out}$ 是电感电流 $I_L$ 的变化率与电感值 $L$ 的乘积。

输入和输出关系

输入 $V_{in}(t)$ 是一个阶跃函数，即在 $t=0$ 时从0跳变到1V（单位阶跃函数）。我们需要求解输出电压 $V_{out}(t)$。
给定的输出表达式为：
$$V_{out}=\frac{1}{\mu RC}{e}^{\lambda t}\sin (\mu t)$$
其中：
$$\lambda =-\frac{1}{2RC},\mu =\sqrt{\frac{1}{LC}-{\lambda }^2}$$

验证和计算

1. 参数计算：
   • $R=100\Omega$
   • $C=20nF=20\times 10^{-9}F$
   • $L=20\mu H=20\times 10^{-6}H$
     计算 $RC$ 和 $LC$：
     $$RC=100\times 20\times 10^{-9}=2\times 10^{-6}\text{ 秒}$$
     $$LC=(20\times 10^{-6})\times (20\times 10^{-9})=4\times 10^{-13}{\text{ 秒}}^2$$
2. 计算 $\lambda$ 和 $\mu$：
   $$\lambda =-\frac{1}{2RC}=-\frac{1}{2\times 2\times 10^{-6}}=-2.5\times 10^5{\text{ s}}^{-1}$$
   $$\frac{1}{LC}=\frac{1}{4\times 10^{-13}}=2.5\times 10^{12}{\text{ s}}^{-2}$$
   $${\lambda }^2=(2.5\times 10^5{)}^2=6.25\times 10^{10}{\text{ s}}^{-2}$$
   $$\mu =\sqrt{\frac{1}{LC}-{\lambda }^2}=\sqrt{2.5\times 10^{12}-6.25\times 10^{10}}=\sqrt{2.4375\times 10^{12}}\approx 1.562\times 10^6\text{ rad/s}$$
3. 计算 $\frac{1}{\mu RC}$：
   $$\mu RC=(1.562\times 10^6)\times (2\times 10^{-6})\approx 3.124$$
   $$\frac{1}{\mu RC}\approx \frac{1}{3.124}\approx 0.32$$
   因此，输出电压表达式为：
   $$V_{out}(t)=0.32e^{-2.5\times 10^{5}t}\sin (1.562\times 10^{6}t)$$

物理意义

• 这是一个二阶系统的典型响应（因为有电感和电容），表现为阻尼振荡。
• $\lambda$ 表示阻尼因子，决定振荡的衰减速度。
• $\mu$ 是振荡的角频率，决定振荡的频率。
• 输出 $V_{out}(t)$ 是一个衰减的正弦波，初始时从0开始振荡，最终由于阻尼而趋于0。

总结

给定的 $V_{out}(t)$ 表达式是正确的，描述了阶跃输入下该电路的输出响应。输出是一个阻尼振荡信号，随着时间衰减。如果你需要进一步分析（例如绘制波形或计算特定时间点的值），可以告诉我！

这个代码是一个直接实现的例子，用于计算 RLC电路（串联电阻、电感和电容）中输出电压在给定时间点的值。

结构体 circuit 的解释

circuit 结构体表示一个 RLC 电路，并且基于电路的电阻（r）、电感（l）和电容（c）来计算输出电压。这个结构体通过其构造函数初始化电路的参数，并通过运算符 () 实现计算给定时间点的输出电压。

关键参数

• r_：电阻（$R$），单位为欧姆（Ω）。
• c_：电容（$C$），单位为法拉（F）。
• lambda_：与系统的阻尼系数相关的参数，公式为：
  $$\lambda =-\frac{1}{2RC}$$
  这个值控制了电路振荡的衰减速度。
• mu_：与电路自然频率相关的参数，公式为：
  $$\mu =\sqrt{\frac{1}{LC}-{\lambda }^2}$$
  这个值控制了电路振荡的频率。

构造函数

circuit(double r, double l, double c)
: r_(r), c_(c) {lambda_ = -1.0 / (2.0 * r * c);mu_ = sqrt((1.0 / (l * c)) - lambda_ * lambda_);
}

• 构造函数接收三个参数：电阻 r，电感 l，电容 c。
• 它根据这些参数计算出阻尼系数 (lambda_) 和角频率 (mu_)，这些值会在后续计算中使用。

operator() 用于计算输出电压

double operator()(double t) {return (1.0 / (mu_ * r_ * c_)) * exp(lambda_ * t) * sin(mu_ * t);
}

• 通过重载 operator()，这个结构体变得像一个函数一样，可以在给定时间 t 时计算输出电压 $V_{\text{out}}(t)$。
• 输出电压的计算公式是：
  $$V_{\text{out}}(t)=\frac{1}{\mu \cdot R\cdot C}\cdot e^{\lambda t}\cdot \sin (\mu t)$$
  其中：
  • $e^{\lambda t}$ 代表电路中由于电阻而导致的衰减，随着时间推移，电压幅度会逐渐减小。
  • $\sin (\mu t)$ 代表电路中的振荡，由电感和电容共同作用形成。它决定了振荡的频率。
  • $\mu$ 和 $\lambda$ 是根据电路组件的值计算得到的常数。

公式的含义

1. 衰减项：指数项 $e^{\lambda t}$ 代表电阻引起的衰减效应。随着时间推移，电压的振幅逐渐减小。
2. 振荡项：正弦项 $\sin (\mu t)$ 代表电感和电容引起的振荡。这是电路的固有振荡频率，受电感 $L$ 和电容 $C$ 的影响。
3. 输出电压：给定时间 t，根据电路的电阻、电感和电容的参数，计算出输出电压 $V_{\text{out}}(t)$。

示例使用

假设我们有一个如下的RLC电路：

• 电阻 $R=100\Omega$,
• 电感 $L=20\mu H$,
• 电容 $C=20nF$。
  我们可以创建一个 circuit 类的对象，并计算在某个给定时间 t 的输出电压 $V_{\text{out}}(t)$。

int main() {circuit myCircuit(100, 20e-6, 20e-9);  // R = 100Ω, L = 20µH, C = 20nFdouble t = 0.01;  // 时间点double vout = myCircuit(t);  // 计算该时刻的输出电压std::cout << "Output voltage at t = " << t << " is " << vout << std::endl;return 0;
}

总结

• 这个 circuit 结构体直接计算了串联RLC电路在给定时间点的输出电压。
• 计算中涉及到电路的衰减（由电阻引起）和振荡（由电感和电容引起）。
• 该实现通过使用标准的物理公式，提供了一种高效模拟RLC电路行为的方式。

Direct Implementation: Simulation and Plotting

这段代码演示了如何使用之前提到的 circuit 类来模拟 RLC 电路的输出电压，并将结果输出到控制台。

代码解析

int main() {// 创建电路对象，电阻100Ω，电感20μH，电容20nFcircuit ckt(100.0, 20e-6, 20e-9);// 循环遍历时间从0到10微秒，步长为0.1微秒for (double t = 0; t < 10e-6; t += 0.1e-6) {// 输出当前时间和该时刻的输出电压std::cout << t << " " << ckt(t) << std::endl;}
}

步骤解释

1. 创建电路对象：

   circuit ckt(100.0, 20e-6, 20e-9);
   

   这里我们实例化了一个 circuit 对象 ckt，其中：
   • 电阻 $R=100\Omega$
   • 电感 $L=20\mu H$
   • 电容 $C=20nF$
     ckt 对象会根据这些参数计算出 lambda_ 和 mu_，用于后续计算输出电压。
2. 时间步进并计算电压：

   for (double t = 0; t < 10e-6; t += 0.1e-6) {std::cout << t << " " << ckt(t) << std::endl;
   }
   

   • 循环的时间从 t = 0 开始，到 t = 10μs 结束。
   • 时间步长是 0.1μs，即每次增加 0.1e-6 秒。
   • 每次循环，调用 ckt(t) 来计算在时间 t 时刻的输出电压 $V_{\text{out}}(t)$，并将当前时间 t 和电压值输出。

输出解释

每次循环，程序会输出一个时间点 t 和对应的电压值 $V_{\text{out}}(t)$，例如：

0.0 0.0
0.1e-6 0.15
0.2e-6 0.28
...

Numeric Integration: A More General Approach

在传统的电路分析中，使用闭式解法（例如求解微分方程）来得到电路的行为是常见的做法。然而，这种方法有其局限性，尤其在复杂电路中，无法总是通过简单的闭式解来解决问题。为了克服这些问题，我们可以使用 数值积分 方法来模拟电路的行为。

问题的挑战：

1. 手动解决微分方程：传统方法通常需要我们手动解决微分方程，这既繁琐又容易出错。
2. 无法一般化：有些方程无法通过闭式解解决。许多复杂的电路包含非线性或变系数的微分方程，无法得到简单的解析解。

Boost.ODEInt 库：

Boost.ODEInt 是 Boost 库中的一个数值积分库。它自 2012 年 9 月被加入 Boost，并提供了一个简单的接口用于数值积分。使用该库，你可以自动地解决微分方程，避免手动推导解析解，同时它也支持复杂方程的求解。

Boost.ODEInt 的基本使用要求：

1. 状态定义：你需要定义一个表示电路状态的容器。这个容器保存电路的所有状态变量，例如电压、电流等。
2. 计算状态变化的可调用对象：这通常是一个函数或者一个类，计算每个时刻电路状态的变化率（例如电流、温度等随时间的变化）。
3. 初始条件：需要提供电路的初始状态。
4. 结果输出：你需要一个可调用对象来接收计算结果，通常用于打印或者绘制图形。

如何在我们的电路中使用 Boost.ODEInt：

假设我们使用 Boost.ODEInt 来求解电路中的微分方程。我们不再需要手动推导电路的每个微分方程，而是通过数值方法直接求解。
步骤：

1. 定义电路的状态：我们需要定义电路的状态向量 $X(t)$，通常包含电压和电流。
2. 定义微分方程：我们需要提供一个计算电路状态变化的函数，通常会根据欧姆定律和基尔霍夫定律来推导微分方程。
3. 初始化电路状态：设置电路的初始条件（例如初始电压、电流）。
4. 调用 Boost.ODEInt 进行数值积分：利用 Boost.ODEInt 库求解微分方程并输出结果。

代码示例

以下是一个使用 Boost.ODEInt 库来模拟电路的基本示例。假设我们有一个简单的 RLC 电路，并且用 Boost 进行数值求解。

#include <boost/numeric/odeint.hpp>
#include <iostream>
#include <vector>
using namespace boost::numeric::odeint;
// 定义电路的状态
typedef std::vector<double> state_type;
// 电路的微分方程
void circuit_equations(const state_type &x, state_type &dxdt, double t) {// 假设电路是一个RLC电路：Vin = L * dIL/dt + Vout// dxdt[0] 表示电压的变化，dxdt[1] 表示电流的变化double R = 100.0;  // 电阻 100Ωdouble L = 20e-6;  // 电感 20μHdouble C = 20e-9;  // 电容 20nFdxdt[0] = x[1] / C;               // 电压的变化率：dVout/dt = I/Cdxdt[1] = -(x[0] / L) - (x[1] / R); // 电流的变化率：dIL/dt = -Vout/L - IL/R
}
// 输出结果
void write_output(const state_type &x, double t) {std::cout << t << " " << x[0] << " " << x[1] << std::endl;
}
int main() {// 初始条件：Vout = 0, IL = 0state_type x(2, 0.0);  // x[0]是电压，x[1]是电流double t_start = 0.0;   // 起始时间double t_end = 10e-6;   // 结束时间double dt = 0.1e-6;     // 时间步长// 使用Runge-Kutta方法进行数值积分integrate(circuit_equations, x, t_start, t_end, dt, write_output);return 0;
}

代码解释：

1. 状态定义：state_type 是一个 std::vector<double>，用来表示电路的状态。在这个例子中，x[0] 是电压 $V_{\text{out}}$，x[1] 是电流 $I_L$。
2. 微分方程：circuit_equations 函数实现了电路的微分方程，根据欧姆定律和基尔霍夫电压定律（KVL）来推导电路状态的变化。
   • 电压的变化率：$\frac{d{V}_{\text{out}}}{dt}=\frac{I_L}{C}$
   • 电流的变化率：$\frac{d{I}_L}{dt}=-\frac{V_{\text{out}}}{L}-\frac{I_L}{R}$
3. 初始条件：在 main 函数中，我们设定了电路的初始条件：电压 $V_{\text{out}}(0)=0$，电流 $I_L(0)=0$。
4. 数值积分：integrate 函数会进行数值积分，使用 Runge-Kutta 方法对电路的状态进行时间步进。在每个时间点，它会调用 write_output 来输出电路状态（电压和电流）。

总结

• Boost.ODEInt 库提供了一个强大的数值积分工具，能够让你轻松求解复杂的微分方程，适用于无法通过解析解法求解的电路。
• 使用该库，你只需要定义电路的微分方程和初始条件，库会自动处理数值积分过程，避免了手动推导和求解微分方程的繁琐。
• 这种方法非常适合用于模拟复杂的电路，尤其是在多种电气元件和动态变化的情况下。

Numeric Integration with ODEInt: Circuit Definition

在电路仿真中，我们通常需要解一个系统的微分方程。对于RLC电路，状态变量可以是电流和电压（例如，电感电流 $I_L$ 和输出电压 $V_{\text{out}}$）。为了模拟电路的行为，我们将其转化为一个微分方程系统，并使用数值积分方法进行求解。Boost.ODEInt 是一个用于求解微分方程的强大工具，它可以帮助我们直接对电路模型进行数值求解。

电路状态定义

在我们给出的代码中，电路的状态是两个独立的变量：

• $I_L$（电感电流）
• $V_{\text{out}}$（输出电压）
  这些变量构成了电路的状态向量。为了将这些变量与其微分方程的变化率结合，我们需要使用数值积分来更新这些状态。

代码解析

typedef std::array<double, 2> state_t;

在这段代码中，state_t 定义为一个大小为2的 std::array<double, 2>，它表示电路的状态。数组的第一个元素（x[0]）是电感电流 $I_L$，第二个元素（x[1]）是输出电压 $V_{\text{out}}$。

struct circuit {circuit(double r, double l, double c): r_(r), l_(l), c_(c) {}void operator()(state_t const& x, state_t& dxdt, double t) {// 计算当前状态的变化率dxdt[0] = ((1 - x[0]) / r_ - x[1]) / c_;dxdt[1] = x[0] / l_;}
private:double r_, l_, c_; // 电路参数，非状态变量
};

这里定义了一个 circuit 结构体，它有以下成员：

• 构造函数：接受电路的参数 $R$（电阻）、$L$（电感）和 $C$（电容），并初始化这些成员变量。
• operator()：该函数定义了电路状态的变化。它接受：
  • x：当前电路状态，其中 x[0] 是电感电流 $I_L$，x[1] 是输出电压 $V_{\text{out}}$。
  • dxdt：用于存储状态变量的变化率（即导数）。
  • t：当前的时间，虽然在这里没有使用，但它通常用于时间相关的计算。

电路的微分方程

根据电路的物理模型，电路的状态变化可以表示为微分方程。假设我们有一个简单的 RLC 电路，其微分方程如下：

1. 电流 $I_L$ 和电压 $V_{\text{out}}$ 的变化：
   • 电感电流的变化率：$\frac{d{I}_L}{dt}=\frac{V_{\text{in}}-V_{\text{out}}}{L}$
   • 输出电压的变化率：$\frac{d{V}_{\text{out}}}{dt}=\frac{I_L}{C}$
     在这个实现中，微分方程被转化成了 dxdt[0] 和 dxdt[1]，表示每个状态变量的变化率。

• dxdt[0] = ((1 - x[0]) / r_ - x[1]) / c_; 表示电感电流 $I_L$ 的变化率，取决于电阻、电容和电感。
• dxdt[1] = x[0] / l_; 表示输出电压 $V_{\text{out}}$ 的变化率，取决于电感电流 $I_L$。

总结

这段代码使用了数值积分的基本概念，将电路的微分方程转化为状态变量的变化率，并通过数值方法（如Boost.ODEInt）来求解系统的时间演化。这种方法可以有效地模拟电路的行为，尤其是在没有简单解析解的情况下。

Numeric Integration with ODEInt - Simulation

int main() {using namespace boost::numeric::odeint;circuit ckt(100.0, 20e-6, 20e-9);state_t x{0.0, 0.0};                   // initial conditionsintegrate(ckt, x, 0.0, 10e-6, 0.1e-6,  // time range// sink for values:[](state_t const& x, double t) { std::cout << t << " " << x[0] << std::endl; });
}

这段代码演示了如何使用 Boost.ODEInt 库进行数值积分模拟。ODEInt 是一个强大的数值微分方程求解器，它能帮助我们通过离散的时间步长来求解电路或物理系统的状态。以下是代码解析：

代码解析

using namespace boost::numeric::odeint;

首先，我们导入了 boost::numeric::odeint 命名空间，这样可以直接使用 ODEInt 库中的函数和类型。

circuit ckt(100.0, 20e-6, 20e-9);

创建了一个 circuit 对象 ckt，并为其提供了三个参数：电阻 $R=100.0$ 欧姆，电感 $L=20e-6$ 亨利，电容 $C=20e-9$ 法拉。这些是电路的基本参数。

state_t x{0.0, 0.0}; // initial conditions

定义了一个状态向量 x，表示电路的初始状态。在这个例子中，x[0] 是电感电流 $I_L$，x[1] 是输出电压 $V_{\text{out}}$。这两个变量都初始化为 0.0，表示在时间 $t=0$ 时电流和电压都为零。

integrate( ckt, x,0.0, 10e-6, 0.1e-6, // time range// sink for values:[](state_t const& x, double t) {std::cout << t << " " << x[0] << std::endl;});

这是 ODEInt 库的核心部分，它用于求解微分方程并对结果进行输出。以下是它的各个部分：

• ckt: 传递给 integrate 函数的是电路对象 ckt，该对象提供了求解微分方程所需的函数（operator()），计算电路状态的变化率。
• x: 这是一个 state_t 类型的变量，存储了电路的初始状态（电感电流和输出电压）。
• 0.0, 10e-6: 这两个参数定义了模拟的时间范围，从 t = 0.0 到 t = 10e-6（即 10 微秒）。
• 0.1e-6: 这是积分的时间步长，即每次积分之间的时间间隔为 0.1e-6（即 0.1 微秒）。通过调节步长，可以控制模拟的精度和运行速度。
• Lambda function (sink for values): 这个匿名函数（lambda 函数）负责在每个时间步打印电感电流 $I_L$ 的值。它接受两个参数：
  • state_t const& x：当前状态的变量，这里是电路的状态向量。
  • double t：当前的时间。
    这个 lambda 函数将每个时间步的电感电流输出到标准输出。

总结

通过这段代码，我们实现了一个简单的数值积分模拟，计算了电路中电感电流 $I_L$ 随时间的变化。ODEInt 库自动处理微分方程的数值求解，并在每个时间步打印计算出的电感电流。

进一步优化

1. 调整步长：可以调整时间步长（例如从 0.1e-6 改为 0.05e-6）以提高精度，或者增加模拟的时间范围来获得更详细的模拟结果。
2. 更多状态变量：除了电流和电压，还可以引入其他状态变量，扩展到更复杂的电路模型。
3. 图形可视化：通过使用图形库（如 matplotlib）将模拟结果绘制成图，方便分析电路的响应。

系统化过程：修改节点分析（MNA）

这个过程的目标是自动化地解决任意电路的微分方程。**修改节点分析（MNA）**是一种广泛使用的技术，能够以数学方式表示电路，它可以应用于包含电阻、电容、电感和电压或电流源等多种类型的电路。

MNA 概述

在 MNA 中，目标是将电路转换成一个 线性方程组，可以通过数值方法求解。关键思路是用一个 状态向量 来表示电路的所有行为信息（电压、电流等）。
电路的通用方程形式为：
$$C\frac{dX}{dt}=GX+u$$
其中：

• $X$ 是 状态向量，包含电路中不同节点的电压（对于电压控制元件如电阻和电容）以及电感器电流或控制源的电流。
• $C$ 是 电容矩阵，代表时间导数项（即电容和电感）。
• $G$ 是 导纳矩阵，代表时间不变的电路元件（电阻、电压源等）。
• $u$ 是 输入向量，包含电路中电压源或电流源的值。

各个组成部分的详细解析：

1. 状态向量（$X$）：
   • 状态向量包含描述电路当前状态的电压和电流。例如，如果我们有节点电压 $V_1,V_2,\dots$ 和电感电流 $I_L$，则状态向量可能是：
     $$X=\left[\begin{array}{c}{V}_1\\ V_2\\ ⋮\\ I_L\end{array}\right]$$
2. 电容矩阵（$C$）：
   • 这个矩阵包含了关于电容和电感的信息，因为它们涉及时间导数。对于电容器，电压随时间变化的导数与电流成正比，类似地，对于电感器，电流随时间变化的导数与电压成正比。
3. 导纳矩阵（$G$）：
   • 这个矩阵包含了 时间不变 的电路元件，例如电阻、电压源和电流源。它们不涉及时间导数。
4. 输入向量（$u$）：
   • 输入向量包含了电压源和电流源的值，这些源驱动着电路的响应。

使用 MNA 的好处

• 自动化：一旦电路被转换为这种标准形式，求解方程组的过程就可以自动化。你不再需要手动重新排列微分方程。
• 可扩展性：MNA 可以适应 任意电路，无论是包含线性元件（电阻、电容、电感）还是非线性元件（如二极管、晶体管）。
• 高效性：使用数值求解器（例如 ODEInt 等库中的矩阵求解器）来求解这组线性方程是计算上高效的。

下一步：

1. 电路表示：将电路元件（电阻、电容、电感等）以矩阵形式表示 $C$、$G$ 和 $u$。电路中的每个元件都可以映射到这些矩阵中的对应项。
2. 求解系统：一旦系统处于标准形式，我们就可以使用 数值积分方法（例如之前讨论过的 Boost.ODEInt）来求解。
3. 用户界面：提供一个接口，用户可以定义任意电路（包括其元件和连接），然后系统自动生成 MNA 所需的矩阵并求解方程。

示例应用场景：

假设用户定义了一个简单的电路，包含一个电阻、电容和电压源：

1. 电阻将贡献到 导纳矩阵 $G$。
2. 电容将贡献到 电容矩阵 $C$。
3. 电压源将被包含在 输入向量 $u$ 中。
   然后，系统将自动生成 MNA 方程组：
   $$C\frac{dX}{dt}=GX+u$$
   并求解这个方程组，模拟电路的行为。
   通过系统化电路仿真过程，我们可以高效地求解和分析电路，而不需要每次手动操作复杂的微分方程。这也为模拟更复杂的电路打开了大门，这些电路是手动处理起来非常困难的。

修改节点分析（Modified Nodal Analysis, MNA）

在之前的讨论中，我们介绍了如何使用 修改节点分析（MNA）来表示电路。接下来，我们将通过一个 示例电路 来更深入地理解这一方法，展示如何将电路的方程转化为标准的矩阵形式，并且通过数值方法求解这些方程。

电路示例

在该示例中，我们考虑一个简单的电路，其中包含：

• 电容 $C$
• 电感 $L$
• 电压源 $Vin(t)$
• 电流源 $Iin$
• 电阻 $R$
  电路的节点矩阵表示如下：
  电路的节点方程（$dX/dt$）将会是：
  $$C\cdot \frac{dX}{dt}=\left[\begin{array}{cccc}0& 0& 0& 0\\ 0& C& 0& 0\\ 0& 0& L& 0\\ 0& 0& 0& 0\end{array}\right]\cdot \frac{d}{dt}\left[\begin{array}{c}{V}_0\\ V_1\\ I_L\\ I_{in}\end{array}\right]$$

$$G\cdot X+u=\left[\begin{array}{cccc}\frac{1}{R}& -\frac{1}{R}& 0& 0\\ -\frac{1}{R}& \frac{1}{R}& 1& 0\\ 0& -1& 0& 0\\ -1& 0& 0& 0\end{array}\right]\cdot \left[\begin{array}{c}{V}_0\\ V_1\\ I_L\\ I_{in}\end{array}\right]+\left[\begin{array}{c}0\\ 0\\ 0\\ -V_{in}(t)\end{array}\right]$$

矩阵形式

• $C$ 是电容矩阵，包含了与电容、电感相关的项。
• $G$ 是导纳矩阵，包含了电阻、电源和电流源的值。
• $u$ 是输入矩阵，包含了电压源的输入 $Vin(t)$。
  这种表示方式使得电路的状态可以通过一个 状态向量 $X=[V_0,V_1,I_L,I_{\text{in}}{]}^T$ 来表示。这种方法使用了矩阵运算来同时处理多个节点和元件，使得我们能够通过数值方法（如Boost.ODEInt）求解电路的行为。

修改节点分析（MNA）的求解过程

我们最终希望得到这样的形式，来满足 Boost.ODEInt 的输入要求：
$$\frac{dX(t)}{dt}=-C^{-1}GX+C^{-1}u$$
这样我们可以通过数值积分方法来求解 $X(t)$，从而获得电路的输出。接下来，我们来看看如何利用矩阵库来实现这些计算。

我们需要一个矩阵库

要表示和计算这个电路模型，我们需要使用一个矩阵库来帮助我们执行矩阵运算。C++ 标准库本身不提供矩阵运算的功能，因此我们需要依赖外部库。
常用的矩阵库：

• Eigen：这是一个非常流行的线性代数库，支持矩阵运算、向量运算和求解方程组。
• Boost uBLAS：Boost 库中的一个模块，支持矩阵和向量的基本操作。
• Armadillo：另一个强大的 C++ 数值线性代数库，常用于机器学习和工程计算。
  在这个例子中，我们可以使用 Boost uBLAS 来构建和操作矩阵。Boost 是 C++ 社区中非常流行的库，它提供了一个高效且简单的接口来执行矩阵运算。

示例代码：MNA 电路仿真

以下是一个简化的代码示例，展示了如何使用 MNA 方程来定义电路，并使用数值方法求解：

#include <boost/numeric/odeint.hpp>
#include <iostream>
#include <array>
// 状态向量类型：包含电压和电流
typedef std::array<double, 4> state_t; // V0, V1, IL, Iin
// 电路结构体
struct circuit {circuit(double r, double l, double c): r_(r), l_(l), c_(c) {}// 求解 dX/dt 方程void operator()(state_t const& x, state_t& dxdt, double t) {// 电流、电压的变化率dxdt[0] = (x[1] - x[0]) / r_ - x[3] / c_; // V0变化dxdt[1] = (x[0] - x[1]) / l_;              // V1变化dxdt[2] = (x[0]) / l_;                     // 电感电流变化dxdt[3] = Vin(t);                          // 输入电流变化}// 电压源函数double Vin(double t) {return 1.0;  // 假设为一个常数电压源，Vin(t) = 1.0V}
private:double r_, l_, c_; // 电阻、电感、电容的参数
};
int main() {using namespace boost::numeric::odeint;circuit ckt(100.0, 20e-6, 20e-9);  // 电阻、电感、电容值state_t x = {0.0, 0.0, 0.0, 0.0};  // 初始状态：电压和电流都为0// 求解时间范围 [0, 10e-6]，时间步长为 0.1e-6integrate(ckt, x, 0.0, 10e-6, 0.1e-6,// 输出每一步的值[](state_t const& x, double t) {std::cout << t << " " << x[0] << std::endl;  // 输出电压 V0});
}

总结：

通过 修改节点分析（MNA），我们将电路的微分方程转化为矩阵形式，并使用数值方法（例如 Boost.ODEInt）求解。这种方法非常强大，因为它可以处理 任意电路，无论电路多么复杂。只要电路的元件和连接信息能被转换为合适的矩阵和状态向量，就可以使用此方法进行模拟。

Eigen 库概述

Eigen 是一个功能强大的矩阵库，使用了 表达式模板（expression templates）来提高运算效率。Eigen 支持以下几个主要特性：

• 稀疏矩阵处理：对于包含大量零元素的矩阵，Eigen 提供了紧凑存储的方式，从而提高内存使用效率。
• 算法支持：例如常见的矩阵分解算法（LU、QR、Cholesky）、特征值计算等。
• 加速功能：利用懒计算（lazy evaluation）和 SIMD（单指令多数据）并行化，Eigen 能有效加速大规模矩阵计算。

Eigen 实现 MNA Stamps

在使用 修改节点分析（MNA）进行电路建模时，我们通常需要生成 MNA 印记（stamps），其中每个元件都能通过一个“印记”来贡献到系统的状态矩阵中。下面我们通过一个电阻器的例子来展示如何使用 Eigen 库来实现这些印记。

电阻器示例

假设我们有一个电阻器连接在电路的两个节点之间。我们可以通过下面的函数为这个电阻器生成印记：

template<int sz>
void stamp_r(Eigen::Matrix<double, sz, sz>& G,Eigen::Matrix<double, sz, sz> const&,int node1, int node2, double r) {G(node1, node1) += 1.0 / r;   // 第一个节点的电流项G(node1, node2) -= 1.0 / r;   // 两个节点之间的电流项G(node2, node2) += 1.0 / r;   // 第二个节点的电流项G(node2, node1) -= 1.0 / r;   // 两个节点之间的电流项
}

在这个函数中，我们将电阻器的贡献添加到 导纳矩阵（$G$） 中。具体而言，电阻器将相应的 导纳（即 $1/R$）添加到电流-电压关系中。
使用示例：

stamp_r(G, C, 0, 1, r);

这表示在矩阵 $G$ 中添加一个连接节点 0 和节点 1 的电阻器，其阻值为 $r$。

处理奇异矩阵 C

在 MNA 中，矩阵 $C$ 代表了电容器和电感器的时间导数项。对于包含电容和电感的电路，矩阵 $C$ 可能是 奇异的，即它的行列式为零。为了使用这些方程，我们需要对矩阵 $C$ 进行正则化。

正则化方法

通过 高斯消元法，我们可以将矩阵 $C$ 转换为一个非奇异的形式，通常我们可以通过 LU 分解来实现：

auto lu = C.fullPivLu(); // 进行 LU 分解
Eigen::MatrixXd Cprime = lu.matrixLU().template triangularView<Eigen::Upper>();  // 获取上三角矩阵

这个操作将矩阵 $C$ 转换为一个更容易求解的形式。接下来，我们可以应用 LU 分解的结果来调整 导纳矩阵 $G$ 和 输入向量 $u$。

简化电路状态

一旦我们对矩阵进行了正则化，我们可以将电路的状态简化为一个包含较少状态变量的系统。具体来说，可以通过消除一些变量（例如电压 $V_0$ 或电流源 $I_{\text{in}}$）来减少问题的维度。
通过以下矩阵操作，我们可以得到一个简化的系统：

// 从 Cprime 和 Gprime 中提取子矩阵，减少系统的维度
Eigen::Matrix2d Cnew = Cprime.topLeftCorner(2, 2);
Eigen::Matrix2d G11 = Gprime.topLeftCorner(2, 2);

通过子矩阵运算，我们可以得到一个简化的导纳矩阵 $G_{\text{new}}$，并继续进行数值求解。

求解简化后的系统

在这个简化的系统中，我们可以通过矩阵反演来求解状态的变化：

Eigen::Matrix2d Gnew = G11 - G12 * G22.inverse() * G21;
Eigen::Matrix2d Cnew_inv = Cnew.inverse();
dxdt = -Cnew_inv * (Gnew * x + Bnew * u);

这个公式将会生成一个简化的 状态变化（$dx/dt$） 方程，从而可以用于数值求解（例如使用 Boost.ODEInt）来模拟电路的行为。

总结

使用 Eigen 库，我们能够高效地处理和求解电路的 修改节点分析（MNA）模型，尤其是在处理 稀疏矩阵 和 奇异矩阵 时，Eigen 提供了强大的工具和算法（如 LU 分解）。通过这些工具，我们可以将电路的状态方程转化为标准的数值积分形式，并通过数值方法进行仿真。这使得我们能够自动化和通用地模拟各种电路，而不需要手动推导和求解每个电路的方程。