電磁気学の授業の準備「gnuplot でベクトル場（電場・磁場）を描く」の Maxima 版。

点電荷がつくる電場 $\boldsymbol{E}$

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

$z=0$ の $xy$ 平面で2次元ベクトル場を描く。

$$\boldsymbol{r} = (x, y, 0), \quad \boldsymbol{r}_1 = (x_1, y_1, 0)$$$$E_x = s \frac{x-x_1}{\left((x-x_1)^2 + (y-y_1)^2\right)^{3/2}}$$$$E_y = s \frac{y-y_1}{\left((x-x_1)^2 + (y-y_1)^2\right)^{3/2}}$$

$s$ は $\mbox{scaling}$ 係数。実際の描画の際にベクトルの長さを微調整する係数。

参考：「点電荷の電荷密度と電場」

正電荷がつくる電場

In [1]:

/* 電場 E */
Ex(x, y):= scaling * x/sqrt(x**2 + y**2)**3$
Ey(x, y):= scaling * y/sqrt(x**2 + y**2)**3$

/* 2次元ベクトルは vector([始点のx, 始点のy], [x成分, y成分]) */
vecE(x, y):= vector([x, y], [Ex(x-x1, y-y1), Ey(x-x1, y-y1)])$

/* (x1, y1) は正電荷の位置 */
[x1, y1]: [0, 1]$

/* -3 <= x <= 3, -3 <= y <= 3 で 1/2 刻みの格子点を始点に */
/* makelist では1次元のリストしか作れないようなので，少し工夫 */
coordx: setify(makelist(-3+i/2, i, 0, 12))$
coordxy: listify(cartesian_product(coordx,coordx))$

/* 正電荷の位置では電場の分母がゼロになるので除く */
coordxy: float(delete([x1, y1], coordxy))$

scaling: 0.2$
/* ベクトル場のデータ */
vecfE: makelist(vecE(k[1], k[2]), k, coordxy)$

/* draw the vector field */
draw2d(
  /* title のフォントサイズの変更例。gnuplot の流儀。 */
  title = "{/=16 正電荷による電場}", 
  xtics=false, ytics=false,

  /* 表示範囲 */
  xrange = [-3, 3], yrange = [-3, 3],
  /* 縦横比 */
  proportional_axes=xy,

  /* 正電荷の位置 */
  color = black,
  point_size = 1.5, 
  point_type = 6, 
  points([[x1, y1]]),
  label(["＋", x1, y1]), 

  /* ベクトル場の表示 */
  color = blue,
  head_length = 0.2,
  head_angle  = 20, 
  vecfE
)$

In [2]:

/* 弘大 JupyterHub では 
   set_draw_defaults(file_name="~/.maxplot",terminal='svg)$
   されているので。*/
system("cp ~/.maxplot.svg ./maxvec-fig01.svg")$

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

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

$$\hat{\boldsymbol{E}} \equiv \frac{\boldsymbol{E}}{\sqrt{\boldsymbol{E}\cdot\boldsymbol{E}}} $$

正電荷がつくる電場の向き

In [3]:

/* 電場 E */
Ex(x, y):= x/sqrt(x**2 + y**2)**3$
Ey(x, y):= y/sqrt(x**2 + y**2)**3$

/* 規格化された電場 E */
E(x, y):= sqrt(Ex(x, y)**2 + Ey(x, y)**2)$
hEx(x, y):= scaling * Ex(x, y)/E(x, y)$
hEy(x, y):= scaling * Ey(x, y)/E(x, y)$

vechE(x, y):= vector([x, y], [hEx(x-x1, y-y1), hEy(x-x1, y-y1)])$

scaling: 0.3$
/* ベクトル場のデータ */
vecfhE: makelist(vechE(k[1], k[2]), k, coordxy)$

draw2d(
  /* title のフォントサイズの変更例。gnuplot の流儀。 */
  title = "{/=16 正電荷による電場の向き}", 
  xtics=false, ytics=false,

  /* 表示範囲 */
  xrange = [-3, 3], yrange = [-3, 3],
  /* 縦横比 */
  proportional_axes=xy,

  /* 正電荷の位置 */
  color = black,
  point_size = 1.5, 
  point_type = 6, 
  points([[x1, y1]]),
  label(["＋", x1, y1]), 

  /* ベクトル場の表示 */
  color = blue,
  head_length = 0.1,
  head_angle  = 20, 
  vecfhE
)$

In [4]:

/* 弘大 JupyterHub では 
   set_draw_defaults(file_name="~/.maxplot",terminal='svg)$
   されているので。*/
