如何使用scipy.integrate.quad求解带电线段电场无量纲量并绘图

1. 无量纲积分表达式推导 #

取带电线段上距离中心位置为$x$的电荷微元$dx$，微元带电量$dq=\lambda dx$，由于对称性仅需保留电场竖直分量：

$dE = \frac{k dq \cdot z}{(x^2 + z²⁾{3/2}}$
其中$k=\frac{1}{4\pi\varepsilon_0}$，总电场为微元电场在整条线段的积分，积分区间$x\in[-L, L]$：
$E = \int_{-L}^{L} \frac{k\lambda z}{(x^2 + z²⁾{3/2}} dx$
引入无量纲变量：

• $\xi = \frac{z}{L}$（即要求的$Z/L$自变量）
• $u = \frac{x}{L}$，替换后$dx = L du$，积分区间变为$u\in[-1,1]$
  代入原式化简后可得无量纲量$\frac{EL}{k\lambda}$的积分表达式：

$\frac{EL}{k\lambda} = \xi \int_{-1}^{1} \frac{1}{(u^2 + \xi²⁾{3/2}} du$
也可利用积分对称性简化为$\frac{EL}{k\lambda} = 2\xi \int_{0}^{1} \frac{1}{(u^2 + \xi²⁾{3/2}} du$，计算效率更高。

2. Python 实现代码 #

import numpy as np
from scipy.integrate import quad
import matplotlib.pyplot as plt

def calc_dimensionless_E(xi):
    # 定义被积函数
    integrand = lambda u: 1 / (u**2 + xi**2)**(3/2)
    # quad返回第一个值为积分结果，第二个为计算误差
    integral_res, _ = quad(integrand, -1, 1)
    return xi * integral_res

# 生成ξ从0.5到10的采样点，共200个
xi_list = np.linspace(0.5, 10.0, 200)
E_list = [calc_dimensionless_E(xi) for xi in xi_list]

# 绘图
plt.figure(figsize=(8,5))
plt.plot(xi_list, E_list, 'b-', linewidth=2)
plt.xlabel(r'$\xi = z/L$', fontsize=12)
plt.ylabel(r'$EL/k\lambda$', fontsize=12)
plt.title('无量纲电场随z/L变化曲线', fontsize=14)
plt.grid(alpha=0.3)
plt.show()

3. 结果验证 #

该积分存在解析解$\frac{EL}{k\lambda} = \frac{2}{\sqrt{\xi^2 + 1}}$，可直接和数值积分结果对比，二者误差在$10^{-10}$量级，符合计算要求。

内容的提问来源于stack exchange，提问作者Antonia Farina