点電荷による電気力線¶

点電荷による電界，電位¶

クーロンの法則は，$N$個の点電荷$Q_i$ [C]が位置ベクトル$\boldsymbol{r}_i$ [m]（$i=1,2, \cdots, N$）にあるとき，観測点$\boldsymbol{r}$ [m]の静電界$\boldsymbol{E}$ [N/C]は， \begin{gather} \boldsymbol{E} = \sum_{i=1}^N \boldsymbol{E}_i, \ \ \ \ \ \boldsymbol{E}_i = \frac{1}{4 \pi \epsilon_0} \frac{Q_i(\boldsymbol{r}-\boldsymbol{r}_i)}{|\boldsymbol{r}-\boldsymbol{r}_i|^3} \end{gather} ここで， \begin{gather} \frac{1}{4 \pi \epsilon_0} = 9 \times 10^9 \end{gather} ただし，$\epsilon_0= 8.842 \times 10^{-12}$ [F/m]は真空中の誘電率を示す．また，電位$V$は， \begin{gather} V = \frac{1}{4 \pi \epsilon_0} \sum_{i=1}^N \frac{Q_i}{R_i}, \ \ \ \ \ R_i = |\boldsymbol{r}-\boldsymbol{r}_i| \end{gather}

In [1]:
from json.encoder import INFINITY
import matplotlib.pyplot as plt  
import numpy as np  
from scipy.integrate import solve_ivp
import scienceplots
plt.style.use(['science', 'notebook'])
#plt.rcParams['font.family'] = 'Times New Roman' # font familyの設定
plt.rcParams['font.family'] = 'serif' # font familyの設定
plt.rcParams['mathtext.fontset'] = 'cm' # math fontの設定
In [2]:
def uvector(aa):# 多次元ベクトルの振幅，単位ベクトル
    norm = np.linalg.norm(aa)
    uv = aa / norm
    return norm, uv
In [3]:
# r: 観測点での位置ベクトルの配列
# rd: 点電荷のある位置ベクトルの配列
# q: 点電荷の電荷量の配列
def p1charge_3d(r, rd, q): # 1つの点電荷による電界，電位
    qk = q/4.0/np.pi
    rr = r-rd
    norm = np.linalg.norm(rr, axis=2, keepdims=True) # 観測点と点電荷との間の距離
    v = qk / norm # 電位の配列
    with np.errstate(divide="ignore", invalid="ignore"):
        e = qk*rr/norm**3 # 電界の配列
    return e, v, norm
In [4]:
def pncharge_3d(r, rdn, q): # 多数の点電荷による電界，電位
    rd = rdn[0].reshape(1, 1, -1)
    e, v, R = p1charge_3d(r, rd, q[0])
    nq = len(q)
    if nq>=2:
        for i in range(nq-1):
            rd = rdn[i+1].reshape(1, 1, -1)
            ei, vi, Ri = p1charge_3d(r, rd, q[i+1])
            e = e + ei
            v = v + vi
            R = np.append(R, Ri)
    return e, v, R
In [5]:
def pncharge_3d_fast(r, rdn, q):
    e = np.zeros_like(r)
    v = np.zeros(r.shape[:2])
    R_list = []
    for i in range(len(q)):
        ei, vi, Ri = p1charge_3d(r, rdn[i].reshape(1,1,-1), q[i])
        e += ei
        v += vi
        R_list.append(Ri)
    R = np.concatenate(R_list, axis=0)
    return e, v, R
In [6]:
# Fig.1.3(b) に対応
q = np.array([1.0, -1.0]) # 点電荷の電荷量
rd_p = np.array([[0.0, 1.0, 0.0]]) # 正の点電荷の位置ベクトル
rd_m = np.array([[0.0, -1.0, 0.0]]) # 負の点電荷の位置ベクトル
rd = np.concatenate([rd_p, rd_m])
rd
Out[6]:
array([[ 0.,  1.,  0.],
       [ 0., -1.,  0.]])
In [7]:
fig = plt.figure() # グラフ領域の作成
plt.axes().set_aspect('equal')
plt.xlim(-3.0, 3.0) # x軸範囲の設定
plt.ylim(-3.0, 3.0) # y軸範囲の設定
plt.xticks(np.linspace(-3.0, 3.0, 7))
plt.yticks(np.linspace(-3.0, 3.0, 7))
plt.xlabel("$x$ [m]") # 横軸のラベル設定
plt.ylabel("$y$ [m]") # 縦軸のラベル設定
plt.title("Electric Field Lines of a Dipole")
plt.plot(rd_p[:,0], rd_p[:,1],'o', label='Positive Charge')
plt.plot(rd_m[:,0], rd_m[:,1],'o', label='Negative Charge')
fig.tight_layout()
plt.legend(ncol=1, loc='best', fancybox=False, frameon = True)
fig.savefig('fig_1_3_b_none.pdf')
No description has been provided for this image

電気力線の方程式¶

