Matplotlib で電場ベクトルの方向場を描く

正負の点電荷がつくる電場の向きを表す方向場を Python の Matplotlib で描く。

正負の点電荷がつくる電場

$\boldsymbol{r}_1$ にある正電荷 $q\ (>0)$ と，$\boldsymbol{r}_2$ にある負電荷 $-q$ がつくる電場ベクトル $\boldsymbol{E}$ は

\begin{eqnarray}
\boldsymbol{E}
&=& \frac{q}{4\pi \varepsilon_0} \frac{\boldsymbol{r} – \boldsymbol{r}_1}{|\boldsymbol{r} – \boldsymbol{r}_1|^3}
– \frac{q}{4\pi \varepsilon_0} \frac{\boldsymbol{r} – \boldsymbol{r}_2}{|\boldsymbol{r} – \boldsymbol{r}_2|^3}
\end{eqnarray}

$z=0$ の $xy$ 平面で2次元ベクトル場を描く。適宜無次元化を行い，
$\displaystyle \frac{q}{4\pi \varepsilon_0}\rightarrow 1$ とする。また，

$$\boldsymbol{r} = (x, y, 0), \quad \boldsymbol{r}_1 = (0, 1, 0),
\quad \boldsymbol{r}_2 = (0, -1, 0)$$

とすると，

$$E_x = \frac{x}{\left(x^2 + (y-1)^2\right)^{3/2}}-\frac{x}{\left(x^2 + (y+1)^2\right)^{3/2}}$$$$E_y = \frac{y-1}{\left(x^2 + (y-1)^2\right)^{3/2}}-
\frac{y+1}{\left(x^2 + (y+1)^2\right)^{3/2}}$$

$\boldsymbol{E}$ の向きを表す方向場（単位ベクトル場）

一般にベクトルは大きさと向きを表す。電場 $\boldsymbol{E}$ の大きさは電荷からの距離の2乗に反比例するので，ベクトルの長さが大きく変化する。見やすくするために，以下のように $\boldsymbol{E}$ の向きを表す方向場，すなわち大きさが $1$ の単位ベクトル場を定義して，ベクトルの向きのみを正しく描画する。

$$\hat{\boldsymbol{E}} \equiv \frac{\boldsymbol{E}}{\sqrt{\boldsymbol{E}\cdot\boldsymbol{E}}}$$$$\hat{E}_x = \frac{E_x}{\sqrt{E_x^2 + E_y^2}}$$$$\hat{E}_y = \frac{E_y}{\sqrt{E_x^2 + E_y^2}}$$

matplotlib.pyplot.quiver() で電場の方向場を描く

グラフのサイズ等のデフォルト設定

In [1]:
# NumPy 
import numpy as np
# Matplotlib 
import matplotlib.pyplot as plt

# グラフを SVG で Notebook にインライン表示させる設定
%config InlineBackend.figure_formats = ['svg']
# グラフのサイズ
plt.rcParams["figure.figsize"] = [6.4, 6.4]
# フォント設定
plt.rcParams['font.family'] = 'Noto Sans CJK JP'
plt.rcParams['font.size'] = 14
plt.rcParams['xtick.labelsize'] = 12
plt.rcParams['ytick.labelsize'] = 12

電場ベクトルの定義

電場ベクトルを x, y の2変数関数として定義。

In [2]:
def Ex(x, y):
    return x/np.sqrt(x**2 + y**2)**3
def Ey(x, y):
    return y/np.sqrt(x**2 + y**2)**3

def E2x(x, y):
    return Ex(x, y-1) - Ex(x, y+1)
def E2y(x, y):
    return Ey(x, y-1) - Ey(x, y+1)

def E(x, y):
    return np.sqrt(E2x(x, y)**2 + E2y(x, y)**2)

def hE2x(x, y):
    return scaling * E2x(x, y)/E(x, y)
def hE2y(x, y):
    return scaling * E2y(x, y)/E(x, y)

ベクトル場を描く

点電荷が $y$ 軸上にあるので，リストを作成する際，ベクトルの始点が $y$ 軸に来ないようにしておく。（点電荷上では電場の分母がゼロになってしまうから。）
以下の例では，$x$ 座標は

$$-5.5, -4.5, \cdots, -0.5, 0.5, \cdots, 5.5$$

となり，$y$ 軸上 ($x = 0$) に始点がこない。

In [3]:
x = y = np.linspace(-5.5, 5.5, 12)
x, y = np.meshgrid(x, y)

plt.quiver() の書式は…

plt.quiver(始点のx座標, 始点のy座標, ベクトルのx成分, ベクトルのy成分)
In [4]:
# ベクトルの長さ調整用
scaling = 0.5

# 表示範囲
plt.xlim(-6.5, 6.5)
plt.ylim(-6.5, 6.5)
# グリッドの表示
plt.grid()
# x の主目盛を 1 刻みに
plt.xticks(np.arange(-5, 6, 1))
# y の主目盛を 1 刻みに
plt.yticks(np.arange(-5, 6, 1))

# ベクトル場の表示
plt.quiver(x, y, hE2x(x, y), hE2y(x, y));

点電荷の位置にラベル

正電荷，負電荷の位置にラベルをつけてみる。

