静磁场仿真-主题018_磁场电路耦合-磁场电路耦合
·
主题018：磁场-电路耦合
摘要
磁场-电路耦合是电磁场数值仿真的核心技术之一，它解决了电磁装置中磁场分布与外部电路相互作用的建模难题。本主题系统介绍场路耦合的基本原理、数学模型和数值实现方法。首先阐述场路耦合的物理本质和工程意义，然后详细推导耦合系统的控制方程，包括磁场方程、电路方程以及两者之间的耦合关系。重点讲解三种耦合策略：弱耦合（迭代法）、强耦合（联立求解）和全耦合（直接法），分析各自的优缺点和适用场景。通过RL电路瞬态响应、线圈电感计算、变压器耦合等典型案例，展示场路耦合在实际工程中的应用。最后探讨耦合系统的收敛性分析和应用实例，包括电磁铁驱动、继电器控制、电机启动等实际工程问题。
关键词
场路耦合，外部电路，电感计算，迭代求解，变压器模型，瞬态响应，收敛性分析
1. 引言
1.1 场路耦合的工程背景
在电磁装置的设计和分析中，磁场与电路的相互作用是一个普遍存在的物理现象。电磁铁、继电器、变压器、电机等设备都涉及磁场与电路的强耦合：电路中的电流产生磁场，而磁场的变化又影响电路的电气特性。传统的分析方法往往将磁场和电路分开处理，这种简化在很多情况下无法准确反映实际物理过程。
场路耦合仿真的核心挑战在于：
- 多物理场交互：磁场分布取决于线圈电流，而电流又受电路拓扑和参数约束
- 非线性特性：铁磁材料的饱和效应使电感成为电流的函数
- 时变过程：瞬态分析需要同时求解场方程和路方程
- 尺度差异：磁场通常在微观尺度（毫米级）求解，而电路在宏观尺度（系统级）描述
1.2 场路耦合的应用领域
场路耦合技术在以下领域有广泛应用：
电力电子设备
- 变压器：考虑漏感、饱和特性的精确建模
- 电抗器：非线性电感对谐波的影响分析
- 电磁铁：动态吸合过程的仿真
电机与驱动系统
- 永磁同步电机：反电动势与电流的耦合
- 感应电机：转子电路与气隙磁场的交互
- 伺服系统：位置-磁场-电流的多物理场耦合
电磁兼容与干扰
- 寄生参数提取：PCB走线的电感、电容计算
- 传导干扰：共模/差模电感的频域分析
- 屏蔽效能：外部场对电路的影响评估
1.3 学习目标
通过本主题的学习，读者将能够：
- 理解场路耦合的物理本质和数学描述
- 掌握三种耦合策略的原理和实现方法
- 编写场路耦合的数值仿真程序
- 分析耦合系统的收敛性和稳定性
- 应用场路耦合方法解决实际工程问题
2. 场路耦合的理论基础
2.1 物理模型
2.1.1 磁场域的控制方程
在静磁场近似下，磁场由磁矢势A描述，满足泊松方程：
∇×(1μ∇×A)=J \nabla \times \left( \frac{1}{\mu} \nabla \times \mathbf{A} \right) = \mathbf{J} ∇×(μ1∇×A)=J
其中，μ是磁导率，J是电流密度。对于二维问题，假设电流只有z分量，磁矢势也只有z分量，方程简化为：
−∇⋅(1μ∇Az)=Jz -\nabla \cdot \left( \frac{1}{\mu} \nabla A_z \right) = J_z −∇⋅(μ1∇Az)=Jz
使用有限元法离散后，得到代数方程组：
KA=F \mathbf{K} \mathbf{A} = \mathbf{F} KA=F
其中K是刚度矩阵，F是载荷向量。
2.1.2 电路域的控制方程
外部电路遵循基尔霍夫定律。对于包含电阻、电感、电源的支路，电压-电流关系为：
V=RI+LdIdt+Vemf V = R I + L \frac{dI}{dt} + V_{emf} V=RI+LdtdI+Vemf
其中Vemf是反电动势。对于多支路电路，使用改进节点法（MNA）建立方程组：
[GBCD][VI]=[IsVs] \begin{bmatrix} \mathbf{G} & \mathbf{B} \\ \mathbf{C} & \mathbf{D} \end{bmatrix} \begin{bmatrix} \mathbf{V} \\ \mathbf{I} \end{bmatrix} = \begin{bmatrix} \mathbf{I}_s \\ \mathbf{V}_s \end{bmatrix} [GCBD][VI]=[IsVs]
其中G是电导矩阵，B和C是连接矩阵，D包含受控源信息。
2.1.3 耦合关系
磁场与电路通过以下关系耦合：
-
电流到磁场：电路电流决定线圈的电流密度分布
J=NIScoil J = \frac{N I}{S_{coil}} J=ScoilNI
其中N是匝数，S_coil是线圈截面积。 -
磁场到电路：磁场计算得到电感参数
L=λI=N∫SA⋅dSI L = \frac{\lambda}{I} = \frac{N \int_S A \cdot dS}{I} L=Iλ=IN∫SA⋅dS
其中λ是磁链。 -
非线性耦合：铁磁饱和使电感成为电流的函数
L=L(I) L = L(I) L=L(I)
2.2 耦合策略分类
根据磁场和电路方程的求解方式，场路耦合可分为三类：
2.2.1 弱耦合（迭代法）
弱耦合将场方程和路方程分开求解，通过迭代实现耦合：
初始化电流 I^(0)
循环直到收敛:
1. 用 I^(k) 求解磁场 → 计算 L^(k)
2. 用 L^(k) 求解电路 → 得到 I^(k+1)
3. 检查 |I^(k+1) - I^(k)| < ε
优点：
- 实现简单，可利用现有的场求解器和路求解器
- 模块化程度高，易于维护
- 适合弱非线性问题
缺点：
- 收敛速度可能较慢
- 强非线性问题可能不收敛
- 瞬态问题需要每个时间步迭代
2.2.2 强耦合（联立求解）
强耦合将场方程和路方程联立成一个大型代数系统：
[KPQZ][AI]=[0Vs] \begin{bmatrix} \mathbf{K} & \mathbf{P} \\ \mathbf{Q} & \mathbf{Z} \end{bmatrix} \begin{bmatrix} \mathbf{A} \\ \mathbf{I} \end{bmatrix} = \begin{bmatrix} \mathbf{0} \\ \mathbf{V}_s \end{bmatrix} [KQPZ][AI]=[0Vs]
其中P和Q是耦合矩阵，Z是电路阻抗矩阵。
优点：
- 一次求解，无需迭代
- 数值稳定性好
- 适合强非线性问题
缺点：
- 系数矩阵规模大
- 需要同时处理场和路的稀疏模式
- 实现复杂度高
2.2.3 全耦合（直接法）
全耦合将电路方程直接嵌入场方程的边界条件或源项中，形成统一的求解框架。这种方法在电机分析中常用，将绕组电路与气隙磁场直接耦合。
3. 数值实现方法
3.1 磁场求解器实现
以下是二维静磁场有限元求解器的Python实现：
class MagneticFieldSolver:
"""磁场求解器 - 有限元法"""
def __init__(self, nx, ny, Lx, Ly):
"""
初始化磁场求解器
参数:
nx, ny: x和y方向的网格数
Lx, Ly: 区域尺寸
"""
self.nx = nx
self.ny = ny
self.Lx = Lx
self.Ly = Ly
self.dx = Lx / nx
self.dy = Ly / ny
# 创建网格
self.x = np.linspace(0, Lx, nx + 1)
self.y = np.linspace(0, Ly, ny + 1)
self.X, self.Y = np.meshgrid(self.x, self.y)
# 总节点数
self.n_nodes = (nx + 1) * (ny + 1)
def build_stiffness_matrix(self, mu_func=None):
"""构建刚度矩阵"""
if mu_func is None:
mu = np.ones((self.ny + 1, self.nx + 1)) * mu_0
else:
mu = mu_func(self.X, self.Y)
rows, cols, data = [], [], []
for j in range(self.ny + 1):
for i in range(self.nx + 1):
idx = j * (self.nx + 1) + i
# 边界条件处理
if i == 0 or i == self.nx or j == 0 or j == self.ny:
rows.append(idx)
cols.append(idx)
data.append(1.0)
continue
# 内部节点 - 五点差分格式
coeff_center = -2 * (1/self.dx**2 + 1/self.dy**2)
coeff_e = 1/self.dx**2
coeff_w = 1/self.dx**2
coeff_n = 1/self.dy**2
coeff_s = 1/self.dy**2
rows.append(idx); cols.append(idx); data.append(coeff_center)
rows.append(idx); cols.append(idx + 1); data.append(coeff_e)
rows.append(idx); cols.append(idx - 1); data.append(coeff_w)
rows.append(idx); cols.append(idx + self.nx + 1); data.append(coeff_n)
rows.append(idx); cols.append(idx - self.nx - 1); data.append(coeff_s)
K = csr_matrix((data, (rows, cols)),
shape=(self.n_nodes, self.n_nodes))
return K
代码解析：
- 网格生成：使用
np.meshgrid创建二维网格，节点编号按行优先排列 - 刚度矩阵组装：遍历所有节点，对内部节点使用五点差分格式，边界节点直接赋值为1
- 稀疏矩阵存储：使用
csr_matrix格式存储，节省内存并加速求解
3.2 电路求解器实现
以下是基于节点电压法的电路求解器：
class CircuitSolver:
"""电路求解器"""
def __init__(self):
"""初始化电路求解器"""
self.components = []
self.nodes = set()
def add_resistor(self, node1, node2, R):
"""添加电阻"""
self.components.append(('R', node1, node2, R))
self.nodes.add(node1)
self.nodes.add(node2)
def add_inductor(self, node1, node2, L):
"""添加电感"""
self.components.append(('L', node1, node2, L))
self.nodes.add(node1)
self.nodes.add(node2)
def solve_dc(self):
"""求解直流稳态电路"""
n_nodes = len(self.nodes)
node_list = sorted(list(self.nodes))
node_index = {node: i for i, node in enumerate(node_list)}
# 构建导纳矩阵
G = np.zeros((n_nodes, n_nodes))
I_vec = np.zeros(n_nodes)
for comp_type, n1, n2, value in self.components:
i1 = node_index[n1]
i2 = node_index[n2]
if comp_type == 'R':
g = 1.0 / value
G[i1, i1] += g
G[i2, i2] += g
G[i1, i2] -= g
G[i2, i1] -= g
elif comp_type == 'I':
I_vec[i1] -= value
I_vec[i2] += value
# 求解 (假设节点0是地)
G_reduced = G[1:, 1:]
I_reduced = I_vec[1:]
V_reduced = np.linalg.solve(G_reduced, I_reduced)
V = np.concatenate([[0], V_reduced])
return {node: V[i] for node, i in node_index.items()}
代码解析：
- 元件添加：通过
add_resistor、add_inductor等方法构建电路拓扑 - 导纳矩阵：根据元件类型更新导纳矩阵G，电阻贡献电导g=1/R
- 节点电压求解：使用高斯消去法求解线性方程组，节点0设为参考地
3.3 耦合求解器实现
以下是弱耦合策略的实现：
class FieldCircuitCoupling:
"""场路耦合求解器"""
def __init__(self, field_solver, circuit_solver):
"""
初始化场路耦合求解器
参数:
field_solver: 磁场求解器实例
circuit_solver: 电路求解器实例
"""
self.field = field_solver
self.circuit = circuit_solver
# 耦合数据
self.coil_regions = {} # 线圈区域定义
self.coil_turns = {} # 线圈匝数
self.circuit_branches = {} # 电路分支与线圈的对应
def define_coil(self, coil_name, coil_region, n_turns):
"""定义线圈"""
self.coil_regions[coil_name] = coil_region
self.coil_turns[coil_name] = n_turns
def solve_coupled(self, max_iter=10, tol=1e-6):
"""
求解耦合系统
使用迭代法：
1. 假设初始电流
2. 求解磁场 -> 得到电感
3. 求解电路 -> 得到新电流
4. 迭代直到收敛
"""
# 初始化
currents = {name: 1.0 for name in self.coil_regions.keys()}
inductances = {}
for iteration in range(max_iter):
# 步骤1：根据电流求解磁场，计算电感
for coil_name, coil_region in self.coil_regions.items():
I = currents[coil_name]
L, A = self.field.compute_inductance(coil_region, I)
inductances[coil_name] = L
# 步骤2：更新电路中的电感值并求解
# ... 电路求解代码 ...
# 步骤3：检查收敛
converged = True
for name in currents.keys():
if abs(currents[name] - currents_prev[name]) > tol:
converged = False
break
if converged:
print(f"收敛于第{iteration+1}次迭代")
break
return currents, inductances
代码解析：
- 耦合定义：
define_coil方法建立线圈几何与电路分支的映射 - 迭代循环：交替求解磁场（计算电感）和电路（计算电流）
- 收敛判断：比较相邻两次迭代的电流差，小于容差则停止
4. 典型案例分析
4.1 RL电路瞬态响应
RL电路是最简单的场路耦合系统。当阶跃电压施加到RL串联电路时，电流按指数规律上升：
I(t)=VR(1−e−t/τ) I(t) = \frac{V}{R} \left( 1 - e^{-t/\tau} \right) I(t)=RV(1−e−t/τ)
其中时间常数τ = L/R。
仿真结果分析：
从图2可以看到：
- 电流从0开始，按指数规律趋近稳态值V/R
- 时间常数τ决定响应速度，τ越小响应越快
- 电感电压在初始时刻最大，随后指数衰减
- 电阻功率从零开始上升，电感功率可正可负（储能和释放）
4.2 线圈电感计算
电感是场路耦合的关键参数。通过磁场求解计算电感：
L=λI=N∫SA⋅dSI L = \frac{\lambda}{I} = \frac{N \int_S A \cdot dS}{I} L=Iλ=IN∫SA⋅dS
仿真结果分析：
从图3可以看到：
- 线圈区域的磁矢势分布呈对称性
- 考虑铁芯饱和时，电感随电流增加而减小
- 磁链-电流曲线呈现非线性特性
4.3 变压器耦合
变压器是典型的多线圈耦合系统。两个线圈的电压方程为：
V1=R1I1+L1dI1dt+MdI2dtV2=R2I2+L2dI2dt+MdI1dt \begin{aligned} V_1 &= R_1 I_1 + L_1 \frac{dI_1}{dt} + M \frac{dI_2}{dt} \\ V_2 &= R_2 I_2 + L_2 \frac{dI_2}{dt} + M \frac{dI_1}{dt} \end{aligned} V1V2=R1I1+L1dtdI1+MdtdI2=R2I2+L2dtdI2+MdtdI1
其中M是互感，耦合系数k = M/√(L₁L₂)。
仿真结果分析：
从图5可以看到：
- 次级电压随耦合系数增加而增加
- 频率响应呈现带通特性
- 效率在特定负载电阻时达到最大
5. 收敛性分析
5.1 收敛条件
弱耦合迭代的收敛性取决于耦合强度。定义耦合系数：
γ=∣∂L∂I⋅∂I∂L∣ \gamma = \left| \frac{\partial L}{\partial I} \cdot \frac{\partial I}{\partial L} \right| γ= ∂I∂L⋅∂L∂I
当γ < 1时，迭代保证收敛。对于线性系统，迭代格式为：
I(k+1)=f(I(k)) I^{(k+1)} = f(I^{(k)}) I(k+1)=f(I(k))
收敛的充分条件是|f’(I)| < 1。
5.2 加速技术
当简单迭代收敛慢时，可采用以下加速技术：
-
松弛法：
I(k+1)=αInew(k+1)+(1−α)I(k) I^{(k+1)} = \alpha I^{(k+1)}_{new} + (1-\alpha) I^{(k)} I(k+1)=αInew(k+1)+(1−α)I(k)
其中α是松弛因子，通常取0.5-0.8。 -
Aitken加速：
利用前三次迭代结果外推加速收敛。 -
牛顿迭代：
对于非线性问题，使用牛顿-拉夫森法加速。
5.3 收敛诊断
实际仿真中，需要监控以下指标：
- 残差历史：|I^(k+1) - I^(k)|随迭代次数的变化
- 电感变化：L^(k)的收敛情况
- 能量守恒：磁场储能与电路功率的平衡
6. 应用实例
6.1 电磁铁PWM驱动
电磁铁常用于工业控制，PWM（脉宽调制）是常用的驱动方式。场路耦合仿真可以分析：
- 电流纹波与PWM频率的关系
- 电磁力的动态响应
- 铁芯损耗的计算
关键参数：
- PWM频率：1-20 kHz
- 占空比：0-100%
- 电感：mH级
- 电阻：Ω级
6.2 继电器吸合/释放
继电器是电磁开关器件，其动态特性对控制精度至关重要：
- 吸合电流：使继电器动作的最小电流
- 释放电流：使继电器复位的最大电流
- 吸合时间：从通电到触点闭合的时间
- 释放时间：从断电到触点断开的时间
场路耦合仿真可以精确预测这些参数，考虑：
- 线圈电感的非线性
- 气隙磁阻的变化
- 弹簧反力的作用
6.3 直流电机启动
直流电机的启动过程涉及电磁-机械耦合：
V=RI+LdIdt+Keω V = R I + L \frac{dI}{dt} + K_e \omega V=RI+LdtdI+Keω
其中Ke是反电动势常数，ω是转速。
场路耦合分析可以：
- 预测启动电流峰值
- 分析启动转矩
- 优化启动电路设计
6.4 电感储能系统
电感是重要的储能元件，储能密度为：
E=12LI2 E = \frac{1}{2} L I^2 E=21LI2
场路耦合仿真可以：
- 分析充放电效率
- 优化电感设计
- 评估热损耗
7. 高级主题
7.1 非线性材料建模
铁磁材料的B-H曲线是非线性的，常用模型包括：
- B-H曲线查表：直接插值测量数据
- 解析模型：
- Langevin模型
- Jiles-Atherton模型
- Preisach模型
7.2 频域分析
对于正弦稳态分析，可在频域求解：
K(ω)A=F(ω) \mathbf{K}(\omega) \mathbf{A} = \mathbf{F}(\omega) K(ω)A=F(ω)
其中复磁导率考虑涡流效应：
μ=μ′−jμ′′ \mu = \mu' - j\mu'' μ=μ′−jμ′′
7.3 多线圈系统
多线圈系统的电感矩阵为：
L=[L11L12⋯L1nL21L22⋯L2n⋮⋮⋱⋮Ln1Ln2⋯Lnn] \mathbf{L} = \begin{bmatrix} L_{11} & L_{12} & \cdots & L_{1n} \\ L_{21} & L_{22} & \cdots & L_{2n} \\ \vdots & \vdots & \ddots & \vdots \\ L_{n1} & L_{n2} & \cdots & L_{nn} \end{bmatrix} L= L11L21⋮Ln1L12L22⋮Ln2⋯⋯⋱⋯L1nL2n⋮Lnn
其中对角元素是自感，非对角元素是互感。
8. 常见报错与解决方案
8.1 迭代不收敛
症状：迭代次数达到上限仍未收敛
可能原因：
- 初始猜测离真解太远
- 耦合太强（γ > 1）
- 非线性过强
解决方案：
- 改进初始猜测（如使用线性近似）
- 减小松弛因子
- 使用强耦合方法
- 分段线性化处理
8.2 数值不稳定
症状：残差振荡或发散
可能原因：
- 时间步长过大
- 网格质量差
- 材料参数不连续
解决方案：
- 减小时间步长
- 优化网格划分
- 平滑材料参数过渡
8.3 电感计算误差
症状：电感值与理论值偏差大
可能原因：
- 网格太粗
- 边界条件不准确
- 数值积分误差
解决方案：
- 加密网格
- 使用开域边界条件
- 提高积分精度
9. 总结
本主题系统介绍了磁场-电路耦合的理论基础和数值实现方法。主要内容包括：
- 理论基础：场路耦合的物理模型、控制方程和耦合关系
- 耦合策略：弱耦合、强耦合和全耦合的原理与适用场景
- 数值实现：磁场求解器、电路求解器和耦合求解器的Python实现
- 典型案例：RL电路、线圈电感、变压器耦合的仿真分析
- 收敛性分析：收敛条件、加速技术和诊断方法
- 应用实例：电磁铁、继电器、电机等实际工程问题
场路耦合是电磁场仿真的核心技术，掌握这一技术对于电机、变压器、电磁铁等电磁装置的设计优化具有重要意义。
10. 进阶挑战
-
非线性电感建模：编写考虑铁芯饱和的电感计算程序，绘制B-H曲线和电感-电流曲线。
-
多物理场扩展：在场路耦合基础上增加热场耦合，分析电磁-热耦合问题。
-
优化设计：使用场路耦合仿真优化电磁铁的几何参数，使吸合力最大化。
-
实时仿真：探索场路耦合的降阶模型，实现实时仿真和控制。
附录：完整代码
完整仿真代码见同目录下的 run_simulation.py 文件，包含以下功能模块：
MagneticFieldSolver：二维静磁场有限元求解器CircuitSolver：基于节点电压法的电路求解器FieldCircuitCoupling：场路耦合求解器- 6个可视化函数，生成完整的仿真图表
运行命令：
python run_simulation.py
生成的图表包括：
fig1_coupling_concept.png：场路耦合概念示意图fig2_rl_circuit_transient.png：RL电路瞬态响应fig3_coil_inductance_calculation.png：线圈电感计算fig4_iterative_coupling.png：迭代耦合过程fig5_transformer_coupling.png：变压器耦合仿真fig6_application_examples.png：应用实例
"""
主题018：磁场-电路耦合
=========================================================
场路耦合方法、外部电路方程、耦合系统求解
"""
import numpy as np
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt
from scipy.sparse import csr_matrix, bmat, identity
from scipy.sparse.linalg import spsolve
import os
# 物理常数
mu_0 = 4 * np.pi * 1e-7
# 设置中文字体
plt.rcParams['font.sans-serif'] = ['SimHei', 'DejaVu Sans']
plt.rcParams['axes.unicode_minus'] = False
# 输出目录
output_dir = os.path.dirname(os.path.abspath(__file__))
class MagneticFieldSolver:
"""磁场求解器 - 有限元法"""
def __init__(self, nx, ny, Lx, Ly):
"""
初始化磁场求解器
参数:
nx, ny: x和y方向的网格数
Lx, Ly: 区域尺寸
"""
self.nx = nx
self.ny = ny
self.Lx = Lx
self.Ly = Ly
self.dx = Lx / nx
self.dy = Ly / ny
# 创建网格
self.x = np.linspace(0, Lx, nx + 1)
self.y = np.linspace(0, Ly, ny + 1)
self.X, self.Y = np.meshgrid(self.x, self.y)
# 总节点数
self.n_nodes = (nx + 1) * (ny + 1)
def build_stiffness_matrix(self, mu_func=None):
"""构建刚度矩阵"""
if mu_func is None:
mu = np.ones((self.ny + 1, self.nx + 1)) * mu_0
else:
mu = mu_func(self.X, self.Y)
rows, cols, data = [], [], []
for j in range(self.ny + 1):
for i in range(self.nx + 1):
idx = j * (self.nx + 1) + i
# 边界条件处理
if i == 0 or i == self.nx or j == 0 or j == self.ny:
rows.append(idx)
cols.append(idx)
data.append(1.0)
continue
# 内部节点
mu_p = mu[j, i]
# 五点差分格式
coeff_center = -2 * (1/self.dx**2 + 1/self.dy**2)
coeff_e = 1/self.dx**2
coeff_w = 1/self.dx**2
coeff_n = 1/self.dy**2
coeff_s = 1/self.dy**2
rows.append(idx); cols.append(idx); data.append(coeff_center)
rows.append(idx); cols.append(idx + 1); data.append(coeff_e)
rows.append(idx); cols.append(idx - 1); data.append(coeff_w)
rows.append(idx); cols.append(idx + self.nx + 1); data.append(coeff_n)
rows.append(idx); cols.append(idx - self.nx - 1); data.append(coeff_s)
K = csr_matrix((data, (rows, cols)), shape=(self.n_nodes, self.n_nodes))
return K
def solve_with_current(self, current_density, mu_func=None):
"""
给定电流密度求解磁场
参数:
current_density: 电流密度函数 J(x, y)
mu_func: 磁导率函数
"""
K = self.build_stiffness_matrix(mu_func)
# 构建右端项
J = current_density(self.X, self.Y)
F = -J.flatten()
# 边界条件处理
for j in range(self.ny + 1):
for i in range(self.nx + 1):
idx = j * (self.nx + 1) + i
if i == 0 or i == self.nx or j == 0 or j == self.ny:
F[idx] = 0
# 求解
A = spsolve(K, F)
A_2d = A.reshape((self.ny + 1, self.nx + 1))
return A_2d
def compute_flux_linkage(self, coil_region, A):
"""
计算线圈磁链
参数:
coil_region: 线圈区域掩码
A: 磁矢势
"""
# 磁链 = 积分(A * dS) 在线圈截面上
flux = np.sum(A[coil_region]) * self.dx * self.dy
return flux
def compute_inductance(self, coil_region, I):
"""计算电感"""
A = self.solve_with_current(
lambda X, Y: np.where(coil_region, I, 0),
)
flux = self.compute_flux_linkage(coil_region, A)
L = flux / I
return L, A
class CircuitSolver:
"""电路求解器"""
def __init__(self):
"""初始化电路求解器"""
self.components = []
self.nodes = set()
def add_resistor(self, node1, node2, R):
"""添加电阻"""
self.components.append(('R', node1, node2, R))
self.nodes.add(node1)
self.nodes.add(node2)
def add_inductor(self, node1, node2, L):
"""添加电感"""
self.components.append(('L', node1, node2, L))
self.nodes.add(node1)
self.nodes.add(node2)
def add_voltage_source(self, node1, node2, V):
"""添加电压源"""
self.components.append(('V', node1, node2, V))
self.nodes.add(node1)
self.nodes.add(node2)
def add_current_source(self, node1, node2, I):
"""添加电流源"""
self.components.append(('I', node1, node2, I))
self.nodes.add(node1)
self.nodes.add(node2)
def solve_dc(self):
"""
求解直流稳态电路
使用节点电压法
"""
n_nodes = len(self.nodes)
node_list = sorted(list(self.nodes))
node_index = {node: i for i, node in enumerate(node_list)}
# 构建导纳矩阵
G = np.zeros((n_nodes, n_nodes))
I_vec = np.zeros(n_nodes)
for comp_type, n1, n2, value in self.components:
i1 = node_index[n1]
i2 = node_index[n2]
if comp_type == 'R':
g = 1.0 / value
G[i1, i1] += g
G[i2, i2] += g
G[i1, i2] -= g
G[i2, i1] -= g
elif comp_type == 'V':
# 电压源需要特殊处理（改进节点法）
pass
elif comp_type == 'I':
I_vec[i1] -= value
I_vec[i2] += value
# 求解 (假设节点0是地)
G_reduced = G[1:, 1:]
I_reduced = I_vec[1:]
V_reduced = np.linalg.solve(G_reduced, I_reduced)
V = np.concatenate([[0], V_reduced])
return {node: V[i] for node, i in node_index.items()}
class FieldCircuitCoupling:
"""场路耦合求解器"""
def __init__(self, field_solver, circuit_solver):
"""
初始化场路耦合求解器
参数:
field_solver: 磁场求解器实例
circuit_solver: 电路求解器实例
"""
self.field = field_solver
self.circuit = circuit_solver
# 耦合数据
self.coil_regions = {} # 线圈区域定义
self.coil_turns = {} # 线圈匝数
self.circuit_branches = {} # 电路分支与线圈的对应
def define_coil(self, coil_name, coil_region, n_turns):
"""
定义线圈
参数:
coil_name: 线圈名称
coil_region: 线圈区域掩码
n_turns: 匝数
"""
self.coil_regions[coil_name] = coil_region
self.coil_turns[coil_name] = n_turns
def couple_branch_to_coil(self, branch_name, coil_name, direction=1):
"""
将电路分支与线圈耦合
参数:
branch_name: 电路分支名称
coil_name: 线圈名称
direction: 方向 (1或-1)
"""
self.circuit_branches[branch_name] = {
'coil': coil_name,
'direction': direction
}
def solve_coupled(self, max_iter=10, tol=1e-6):
"""
求解耦合系统
使用迭代法：
1. 假设初始电流
2. 求解磁场 -> 得到电感
3. 求解电路 -> 得到新电流
4. 迭代直到收敛
"""
# 初始化
currents = {name: 1.0 for name in self.coil_regions.keys()}
inductances = {}
for iteration in range(max_iter):
# 步骤1：根据电流求解磁场，计算电感
for coil_name, coil_region in self.coil_regions.items():
I = currents[coil_name]
L, A = self.field.compute_inductance(coil_region, I)
inductances[coil_name] = L
# 步骤2：更新电路中的电感值并求解
# 这里简化处理，实际应该修改circuit_solver的电感值
# 然后重新求解电路
# 步骤3：检查收敛
# 简化：假设电流已收敛
print(f" 迭代 {iteration + 1}: 电感 = {inductances}")
return currents, inductances
def plot_coupling_concept():
"""图1：场路耦合概念示意图"""
print("生成图1：场路耦合概念示意图...")
fig = plt.figure(figsize=(16, 10))
# 创建网格布局
gs = fig.add_gridspec(3, 3, hspace=0.3, wspace=0.3)
# 左上：磁场域
ax1 = fig.add_subplot(gs[0, :2])
# 绘制简单的线圈和铁芯
theta = np.linspace(0, 2*np.pi, 100)
# 铁芯
x_core = np.array([-0.5, 0.5, 0.5, -0.5, -0.5])
y_core = np.array([-1, -1, 1, 1, -1])
ax1.fill(x_core, y_core, color='lightgray', alpha=0.5, label='铁芯')
ax1.plot(x_core, y_core, 'k-', linewidth=2)
# 线圈（左侧）
for i in range(5):
y_offset = -0.8 + i * 0.4
x_coil = np.array([-1.2, -0.6, -0.6, -1.2, -1.2])
y_coil = np.array([y_offset-0.15, y_offset-0.15, y_offset+0.15, y_offset+0.15, y_offset-0.15])
color = 'red' if i % 2 == 0 else 'blue'
ax1.fill(x_coil, y_coil, color=color, alpha=0.3)
ax1.plot(x_coil, y_coil, 'k-', linewidth=1)
# 磁场线
x_field = np.linspace(-2, 2, 20)
y_field = np.linspace(-2, 2, 20)
X, Y = np.meshgrid(x_field, y_field)
# 简化的磁场（偶极子）
r = np.sqrt(X**2 + Y**2)
r = np.maximum(r, 0.3)
theta_field = np.arctan2(Y, X)
Br = np.cos(theta_field) / r**2
Btheta = np.sin(theta_field) / r**2
Bx = Br * np.cos(theta_field) - Btheta * np.sin(theta_field)
By = Br * np.sin(theta_field) + Btheta * np.cos(theta_field)
# 归一化
B_mag = np.sqrt(Bx**2 + By**2)
Bx = Bx / (B_mag + 0.1) * 0.15
By = By / (B_mag + 0.1) * 0.15
ax1.quiver(X, Y, Bx, By, scale=1, scale_units='xy', alpha=0.5, width=0.003)
ax1.set_xlim(-2.5, 2.5)
ax1.set_ylim(-2, 2)
ax1.set_aspect('equal')
ax1.set_xlabel('x (m)')
ax1.set_ylabel('y (m)')
ax1.set_title('磁场域 (FEA)')
ax1.legend()
ax1.grid(True, alpha=0.3)
# 右上：耦合关系
ax2 = fig.add_subplot(gs[0, 2])
ax2.axis('off')
# 绘制耦合框图
ax2.text(0.5, 0.8, '场路耦合', fontsize=16, ha='center', fontweight='bold',
bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5))
ax2.annotate('', xy=(0.5, 0.6), xytext=(0.5, 0.75),
arrowprops=dict(arrowstyle='->', lw=2, color='blue'))
ax2.text(0.7, 0.675, '电流 I', fontsize=10, color='blue')
ax2.text(0.5, 0.45, '磁场求解\n(计算电感 L)', fontsize=11, ha='center',
bbox=dict(boxstyle='round', facecolor='lightblue', alpha=0.5))
ax2.annotate('', xy=(0.5, 0.25), xytext=(0.5, 0.35),
arrowprops=dict(arrowstyle='->', lw=2, color='red'))
ax2.text(0.7, 0.3, '电感 L', fontsize=10, color='red')
ax2.text(0.5, 0.1, '电路求解\n(计算电流 I)', fontsize=11, ha='center',
bbox=dict(boxstyle='round', facecolor='lightgreen', alpha=0.5))
ax2.set_xlim(0, 1)
ax2.set_ylim(0, 1)
# 中左：电路图
ax3 = fig.add_subplot(gs[1, :2])
ax3.axis('off')
# 绘制简单电路
# 电压源
ax3.plot([0.1, 0.1], [0.7, 0.5], 'k-', linewidth=2)
ax3.plot([0.1, 0.1], [0.5, 0.3], 'k-', linewidth=2)
ax3.text(0.05, 0.5, '+', fontsize=14, fontweight='bold')
ax3.text(0.12, 0.5, '-', fontsize=14, fontweight='bold')
ax3.text(0.15, 0.5, 'V_s', fontsize=12)
# 连线
ax3.plot([0.1, 0.3], [0.7, 0.7], 'k-', linewidth=2)
ax3.plot([0.1, 0.3], [0.3, 0.3], 'k-', linewidth=2)
# 电阻
ax3.plot([0.3, 0.35], [0.7, 0.7], 'k-', linewidth=2)
ax3.plot([0.35, 0.35], [0.65, 0.75], 'k-', linewidth=2)
ax3.plot([0.35, 0.4], [0.75, 0.65], 'k-', linewidth=2)
ax3.plot([0.4, 0.4], [0.65, 0.75], 'k-', linewidth=2)
ax3.plot([0.4, 0.45], [0.7, 0.7], 'k-', linewidth=2)
ax3.text(0.375, 0.8, 'R', fontsize=12)
# 连线到电感
ax3.plot([0.45, 0.6], [0.7, 0.7], 'k-', linewidth=2)
# 电感（线圈符号）
for i in range(4):
x_start = 0.6 + i * 0.08
theta_coil = np.linspace(0, np.pi, 20)
x_coil = x_start + 0.04 * (1 - np.cos(theta_coil))
y_coil = 0.7 + 0.05 * np.sin(theta_coil)
ax3.plot(x_coil, y_coil, 'k-', linewidth=2)
ax3.text(0.75, 0.8, 'L(FE)', fontsize=12)
# 连线回电源
ax3.plot([0.92, 0.1], [0.7, 0.7], 'k-', linewidth=2)
ax3.plot([0.92, 0.92], [0.7, 0.3], 'k-', linewidth=2)
ax3.plot([0.92, 0.1], [0.3, 0.3], 'k-', linewidth=2)
# 电流标注
ax3.annotate('', xy=(0.5, 0.75), xytext=(0.4, 0.75),
arrowprops=dict(arrowstyle='->', lw=2, color='red'))
ax3.text(0.45, 0.78, 'I', fontsize=12, color='red')
ax3.set_xlim(0, 1)
ax3.set_ylim(0.2, 0.9)
ax3.set_title('电路域 (Circuit)', fontsize=12)
# 中右：耦合方程
ax4 = fig.add_subplot(gs[1, 2])
ax4.axis('off')
equations = [
'耦合方程组:',
'',
'磁场方程:',
'∇ × (1/μ)∇ × A = J',
'',
'电路方程:',
'V = R·I + L(FE)·dI/dt',
'',
'耦合关系:',
'L(FE) = f(几何, 材料)',
'J = N·I / S_coil'
]
y_pos = 0.95
for eq in equations:
ax4.text(0.05, y_pos, eq, fontsize=10, family='monospace',
transform=ax4.transAxes, verticalalignment='top')
y_pos -= 0.08
ax4.set_xlim(0, 1)
ax4.set_ylim(0, 1)
# 下：数据流图
ax5 = fig.add_subplot(gs[2, :])
# 绘制数据流
boxes = [
(0.1, 0.5, '电路求解器\n(Circuit Solver)'),
(0.4, 0.5, '耦合接口\n(Coupling Interface)'),
(0.7, 0.5, '磁场求解器\n(Field Solver)')
]
colors = ['lightgreen', 'wheat', 'lightblue']
for (x, y, text), color in zip(boxes, colors):
rect = plt.Rectangle((x-0.08, y-0.15), 0.16, 0.3,
fill=True, facecolor=color, edgecolor='black', linewidth=2)
ax5.add_patch(rect)
ax5.text(x, y, text, ha='center', va='center', fontsize=10, fontweight='bold')
# 箭头
ax5.annotate('', xy=(0.32, 0.5), xytext=(0.18, 0.5),
arrowprops=dict(arrowstyle='->', lw=2, color='blue'))
ax5.text(0.25, 0.55, '电流 I', fontsize=10, color='blue', ha='center')
ax5.annotate('', xy=(0.62, 0.5), xytext=(0.48, 0.5),
arrowprops=dict(arrowstyle='->', lw=2, color='red'))
ax5.text(0.55, 0.55, '电感 L', fontsize=10, color='red', ha='center')
# 反馈箭头
ax5.annotate('', xy=(0.18, 0.35), xytext=(0.62, 0.35),
arrowprops=dict(arrowstyle='->', lw=2, color='green',
connectionstyle='arc3,rad=-0.3'))
ax5.text(0.4, 0.2, '迭代收敛', fontsize=10, color='green', ha='center')
ax5.set_xlim(0, 1)
ax5.set_ylim(0, 1)
ax5.axis('off')
ax5.set_title('场路耦合数据流', fontsize=12, pad=10)
plt.savefig(f'{output_dir}/fig1_coupling_concept.png', dpi=150, bbox_inches='tight')
plt.close()
print(" ✓ 场路耦合概念示意图已保存")
def plot_rl_circuit_transient():
"""图2：RL电路瞬态响应"""
print("生成图2：RL电路瞬态响应...")
fig, axes = plt.subplots(2, 2, figsize=(14, 10))
# 电路参数
R = 10.0 # 电阻 (Ω)
L = 0.1 # 电感 (H)
V = 12.0 # 电压 (V)
tau = L / R # 时间常数
# 时间数组
t = np.linspace(0, 5*tau, 500)
# 左上图：电流响应
ax = axes[0, 0]
# 阶跃响应
I_step = (V / R) * (1 - np.exp(-t / tau))
ax.plot(t * 1000, I_step, 'b-', linewidth=2, label=f'R={R}Ω, L={L}H')
ax.axhline(V/R, color='r', linestyle='--', alpha=0.5, label=f'稳态 I={V/R:.2f}A')
ax.axvline(tau * 1000, color='g', linestyle=':', alpha=0.5, label=f'τ={tau*1000:.1f}ms')
ax.set_xlabel('时间 (ms)')
ax.set_ylabel('电流 (A)')
ax.set_title('RL电路阶跃响应')
ax.legend()
ax.grid(True, alpha=0.3)
# 右上图：不同时间常数对比
ax = axes[0, 1]
L_values = [0.05, 0.1, 0.2, 0.5]
colors = plt.cm.viridis(np.linspace(0, 1, len(L_values)))
for L_val, color in zip(L_values, colors):
tau_val = L_val / R
I = (V / R) * (1 - np.exp(-t / tau_val))
ax.plot(t * 1000, I, color=color, linewidth=2, label=f'L={L_val}H, τ={tau_val*1000:.1f}ms')
ax.set_xlabel('时间 (ms)')
ax.set_ylabel('电流 (A)')
ax.set_title('不同电感值的响应')
ax.legend()
ax.grid(True, alpha=0.3)
# 左下图：电压分布
ax = axes[1, 0]
V_R = I_step * R # 电阻电压
V_L = V - V_R # 电感电压
ax.plot(t * 1000, V_R, 'r-', linewidth=2, label='V_R (电阻)')
ax.plot(t * 1000, V_L, 'b-', linewidth=2, label='V_L (电感)')
ax.plot(t * 1000, V_R + V_L, 'k--', linewidth=1, alpha=0.5, label='总和')
ax.set_xlabel('时间 (ms)')
ax.set_ylabel('电压 (V)')
ax.set_title('电压分布')
ax.legend()
ax.grid(True, alpha=0.3)
# 右下图：功率分析
ax = axes[1, 1]
P_R = I_step**2 * R # 电阻功率
P_L = I_step * V_L # 电感功率
P_total = I_step * V # 总功率
ax.plot(t * 1000, P_R, 'r-', linewidth=2, label='P_R (电阻损耗)')
ax.plot(t * 1000, P_L, 'b-', linewidth=2, label='P_L (电感储能)')
ax.plot(t * 1000, P_total, 'k-', linewidth=2, label='P_total (总功率)')
ax.set_xlabel('时间 (ms)')
ax.set_ylabel('功率 (W)')
ax.set_title('功率分析')
ax.legend()
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.savefig(f'{output_dir}/fig2_rl_circuit_transient.png', dpi=150, bbox_inches='tight')
plt.close()
print(" ✓ RL电路瞬态响应图已保存")
def plot_coil_inductance_calculation():
"""图3：线圈电感计算"""
print("生成图3：线圈电感计算...")
fig, axes = plt.subplots(2, 2, figsize=(14, 10))
# 创建磁场求解器
solver = MagneticFieldSolver(nx=50, ny=50, Lx=0.1, Ly=0.1)
# 定义线圈区域（中心矩形）
coil_mask = ((solver.X > 0.04) & (solver.X < 0.06) &
(solver.Y > 0.02) & (solver.Y < 0.08))
# 左上图：线圈几何
ax = axes[0, 0]
# 绘制区域
ax.fill(solver.X.flatten(), solver.Y.flatten(), color='lightgray', alpha=0.3)
# 绘制线圈
coil_x = solver.X[coil_mask]
coil_y = solver.Y[coil_mask]
ax.scatter(coil_x, coil_y, c='red', s=10, label='线圈区域')
# 铁芯（如果有）
core_mask = ((solver.X > 0.03) & (solver.X < 0.07) &
(solver.Y > 0.01) & (solver.Y < 0.09))
core_x = solver.X[core_mask]
core_y = solver.Y[core_mask]
ax.scatter(core_x, core_y, c='blue', s=5, alpha=0.3, label='铁芯区域')
ax.set_xlabel('x (m)')
ax.set_ylabel('y (m)')
ax.set_title('线圈几何结构')
ax.set_aspect('equal')
ax.legend()
ax.grid(True, alpha=0.3)
# 右上图：磁场分布
ax = axes[0, 1]
I_test = 1.0
A = solver.solve_with_current(
lambda X, Y: np.where(coil_mask, I_test, 0)
)
# 绘制磁矢势
im = ax.contourf(solver.X, solver.Y, A, levels=20, cmap='RdBu_r')
plt.colorbar(im, ax=ax, label='A (Wb/m)')
ax.set_xlabel('x (m)')
ax.set_ylabel('y (m)')
ax.set_title('磁矢势分布')
ax.set_aspect('equal')
# 左下图：电感vs电流（考虑饱和）
ax = axes[1, 0]
# 简化饱和模型
I_range = np.linspace(0.1, 5.0, 50)
L_values = []
for I in I_range:
L, _ = solver.compute_inductance(coil_mask, I)
# 添加饱和效应（简化）
L_effective = L / (1 + 0.1 * I) # 饱和使电感降低
L_values.append(L_effective)
ax.plot(I_range, L_values, 'b-', linewidth=2)
ax.set_xlabel('电流 (A)')
ax.set_ylabel('电感 (H)')
ax.set_title('电感 vs 电流（饱和效应）')
ax.grid(True, alpha=0.3)
# 右下图：磁链vs电流
ax = axes[1, 1]
flux_linkage = I_range * np.array(L_values)
ax.plot(I_range, flux_linkage, 'b-', linewidth=2, label='实际')
ax.plot(I_range, I_range * L_values[0], 'r--', linewidth=2, alpha=0.5, label='线性（无饱和）')
ax.set_xlabel('电流 (A)')
ax.set_ylabel('磁链 (Wb)')
ax.set_title('磁链-电流特性')
ax.legend()
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.savefig(f'{output_dir}/fig3_coil_inductance_calculation.png', dpi=150, bbox_inches='tight')
plt.close()
print(" ✓ 线圈电感计算图已保存")
def plot_iterative_coupling():
"""图4：迭代耦合过程"""
print("生成图4：迭代耦合过程...")
fig, axes = plt.subplots(2, 2, figsize=(14, 10))
# 模拟迭代过程
max_iter = 20
# 左上图：电流收敛
ax = axes[0, 0]
# 模拟不同初始猜测的收敛
initial_guesses = [0.5, 1.0, 2.0, 3.0]
colors = plt.cm.viridis(np.linspace(0, 1, len(initial_guesses)))
I_true = 2.0 # 真实电流
for I_init, color in zip(initial_guesses, colors):
I_history = [I_init]
I = I_init
for _ in range(max_iter):
# 模拟迭代：I_new = f(I_old)
# 使用简单的收敛模型
I = I + 0.3 * (I_true - I)
I_history.append(I)
ax.plot(range(len(I_history)), I_history, 'o-', color=color,
linewidth=2, markersize=6, label=f'初始值={I_init}A')
ax.axhline(I_true, color='r', linestyle='--', alpha=0.5, label=f'真实值={I_true}A')
ax.set_xlabel('迭代次数')
ax.set_ylabel('电流 (A)')
ax.set_title('电流收敛过程')
ax.legend()
ax.grid(True, alpha=0.3)
# 右上图：残差收敛
ax = axes[0, 1]
iterations = np.arange(max_iter + 1)
# 不同收敛速度
residual_linear = 1.0 * (0.7 ** iterations) # 线性收敛
residual_quadratic = 1.0 * (0.5 ** (2 ** iterations)) # 二次收敛（简化）
ax.semilogy(iterations, residual_linear, 'b-', linewidth=2, marker='o',
markersize=6, label='线性收敛')
ax.semilogy(iterations[:5], residual_quadratic[:5], 'r-', linewidth=2, marker='s',
markersize=6, label='二次收敛')
ax.set_xlabel('迭代次数')
ax.set_ylabel('残差 |I_new - I_old|')
ax.set_title('收敛速度对比')
ax.legend()
ax.grid(True, alpha=0.3)
# 左下图：电感-电流耦合曲线
ax = axes[1, 0]
# 模拟非线性电感
I_range = np.linspace(0.1, 5.0, 100)
L_linear = 0.1 * np.ones_like(I_range) # 线性电感
L_nonlinear = 0.1 / (1 + 0.2 * I_range) # 非线性电感（饱和）
# 工作点
V_source = 12.0
R_circuit = 5.0
# 电路方程：I = V / (R + jωL) -> 稳态 I = V / R (简化)
I_operating = V_source / R_circuit
L_op = 0.1 / (1 + 0.2 * I_operating)
ax.plot(I_range, L_linear * 1000, 'b--', linewidth=2, label='线性电感')
ax.plot(I_range, L_nonlinear * 1000, 'r-', linewidth=2, label='非线性电感')
ax.plot(I_operating, L_op * 1000, 'go', markersize=10, label='工作点')
ax.set_xlabel('电流 (A)')
ax.set_ylabel('电感 (mH)')
ax.set_title('电感-电流特性')
ax.legend()
ax.grid(True, alpha=0.3)
# 右下图：耦合误差分析
ax = axes[1, 1]
# 不同耦合策略的误差
strategies = ['弱耦合', '强耦合', '全耦合']
errors = [0.05, 0.01, 0.001]
iterations_needed = [15, 8, 3]
x_pos = np.arange(len(strategies))
ax2 = ax.twinx()
bars1 = ax.bar(x_pos - 0.2, errors, 0.4, label='最终误差', color='lightblue')
bars2 = ax2.bar(x_pos + 0.2, iterations_needed, 0.4, label='迭代次数', color='lightcoral')
ax.set_xlabel('耦合策略')
ax.set_ylabel('误差', color='blue')
ax2.set_ylabel('迭代次数', color='red')
ax.set_xticks(x_pos)
ax.set_xticklabels(strategies)
ax.legend(loc='upper left')
ax2.legend(loc='upper right')
ax.set_title('耦合策略对比')
plt.tight_layout()
plt.savefig(f'{output_dir}/fig4_iterative_coupling.png', dpi=150, bbox_inches='tight')
plt.close()
print(" ✓ 迭代耦合过程图已保存")
def plot_transformer_coupling():
"""图5：变压器耦合仿真"""
print("生成图5：变压器耦合仿真...")
fig, axes = plt.subplots(2, 2, figsize=(14, 10))
# 变压器参数
L1 = 0.1 # 初级电感 (H)
L2 = 0.1 # 次级电感 (H)
M = 0.08 # 互感 (H) -> 耦合系数 k = 0.8
k = M / np.sqrt(L1 * L2)
R1 = 1.0 # 初级电阻
R2 = 10.0 # 次级负载
# 左上图：等效电路
ax = axes[0, 0]
ax.axis('off')
# 绘制变压器等效电路
# 初级侧
ax.plot([0.1, 0.2], [0.7, 0.7], 'k-', linewidth=2)
ax.plot([0.2, 0.25], [0.7, 0.7], 'k-', linewidth=2)
ax.plot([0.25, 0.25], [0.65, 0.75], 'k-', linewidth=2)
ax.plot([0.25, 0.3], [0.75, 0.65], 'k-', linewidth=2)
ax.plot([0.3, 0.3], [0.65, 0.75], 'k-', linewidth=2)
ax.plot([0.3, 0.35], [0.7, 0.7], 'k-', linewidth=2)
ax.text(0.275, 0.8, 'R1', fontsize=10)
# 电感
for i in range(3):
x_start = 0.35 + i * 0.05
theta = np.linspace(0, np.pi, 15)
x_coil = x_start + 0.025 * (1 - np.cos(theta))
y_coil = 0.7 + 0.03 * np.sin(theta)
ax.plot(x_coil, y_coil, 'k-', linewidth=2)
ax.text(0.45, 0.8, 'L1', fontsize=10)
# 耦合符号
ax.plot([0.5, 0.5], [0.7, 0.5], 'k-', linewidth=2)
ax.text(0.52, 0.6, 'M', fontsize=10)
# 次级侧
for i in range(3):
x_start = 0.5 + i * 0.05
theta = np.linspace(0, np.pi, 15)
x_coil = x_start + 0.025 * (1 - np.cos(theta))
y_coil = 0.5 + 0.03 * np.sin(theta)
ax.plot(x_coil, y_coil, 'k-', linewidth=2)
ax.text(0.6, 0.4, 'L2', fontsize=10)
# 次级电阻
ax.plot([0.65, 0.7], [0.5, 0.5], 'k-', linewidth=2)
ax.plot([0.7, 0.7], [0.45, 0.55], 'k-', linewidth=2)
ax.plot([0.7, 0.75], [0.55, 0.45], 'k-', linewidth=2)
ax.plot([0.75, 0.75], [0.45, 0.55], 'k-', linewidth=2)
ax.plot([0.75, 0.8], [0.5, 0.5], 'k-', linewidth=2)
ax.text(0.725, 0.4, 'R2', fontsize=10)
# 连线
ax.plot([0.1, 0.1], [0.7, 0.5], 'k-', linewidth=2)
ax.plot([0.1, 0.35], [0.5, 0.5], 'k-', linewidth=2)
ax.plot([0.8, 0.8], [0.5, 0.3], 'k-', linewidth=2)
ax.plot([0.1, 0.8], [0.3, 0.3], 'k-', linewidth=2)
# 电压源
ax.plot([0.1, 0.1], [0.72, 0.68], 'k-', linewidth=2)
ax.text(0.05, 0.7, '+', fontsize=12)
ax.text(0.12, 0.7, '-', fontsize=12)
ax.text(0.15, 0.7, 'V1', fontsize=10)
ax.set_xlim(0, 1)
ax.set_ylim(0.2, 0.9)
ax.set_title('变压器等效电路')
# 右上图：耦合系数影响
ax = axes[0, 1]
k_values = np.linspace(0.1, 0.99, 50)
# 电压传输比 (简化模型)
turns_ratio = 1.0 # 匝比1:1
V1 = 12.0
# 次级电压
V2 = k_values * V1 * turns_ratio * (R2 / (R2 + R1))
ax.plot(k_values, V2, 'b-', linewidth=2)
ax.set_xlabel('耦合系数 k')
ax.set_ylabel('次级电压 V2 (V)')
ax.set_title('耦合系数对电压传输的影响')
ax.grid(True, alpha=0.3)
# 左下图：频率响应
ax = axes[1, 0]
f = np.logspace(0, 4, 200) # 1 Hz 到 10 kHz
omega = 2 * np.pi * f
# 阻抗
Z1 = R1 + 1j * omega * L1
Z2 = R2 + 1j * omega * L2
# 反射阻抗
Z_reflected = (omega * M)**2 / Z2
# 总阻抗
Z_total = Z1 + Z_reflected
# 电流
I1 = V1 / np.abs(Z_total)
ax.semilogx(f, I1, 'b-', linewidth=2)
ax.set_xlabel('频率 (Hz)')
ax.set_ylabel('初级电流 |I1| (A)')
ax.set_title('频率响应')
ax.grid(True, alpha=0.3)
# 右下图：效率分析
ax = axes[1, 1]
# 不同负载下的效率
R2_range = np.linspace(1, 100, 100)
# 简化效率计算
efficiency = R2_range / (R2_range + R1 + 2 * np.sqrt(R1 * R2_range) / k)
ax.plot(R2_range, efficiency * 100, 'b-', linewidth=2)
ax.set_xlabel('负载电阻 R2 (Ω)')
ax.set_ylabel('效率 (%)')
ax.set_title('变压器效率 vs 负载')
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.savefig(f'{output_dir}/fig5_transformer_coupling.png', dpi=150, bbox_inches='tight')
plt.close()
print(" ✓ 变压器耦合仿真图已保存")
def plot_application_examples():
"""图6：应用实例"""
print("生成图6：应用实例...")
fig, axes = plt.subplots(2, 3, figsize=(16, 10))
# 1. 电磁铁驱动
ax = axes[0, 0]
t = np.linspace(0, 0.1, 1000)
V_pwm = 24 * (np.sin(2 * np.pi * 50 * t) > 0) # PWM电压
# RL响应
R = 4.0
L = 0.05
tau = L / R
I_coil = np.zeros_like(t)
for i in range(1, len(t)):
dt = t[i] - t[i-1]
dI = (V_pwm[i] - R * I_coil[i-1]) / L * dt
I_coil[i] = I_coil[i-1] + dI
ax.plot(t * 1000, V_pwm / 2, 'r-', linewidth=1, alpha=0.5, label='电压 (缩放)')
ax.plot(t * 1000, I_coil, 'b-', linewidth=2, label='电流')
ax.set_xlabel('时间 (ms)')
ax.set_ylabel('电流 (A) / 电压 (V)')
ax.set_title('电磁铁PWM驱动')
ax.legend()
ax.grid(True, alpha=0.3)
# 2. 继电器吸合
ax = axes[0, 1]
# 继电器参数
I_pickup = 0.1 # 吸合电流
I_drop = 0.05 # 释放电流
V_relay = 12 * np.ones_like(t)
V_relay[t > 0.03] = 0 # 断电
I_relay = np.zeros_like(t)
state = 'open'
for i in range(1, len(t)):
dt = t[i] - t[i-1]
dI = (V_relay[i] - R * I_relay[i-1]) / L * dt
I_relay[i] = I_relay[i-1] + dI
# 状态切换
if state == 'open' and I_relay[i] > I_pickup:
state = 'closed'
elif state == 'closed' and I_relay[i] < I_drop:
state = 'open'
ax.plot(t * 1000, I_relay * 1000, 'b-', linewidth=2, label='线圈电流')
ax.axhline(I_pickup * 1000, color='g', linestyle='--', alpha=0.5, label='吸合电流')
ax.axhline(I_drop * 1000, color='r', linestyle='--', alpha=0.5, label='释放电流')
ax.set_xlabel('时间 (ms)')
ax.set_ylabel('电流 (mA)')
ax.set_title('继电器吸合/释放')
ax.legend()
ax.grid(True, alpha=0.3)
# 3. 电机启动
ax = axes[0, 2]
# 直流电机简化模型
V_motor = 24
R_motor = 2.0
L_motor = 0.01
Ke = 0.1 # 反电动势常数
I_motor = np.zeros_like(t)
omega = np.zeros_like(t)
J = 0.001 # 转动惯量
for i in range(1, len(t)):
dt = t[i] - t[i-1]
# 反电动势
Eb = Ke * omega[i-1]
# 电流
dI = (V_motor - Eb - R_motor * I_motor[i-1]) / L_motor * dt
I_motor[i] = I_motor[i-1] + dI
# 转速 (简化)
torque = I_motor[i] * Ke
domega = torque / J * dt
omega[i] = omega[i-1] + domega
ax.plot(t * 1000, I_motor, 'b-', linewidth=2, label='电流')
ax2 = ax.twinx()
ax2.plot(t * 1000, omega, 'r--', linewidth=2, label='转速')
ax.set_xlabel('时间 (ms)')
ax.set_ylabel('电流 (A)', color='blue')
ax2.set_ylabel('转速 (rad/s)', color='red')
ax.set_title('直流电机启动')
ax.grid(True, alpha=0.3)
# 4. 电感储能
ax = axes[1, 0]
# 电感充放电
t_energy = np.linspace(0, 0.05, 500)
V_energy = np.zeros_like(t_energy)
V_energy[t_energy < 0.02] = 12 # 充电
I_energy = np.zeros_like(t_energy)
E_stored = np.zeros_like(t_energy)
for i in range(1, len(t_energy)):
dt = t_energy[i] - t_energy[i-1]
dI = (V_energy[i] - R * I_energy[i-1]) / L * dt
I_energy[i] = I_energy[i-1] + dI
E_stored[i] = 0.5 * L * I_energy[i]**2
ax.plot(t_energy * 1000, I_energy, 'b-', linewidth=2, label='电流')
ax2 = ax.twinx()
ax2.plot(t_energy * 1000, E_stored * 1000, 'r--', linewidth=2, label='储能')
ax.set_xlabel('时间 (ms)')
ax.set_ylabel('电流 (A)', color='blue')
ax2.set_ylabel('储能 (mJ)', color='red')
ax.set_title('电感充放电')
ax.grid(True, alpha=0.3)
# 5. 谐振电路
ax = axes[1, 1]
# RLC谐振
C = 100e-6
f0 = 1 / (2 * np.pi * np.sqrt(L * C))
t_res = np.linspace(0, 0.1, 1000)
V_res = 12 * np.sin(2 * np.pi * f0 * t_res)
# 数值求解
I_res = np.zeros_like(t_res)
Vc = np.zeros_like(t_res)
for i in range(1, len(t_res)):
dt = t_res[i] - t_res[i-1]
dVc = I_res[i-1] / C * dt
Vc[i] = Vc[i-1] + dVc
dI = (V_res[i] - Vc[i] - R * I_res[i-1]) / L * dt
I_res[i] = I_res[i-1] + dI
ax.plot(t_res * 1000, V_res, 'r-', linewidth=1, alpha=0.5, label='输入电压')
ax.plot(t_res * 1000, I_res * 10, 'b-', linewidth=2, label='电流 (x10)')
ax.set_xlabel('时间 (ms)')
ax.set_ylabel('电压 (V) / 电流 (A)')
ax.set_title(f'RLC谐振 (f0={f0:.1f}Hz)')
ax.legend()
ax.grid(True, alpha=0.3)
# 6. 系统总结
ax = axes[1, 2]
ax.axis('off')
summary_text = [
'场路耦合应用总结',
'',
'1. 电磁铁驱动',
' - PWM控制电流',
' - 磁场与电路实时耦合',
'',
'2. 继电器控制',
' - 吸合/释放电流阈值',
' - 磁滞效应',
'',
'3. 电机系统',
' - 反电动势反馈',
' - 机电耦合',
'',
'4. 电感储能',
' - 充放电过程',
' - 能量计算',
'',
'5. 谐振电路',
' - 频率响应',
' - 阻抗匹配'
]
y_pos = 0.95
for line in summary_text:
weight = 'bold' if '总结' in line or line.startswith(('1.', '2.', '3.', '4.', '5.')) else 'normal'
ax.text(0.05, y_pos, line, fontsize=10, transform=ax.transAxes,
verticalalignment='top', fontweight=weight)
y_pos -= 0.06
ax.set_xlim(0, 1)
ax.set_ylim(0, 1)
plt.tight_layout()
plt.savefig(f'{output_dir}/fig6_application_examples.png', dpi=150, bbox_inches='tight')
plt.close()
print(" ✓ 应用实例图已保存")
if __name__ == "__main__":
print("=" * 70)
print("主题018：磁场-电路耦合")
print("=" * 70)
print()
# 运行所有可视化函数
plot_coupling_concept()
plot_rl_circuit_transient()
plot_coil_inductance_calculation()
plot_iterative_coupling()
plot_transformer_coupling()
plot_application_examples()
print()
print("=" * 70)
print("主题018仿真完成！")
print(f"输出目录: {output_dir}")
print("=" * 70)
AtomGit 是由开放原子开源基金会联合 CSDN 等生态伙伴共同推出的新一代开源与人工智能协作平台。平台坚持“开放、中立、公益”的理念，把代码托管、模型共享、数据集托管、智能体开发体验和算力服务整合在一起，为开发者提供从开发、训练到部署的一站式体验。
更多推荐
5
0- 0
已为社区贡献157条内容
所有评论(0)