system("cp ~/.maxplot.svg ./maxvec-fig02.svg")$

負電荷がつくる電場の向き

In [5]:

/* 電場 E */
Ex(x, y):= -x/sqrt(x**2 + y**2)**3$
Ey(x, y):= -y/sqrt(x**2 + y**2)**3$

/* 規格化された電場 E */
E(x, y):= sqrt(Ex(x, y)**2 + Ey(x, y)**2)$
hEx(x, y):= scaling * Ex(x, y)/E(x, y)$
hEy(x, y):= scaling * Ey(x, y)/E(x, y)$

/* 2次元ベクトルは vector([始点のx, 始点のy], [x成分, y成分]) */
vechE(x, y):= vector([x, y], [hEx(x-x2, y-y2), hEy(x-x2, y-y2)])$

/* (x2, y2) は負電荷の位置 */
[x2, y2]: [0, -1]$

/* -3 <= x <= 3, -3 <= y <= 3 で 1/2 刻みの格子点を始点に */
/* makelist では1次元のリストしか作れないようなので，少し工夫 */
/* 別解としてこういうやり方もあるということで。 */
/* makelist() で空のリストをつくり，append() で要素を追加 */
coordxy: makelist()$
for i: 0 thru 12 do(
  for j: 0 thru 12 do(
    coordxy: append(coordxy, [[-3+i/2,-3+j/2]])
  )
)$

/* 負電荷の位置では電場の分母がゼロになるので除く */
coordxy: float(delete([x2, y2], coordxy))$

scaling: 0.3$
/* ベクトル場のデータ */
vecfhE: makelist(vechE(k[1], k[2]), k, coordxy)$

draw2d(
  /* title のフォントサイズの変更例。gnuplot の流儀。 */
  title = "{/=16 負電荷による電場の向き}", 
  xtics=false, ytics=false,

  /* 表示範囲 */
  xrange = [-3, 3], yrange = [-3, 3],
  /* 縦横比 */
  proportional_axes=xy,

  /* 負電荷の位置 */
  color = red, 
  point_size = 1.5, 
  point_type = 6, 
  points([[x2, y2]]),
  label(["ー", x2, y2]), 

  /* ベクトル場の表示 */
  color = blue,
  head_length = 0.1,
  head_angle  = 20, 
  vecfhE
)$

In [6]:

/* 弘大 JupyterHub では 
   set_draw_defaults(file_name="~/.maxplot",terminal='svg)$
   されているので。*/
system("cp ~/.maxplot.svg ./maxvec-fig03.svg")$

正電荷と負電荷がつくる電場の向き

参考：「2つの点電荷による電場」

In [7]:

/* 電場 E */
Ex(x, y):= x/sqrt(x**2 + y**2)**3$
Ey(x, y):= y/sqrt(x**2 + y**2)**3$

/* 正電荷と負電荷がつくる電場 */
E2x(x, y):= Ex(x-x1, y-y1) - Ex(x-x2, y-y2)$
E2y(x, y):= Ey(x-x1, y-y1) - Ey(x-x2, y-y2)$

/* 規格化された電場 hat E */
E(x, y):= sqrt(E2x(x, y)**2 + E2y(x, y)**2)$
hE2x(x, y):= scaling * E2x(x, y)/E(x, y)$
hE2y(x, y):= scaling * E2y(x, y)/E(x, y)$

/* 2次元ベクトルは vector([始点のx, 始点のy], [x成分, y成分]) */
vechE(x, y):= vector([x, y], [hE2x(x, y), hE2y(x, y)])$

/* (x1, y1) は正電荷の位置 */
[x1, y1]: [0, 1]$
/* (x2, y2) は負電荷の位置 */
[x2, y2]: [0, -1]$

/* -3 <= x <= 3, -3 <= y <= 3 で 1/3 刻みの格子点を始点に */
coordxy: makelist()$
for i: 0 thru 18 do(
  for j: 0 thru 18 do(
    coordxy: append(coordxy, [[-3+i/3,-3+j/3]])
  )
)$

/* 電荷の位置では電場の分母がゼロになるので除く */
coordxy: float(delete([x2, y2], delete([x1, y1], coordxy)))$

scaling: 0.2$
/* ベクトル場のデータ */
vecfhE: makelist(vechE(k[1], k[2]), k, coordxy)$