In [5]:
# ベクトルの長さ調整用
scaling = 0.5

# 表示範囲
plt.xlim(-6.5, 6.5)
plt.ylim(-6.5, 6.5)
# グリッドの表示
plt.grid()
# x の主目盛を 1 刻みに
plt.xticks(np.arange(-6, 7, 1))
# y の主目盛を 1 刻みに
plt.yticks(np.arange(-6, 7, 1))

# 
plt.quiver(x, y, hE2x(x, y), hE2y(x, y))

# 正電荷の位置に青丸と +
plt.scatter(0, 1, marker='o', 
            # size, color
            s = 300, c = 'white',
            edgecolors = 'blue', linewidths = 2)
plt.text(0, 1, "＋", color = "blue",
         horizontalalignment = "center", 
         verticalalignment = "center")
# 負電荷の位置に赤丸と -
plt.scatter(0, -1, marker='o', 
            # size, color
            s = 300, c = 'white',
            edgecolors = 'red', linewidths = 2)
plt.text(0, -1, "ー", color = "red",
         horizontalalignment = "center", 
         verticalalignment = "center");

Colormap (cmap)

plt.quiver() の書式で5番目の引数を設定すると，その大きさに応じたカラーマップが使える。

plt.quiver(x, y, vx, vy, C, cmap = ...)

以下の例では電場の大きさ E(x, y) に応じて矢印の色が変化する。
カラーマップ cmap = "Blues" を指定しているので，電場の強い（E(x, y) の値が大きい）ところは濃い青で，弱い（値が小さい）ところほど薄い青になる。

In [6]:
# ベクトルの長さ調整用
scaling = 0.5

# 表示範囲
plt.xlim(-6.5, 6.5)
plt.ylim(-6.5, 6.5)
# グリッドの表示
plt.grid()
# x の主目盛を 1 刻みに
plt.xticks(np.arange(-6, 7, 1))
# y の主目盛を 1 刻みに
plt.yticks(np.arange(-6, 7, 1))

# 電場の大きさに応じたカラーマップの設定
plt.quiver(x, y, hE2x(x, y), hE2y(x, y), E(x, y), cmap = "Blues");

cmap の色味，ベクトルの大きさ・太さの設定

電荷の位置から離れると，電場の大きさは急激に小さくなるので，このままでは薄くてほとんど見えないベクトルばかりになる。

これだと少し寂しいのでカラーマップの範囲を設定する clim を調整して，離れた位置のベクトルもそこそこの青色になるように調整する。

上記の例では clim = (0, 2) くらいであろうから，ゼロに近い大きさでも色がつくように手探りで数値をいれてみます。

また，ベクトルそのものの大きさ・太さも少し変更してみる。書式はまず units='xy' などと設定後，ベクトルの長さを伸ばしたいなら scale により小さい値を，太くしたいなら width により大きい値をいれます。これも手探りでいれてみました。

In [7]:
# ベクトルの長さ調整用
scaling = 0.5

# 表示範囲
plt.xlim(-6.5, 6.5)
plt.ylim(-6.5, 6.5)
# グリッドの表示
plt.grid()
# x の主目盛を 1 刻みに
plt.xticks(np.arange(-6, 7, 1))
# y の主目盛を 1 刻みに
plt.yticks(np.arange(-6, 7, 1))

# 
plt.quiver(x, y, hE2x(x, y), hE2y(x, y), E(x, y), cmap = "Blues", 
           # 数値は手探りで
           clim=(-0.1, 0.8),
           # units 設定後，長さは scale, 太さは width で
           units='xy', scale = 0.5, width=0.15);

一応の完成作

以上の設定をまとめてみたのが以下の例。少し短めのベクトルをたくさん描いてみた。

In [8]:
x = y = np.linspace(-5.25, 5.25, 22)
x, y = np.meshgrid(x, y)

# ベクトルの長さ調整用
scaling = 0.23

# 枠と縦軸目盛を非表示に
fig, ax = plt.subplots()
ax.spines['right'].set_visible(False)
ax.spines['left'].set_visible(False)
ax.spines['top'].set_visible(False)
ax.spines['bottom'].set_visible(False)
ax.axes.xaxis.set_visible(False)
ax.axes.yaxis.set_visible(False)

# 表示範囲
plt.xlim(-6, 6)
plt.ylim(-6, 6)

plt.quiver(x, y, hE2x(x, y), hE2y(x, y), E(x, y), cmap = "Blues", 
           # 数値は手探りで
           clim=(-0.1, 0.8), 
           # units 設定後，長さは scale, 太さは width で
           units='xy', scale = 0.4, width=0.08)

# 正電荷の位置に青丸と +
plt.scatter(0, 1, marker='o', 
            # size, color
            s = 300, c = 'white',
            edgecolors = 'blue', linewidths = 2)
plt.text(0, 1, "＋", color = "blue",
         horizontalalignment = "center", 
         verticalalignment = "center")
# 負電荷の位置に赤丸と -
plt.scatter(0, -1, marker='o', 
            # size, color
            s = 300, c = 'white',
            edgecolors = 'red', linewidths = 2)
plt.text(0, -1, "ー", color = "red",
         horizontalalignment = "center", 
         verticalalignment = "center");