電界$\boldsymbol{E} $を \begin{gather} \boldsymbol{E} = E \boldsymbol{a} = E_x \boldsymbol{a}_x + E_y \boldsymbol{a}_y + E_z \boldsymbol{a}_z \end{gather} より，電界に沿う単位ベクトル$\boldsymbol{a}$は， \begin{gather} \boldsymbol{a} = \frac{E_x}{E} \boldsymbol{a}_x + \frac{E_y}{E} \boldsymbol{a}_y + \frac{E_y}{E} \boldsymbol{a}_z \end{gather} ただし， $\boldsymbol{a}_x$，$\boldsymbol{a}_y$，$\boldsymbol{a}_z$は$x, y, z$方向の単位ベクトルを示す． この電界に沿うベクトル線要素$d\boldsymbol{s}$を， \begin{gather} d\boldsymbol{s} = \boldsymbol{a} ds = \boldsymbol{a}_x dx + \boldsymbol{a}_y dy + \boldsymbol{a}_z dz \end{gather} とおくと，単位ベクトル$\boldsymbol{a}$の成分比較より次のような常微分方程式が得られる． \begin{gather} \frac{dx}{ds} = \frac{E_x}{E} = \boldsymbol{a} \cdot \boldsymbol{a}_x\\ \frac{dy}{ds} = \frac{E_y}{E} = \boldsymbol{a} \cdot \boldsymbol{a}_y\\ \frac{dz}{ds} = \frac{E_z}{E} = \boldsymbol{a} \cdot \boldsymbol{a}_z \end{gather} 初期値$s=s_0$のとき， $x=x_0(s_0)$, $y=y_0(s_0)$， $z=z_0(s_0)$として， 積分形で表せば次のようになる． \begin{gather} x(s) = x_0 + \int_{s_0}^s \frac{E_x(s')}{E(s')} ds'\\ y(s) = y_0 + \int_{s_0}^s \frac{E_y(s')}{E(s')} ds'\\ z(s) = z_0 + \int_{s_0}^s \frac{E_z(s')}{E(s')} ds' \end{gather}

In [8]:
# q: 点電荷の電荷量のndarray
# rd: 点電荷の位置ベクトル(x,y)のndarray
def func_eline_3d(s, r, q, rd):
    e, v, R = pncharge_3d(r, rd, q)
    eamp, vs = uvector(e)
    R_min = np.min(R)
    if(R_min < 5.0e-4):
        vs = vs*0.0
    vs0 = vs[0,0,:]
    return vs0
In [9]:
s_ini, s_end, ns = 0.0, 4, 151
s = np.linspace(s_ini, s_end, ns)
nphi = 17
phi = np.linspace(0.0, 2.0*np.pi, nphi)
d_eps = 0.01
d_x, d_y = d_eps*np.cos(phi), d_eps*np.sin(phi)
In [10]:
def e_field_lines_2(q, rd, s, d_x, d_y, nphi):
    nq = len(q)
    ns = len(s)
    ode = np.empty([ns, 3, nphi*nq])
    ij = 0
    for j in range(nq):
        y0 = rd[j,:]
        for i in range(nphi):
            yeps = np.array([d_x[i], d_y[i], 0.0])
            y_ini = y0 + yeps
            if q[j] > 0:
                ode2 = solve_ivp(func_eline_3d, [s_ini, s_end], y_ini, dense_output=True, args=(q, rd,))
                y = ode2.sol(s)
            else:
                ode2 = solve_ivp(func_eline_3d, [s_ini, -s_end], y_ini, dense_output=True, args=(q, rd,))
                y = ode2.sol(-s)

            yt = np.transpose(y)
            plt.plot(yt[:,0], yt[:,1], linewidth=0.5, color="RebeccaPurple")

            # 矢印の追加
            arrow_idx = 7*ns // 20
            vec = yt[arrow_idx+1] - yt[arrow_idx-1]
            if q[j] < 0:
                vec = -vec
            dx, dy = vec[0], vec[1]
            plt.arrow(yt[arrow_idx, 0], yt[arrow_idx, 1], dx, dy,
                      shape='full', head_width=0.1, head_length=0.2,
                      fc="RebeccaPurple", ec="RebeccaPurple", linewidth=0)

            ode[:,:,ij] = yt
            ij += 1
    return ode
In [11]:
fig = plt.figure() # グラフ領域の作成
plt.axes().set_aspect('equal')
plt.xlim(-3.0, 3.0) # x軸範囲の設定
plt.ylim(-3.0, 3.0) # y軸範囲の設定
plt.xticks(np.linspace(-3.0, 3.0, 7))
plt.yticks(np.linspace(-3.0, 3.0, 7))
ode = e_field_lines_2(q, rd, s, d_x, d_y, nphi) # solve_ivp
plt.plot(rd_p[:,0], rd_p[:,1],'o', label='Positive Charge')
plt.plot(rd_m[:,0], rd_m[:,1],'o', label='Negative Charge')
fig.tight_layout()
plt.legend(ncol=1, loc='best', fancybox=False, frameon = True)
fig.savefig('fig_1_3_b.pdf')
plt.show()
No description has been provided for this image
In [12]:
q = np.array([1.0, 1.0, -1.0, -1.0])
rd_p = np.array([[2.0, 0.0, 0.0],[1.0, 0.0, 0.0]])
rd_m = np.array([[-1.0, 0.0, 0.0],[-2.0, 0.0, 0.0]])
rd = np.concatenate([rd_p, rd_m])
rd
Out[12]:
array([[ 2.,  0.,  0.],
       [ 1.,  0.,  0.],
       [-1.,  0.,  0.],
       [-2.,  0.,  0.]])
In [13]:
fig = plt.figure() # グラフ領域の作成
plt.axes().set_aspect('equal')
plt.xlim(-3.0, 3.0) # x軸範囲の設定
plt.ylim(-3.0, 3.0) # y軸範囲の設定
plt.xticks(np.linspace(-3.0, 3.0, 7))
plt.yticks(np.linspace(-3.0, 3.0, 7))
ode = e_field_lines_2(q, rd, s, d_x, d_y, nphi) # solve_ivp
plt.plot(rd_p[:,0], rd_p[:,1],'o', label='Positive Charge')
plt.plot(rd_m[:,0], rd_m[:,1],'o', label='Negative Charge')
fig.tight_layout()
plt.legend(ncol=1, loc='best', fancybox=False, frameon = True)
fig.savefig('fig_1_4_b.pdf')
plt.show()
No description has been provided for this image