draw2d(
  /* title のフォントサイズの変更例。gnuplot の流儀。 */
  title = "{/=16 正電荷と負電荷がつくる電場の向き}", 
  xtics=false, ytics=false,

  /* 表示範囲 */
  xrange = [-3, 3], yrange = [-3, 3],
  /* 縦横比 */
  proportional_axes=xy,
  
  /* 正電荷の位置 */
  color = black,
  point_size = 1.5, 
  point_type = 6, 
  points([[x1, y1]]),
  label(["＋", x1, y1]), 
  
  /* 負電荷の位置 */
  color = red, 
  point_size = 1.5, 
  point_type = 6, 
  points([[x2, y2]]),
  label(["ー", x2, y2]), 
  
  /* ベクトル場の表示 */
  color = blue,
  head_length = 0.1,
  head_angle  = 20, 
  vecfhE
)$

In [8]:

/* 弘大 JupyterHub では 
   set_draw_defaults(file_name="~/.maxplot",terminal='svg)$
   されているので。*/
system("cp ~/.maxplot.svg ./maxvec-fig04.svg")$

電気双極子がつくる電場

参考：「電気双極子による電場」

$$\boldsymbol{E} = 3 \frac{\boldsymbol{r}\cdot\boldsymbol{p}}{r^5} \boldsymbol{r}
– \frac{\boldsymbol{p}}{r^3}$$

$\boldsymbol{p} = (0, 1, 0), \quad z = 0$ として

\begin{eqnarray}
E_x &=& \frac{3}{\left( x^2 + y^2 \right)^{5/2}} x y \\
E_y &=& \frac{3}{\left( x^2 + y^2 \right)^{5/2}} y^2 – \frac{1}{\left( x^2 + y^2 \right)^{3/2}} \\
\end{eqnarray}

電気双極子がつくる電場の向き

In [9]:

/* 電場 E */
Ex(x, y):= 3* x*y/sqrt(x**2 + y**2)**5$
Ey(x, y):= 3* y*y/sqrt(x**2 + y**2)**5 - 1/sqrt(x**2 + y**2)**3$

/* 規格化された電場 E */
E(x, y):= sqrt(Ex(x, y)**2 + Ey(x, y)**2)$
hEx(x, y):= scaling * Ex(x, y)/E(x, y)$
hEy(x, y):= scaling * Ey(x, y)/E(x, y)$

/* 2次元ベクトルは vector([始点のx, 始点のy], [x成分, y成分]) */
vechE(x, y):= vector([x, y], [hEx(x, y), hEy(x, y)])$

/* (x0, y0) は電気双極子の位置 */
[x0, y0]: [0, 0]$

/* -3 <= x <= 3, -3 <= y <= 3 で 1/3 刻みの格子点を始点に */
coordxy: makelist()$
for i: 0 thru 18 do(
  for j: 0 thru 18 do(
    coordxy: append(coordxy, [[-3+i/3,-3+j/3]])
  )
)$

/* 電気双極子の位置では電場の分母がゼロになるので除く */
coordxy: float(delete([x0, y0], coordxy))$

scaling: 0.25$
/* ベクトル場のデータ */
vecfhE: makelist(vechE(k[1], k[2]), k, coordxy)$

draw2d(
  xaxis = true, yaxis = true,
  /* title のフォントサイズの変更例。gnuplot の流儀。 */
  title = "{/=16 電気双極子による電場の向き}", 
  xtics=false, ytics=false,
  
  /* 表示範囲 */
  xrange = [-3, 3], yrange = [-3, 3],
  /* 縦横比 */
  proportional_axes=xy,
  
  /* 電気双極子の位置 */
  color = red, 
  line_width = 5,
  head_length = 0.1,
  head_angle  = 30, 
  vector([0,-0.15],[0,0.3]), 
  
  /* ベクトル場の表示 */
  color = blue,
  line_width = 1,
  head_length = 0.1,
  head_angle  = 20, 
  vecfhE
)$

In [10]:

/* 弘大 JupyterHub では 
   set_draw_defaults(file_name="~/.maxplot",terminal='svg)$
   されているので。*/
system("cp ~/.maxplot.svg ./maxvec-fig05.svg")$

直線電流がつくる磁場 $\boldsymbol{B}$

参考：「直線電流による磁場」

$$\boldsymbol{B} = \frac{\boldsymbol{I}\times\boldsymbol{\rho}}{\rho^2}$$$$ \boldsymbol{I} = (0, 0, I), \quad \boldsymbol{\rho} = (x, y, 0)$$\begin{eqnarray}
B_x &=& – \frac{ I y}{x^2 + y^2} \\
B_y &=& \frac{ I x}{x^2 + y^2}
\end{eqnarray}

直線電流がつくる磁場

In [11]:

/* 磁場 B */
Bx(x, y):= - scaling * y/(x**2 + y**2)$
By(x, y):=   scaling * x/(x**2 + y**2)$

/* 2次元ベクトルは vector([始点のx, 始点のy], [x成分, y成分]) */
vecB(x, y):= vector([x, y], [Bx(x, y), By(x, y)])$

/* -3 <= x <= 3, -3 <= y <= 3 で 1/2 刻みの格子点を始点に */
coordxy: makelist()$
for i: 0 thru 12 do(
  for j: 0 thru 12 do(
    coordxy: append(coordxy, [[-3+i/2,-3+j/2]])
  )
)$

/* 原点では磁場の分母がゼロになるので除く */
coordxy: float(delete([0, 0], coordxy))$

scaling: 0.25$
/* ベクトル場のデータ */
vecfB: makelist(vecB(k[1], k[2]), k, coordxy)$

draw2d(
  xaxis = true, yaxis = true,
  /* title のフォントサイズの変更例。gnuplot の流儀。 */
  title = "{/=16 直線電流による磁場}", 
  xtics=false, ytics=false,
  
  /* 表示範囲 */
  xrange = [-3, 3], yrange = [-3, 3],
  /* 縦横比 */
  proportional_axes=xy,
  
  /* 電線の位置 */
  color = red,
  point_size = 1.5, 
  point_type = 6, 
  points([[0, 0]]),
  label(["・", 0, 0]), 
  
  /* ベクトル場の表示 */
  color = green,
  line_width = 1.5,
  head_length = 0.1,
  head_angle  = 20, 
  vecfB
)$

In [12]:

/* 弘大 JupyterHub では 
   set_draw_defaults(file_name="~/.maxplot",terminal='svg)$
   されているので。*/
system("cp ~/.maxplot.svg ./maxvec-fig06.svg")$

直線電流がつくる磁場の向き

磁場 $\boldsymbol{B}$ の大きさは電線からの距離に反比例するので，上図のようにベクトルの長さが大きく変化する。これはこれで正しいのだが，見やすくするために，以下のように $\boldsymbol{B}$ の向きを表す単位ベクトル場（大きさが $1$）を定義して，ベクトルの向きのみを正しく描画する。

$$\hat{\boldsymbol{B}} \equiv \frac{\boldsymbol{B}}{\sqrt{\boldsymbol{B}\cdot\boldsymbol{B}}}$$

In [13]:

/* 磁場 B */
Bx(x, y):= - y/(x**2 + y**2)$
By(x, y):=   x/(x**2 + y**2)$

/* 規格化された hat B */
B(x, y):= sqrt(Bx(x, y)**2 + By(x, y)**2)$
hBx(x, y):= scaling * Bx(x, y)/B(x, y)$
hBy(x, y):= scaling * By(x, y)/B(x, y)$

/* 2次元ベクトルは vector([始点のx, 始点のy], [x成分, y成分]) */
vechB(x, y):= vector([x, y], [hBx(x, y), hBy(x, y)])$

/* -3 <= x <= 3, -3 <= y <= 3 で 1/3 刻みの格子点を始点に */
coordxy: makelist()$
for i: 0 thru 18 do(
  for j: 0 thru 18 do(
    coordxy: append(coordxy, [[-3+i/3,-3+j/3]])
  )
)$

/* 原点では磁場の分母がゼロになるので除く */
coordxy: float(delete([0, 0], coordxy))$

scaling: 0.25$
/* ベクトル場のデータ */
vecfhB: makelist(vechB(k[1], k[2]), k, coordxy)$

draw2d(
  xaxis = true, yaxis = true,
  /* title のフォントサイズの変更例。gnuplot の流儀。 */
  title = "{/=16 直線電流による磁場の向き}", 
  xtics=false, ytics=false,
  
  /* 表示範囲 */
  xrange = [-3, 3], yrange = [-3, 3],
  /* 縦横比 */
  proportional_axes=xy,
  
  /* 電線の位置 */
  color = red,
  point_size = 1.5, 
  point_type = 6, 
  points([[0, 0]]),
  label(["・", 0, 0]), 
  
  /* ベクトル場の表示 */
  color = green,
  line_width = 1.5,
  head_length = 0.1,
  head_angle  = 20, 
  vecfhB
)$

In [14]:

/* 弘大 JupyterHub では 
   set_draw_defaults(file_name="~/.maxplot",terminal='svg)$
   されているので。*/
system("cp ~/.maxplot.svg ./maxvec-fig07.svg")$

Maxima でベクトル場（電場・磁場）を描く