gnuplot で座標をファイルから読み込んで電気力線・磁力線を描く

gnuplot でベクトル場(電場・磁場)を描く」や「palette でベクトルの色を変えて gnuplot で電場・磁場を描く」で描いてみたが,もう少し。

gnuplot だけでベクトルの始点の座標ファイルを作成し,それを読み込んで電気力線・磁力線を描く例。

点電荷による電場 $\boldsymbol{E}$

正の点電荷による電場

In [20]:
%gnuplot inline svg size 480,480 fixed enhanced font 'Noto Sans CJK JP,16'
In [27]:
reset
unset xtics
unset ytics
unset border

# 正電荷の位置
x1 = 0
y1 = 1

# 電気力線用座標ファイル
set print "xy-line.txt"
do for [j=0: 15] {
  do for [i=1:60] {
    r = 0.2 * i
    theta = pi/8 * j
    x = x1 + r * cos(theta)
    y = y1 + r * sin(theta)
    print sprintf("%8.4f  %8.4f", x, y)
  }
  # with lines するので適宜空行をいれる
  print ""
}
set print

# ベクトルの始点用座標ファイル
# 矢印の混雑防止のため適宜間引く
set print "xy.txt"
do for [j=0: 15] {
  do for [i=2:8:3] {
    r = 0.5 * i
    theta = pi/8 * j
    x = x1 + r * cos(theta)
    y = y1 + r * sin(theta)
    print sprintf("%8.4f  %8.4f", x, y)
  }
}
set print

# 電場 E
Ex(x, y) = (x-x1)/(sqrt((x-x1)**2 + (y-y1)**2)**3)
Ey(x, y) = (y-y1)/(sqrt((x-x1)**2 + (y-y1)**2)**3)
# 電場 E の大きさ
E(x, y) = sqrt(Ex(x, y)**2 + Ey(x, y)**2)
# 規格化された hat E
hEx(x, y) = Ex(x, y)/E(x, y)
hEy(x, y) = Ey(x, y)/E(x, y)

set xrange [-3.5:3.5]
set yrange [-3.5:3.5]
set size ratio 1

# 矢印の長さ調整用
scaling = 0.1

# 正電荷の位置
set label 1 point pt 6 ps 1.6 lc "black" at x1, y1
set label 2 center at first x1, y1 "+" tc "black" font ",22"

set title "正電荷による電場"
unset colorbox
set palette cubehelix start 0 cycles 0 saturation 3 negative

plot \
  "xy-line.txt" \
    u 1:2:(log(E($1,$2))) w l lw 3 lc palette notitle,\
  "xy.txt" \
    u 1:2:(scaling*hEx($1,$2)):(scaling*hEy($1,$2)):(log(E($1,$2))) \
    w vec notitle lw 8 lc palette filled head

In [3]:
set output "./vec-fig01d.svg"
replot
set output

負の点電荷による電場

In [28]:
reset
unset xtics
unset ytics
unset border
set xrange [-3.5:3.5]
set yrange [-3.5:3.5]
set size ratio 1

# 極座標で x, y のデータファイル作成
# 負電荷の位置
x2 = 0
y2 = -1

# 電気力線用座標ファイル
set print "xy-line.txt"
do for [j=0:15] {
  do for [i=1:60] {
    r = 0.2 * i
    theta = pi/8 * j
    x = x2 + r * cos(theta)
    y = y2 + r * sin(theta)
    print sprintf("%8.4f  %8.4f", x, y)
  }
  # with lines するので適宜空行をいれる
  print ""
}
set print

# ベクトルの始点用座標ファイル
# 矢印の混雑防止のため適宜間引く
set print "xy.txt"
do for [j=0:15] {
  do for [i=3:9:3] {
    r = 0.5 * i
    theta = pi/8 * j
    x = x2 + r * cos(theta)
    y = y2 + r * sin(theta)
    print sprintf("%8.4f  %8.4f", x, y)
  }
}
set print

# 電場 E
Ex(x, y) = -(x-x2)/(sqrt((x-x2)**2 + (y-y2)**2)**3)
Ey(x, y) = -(y-y2)/(sqrt((x-x2)**2 + (y-y2)**2)**3)
# 電場 E の大きさ
E(x, y) = sqrt(Ex(x, y)**2 + Ey(x, y)**2)
# 規格化された hat E
hEx(x, y) = Ex(x, y)/E(x, y)
hEy(x, y) = Ey(x, y)/E(x, y)

set xrange [-3.5:3.5]
set yrange [-3.5:3.5]
set size ratio 1

# 矢印の長さ調整用
scaling = 0.1

# 負電荷の位置
set label 3 point pt 6 ps 1.6 lc "black" at x2, y2
set label 4 center at first x2, y2 "ー" tc "black" font ",20"

set title "負電荷による電場"
unset colorbox
set palette cubehelix start 0 cycles 0 saturation 3 negative

plot "xy-line.txt" u 1:2:(log(E($1,$2))) w l lw 3 lc palette notitle,\
  "xy.txt" \
  u 1:2:(scaling*hEx($1,$2)):(scaling*hEy($1,$2)):(log(E($1,$2))) w vec \
  notitle lw 8 lc palette filled head

In [5]:
set output "./vec-fig02d.svg"
replot
set output

一様な線電荷による電場

$z$ 軸上の一様な線電荷密度 $\lambda (\mbox{C}/\mbox{m})$ による電場は,$\boldsymbol{r} \equiv (x, y, 0), \ r^2 \equiv \boldsymbol{r} \cdot\boldsymbol{r} $ としてまとめると

$$\boldsymbol{E} = \frac{\lambda}{2\pi\varepsilon_0} \frac{\boldsymbol{r}}{r^2}$$

定数部分を適宜規格化して

\begin{eqnarray}
E_x &\Rightarrow& \frac{x}{x^2 + y^2} \\
E_x &\Rightarrow& \frac{y}{x^2 + y^2} \\
\end{eqnarray}

In [40]:
# 電気力線
set print "xyz-line.txt"
do for [z=0:2]{
  do for [i=0:11]{
    theta = 2*pi/12 * i
    do for [j=1:20] {
      r = 0.2*j
      x = r * cos(theta)
      y = r * sin(theta)
      print sprintf(" %8.4f %8.4f %8.4f", x, y, z)
    }
    print ""
    print ""
  }
    print ""
    print ""
}
set print

# 電場ベクトルの始点
set print "xyz.txt"
do for [z=0:2]{
  do for [i=0:11:1]{
    theta = 2*pi/12 * i
    do for [j=3:6:3] {
      r = 0.5*j
      x = r * cos(theta)
      y = r * sin(theta)
      print sprintf(" %8.4f %8.4f %8.4f", x, y, z)
    }
    print ""
  }
    print ""
}
set print
In [41]:
%gnuplot inline svg size 640,640 fixed enhanced font 'Noto Sans CJK JP,16'
In [45]:
reset

Ex(x,y) = x/(x**2 + y**2)
Ey(x,y) = y/(x**2 + y**2)
E(x, y) = sqrt(Ex(x,y)**2 + Ey(x,y)**2)
ex(x,y) = scaling*Ex(x,y)/E(x, y)
ey(x,y) = scaling*Ey(x,y)/E(x, y)

unset xtics
unset ytics
unset ztics
# unset border
set xyplane 0
set zeroaxis
set xrange [-4:4]
set yrange [-4:4]
unset colorbox
set palette cubehelix start 0 cycles 0 saturation 3 negative

# 線電荷
set arrow 1 from 0,0,-0.2 to 0,0,2.2 lw 4 lc "black" lt 3  nohead
set key inside
set title "一様な線電荷による電場"
scaling = 0.1
set view 65,15,,
splot "xyz-line.txt" u 1:2:3:(log(E($1,$2))) w l lw 3 lc palette notitle, \
      "xyz.txt" u 1:2:3:(ex($1,$2)):(ey($1,$2)):(0):(log(E($1,$2))) \
        w vec filled head lw 6 lc palette title "電場 E"

In [9]:
set output "./vec-fig-sen.svg"
replot
set output

一様な面電荷による電場

In [10]:
%gnuplot inline svg size 640,480 fixed enhanced font 'Noto Sans CJK JP,16'
In [11]:
reset
set parametric
set isosamples 30
set hidden3d
set xrange [-3:3]
set xrange [-3:3]
set zrange [-3:3]
unset xtics
unset ytics
unset ztics
unset border
set zeroaxis
set xyplane 0
set isosamples 60
set view ,20,,

# 一様な電場ベクトルを arrow で描く
do for [i=1:5] {
  zi = -3 + i
  # x > 0 側
  set arrow i from 0.5,0,zi to 3,0,zi filled head lc "blue" lw 4
  # x < 0 側
  set arrow i+10 from -1.2,0,zi to -3,0,zi filled head lc "blue" lw 4
}

set title "一様な面電荷による電場"
splot 0, u, v lc "black" lw 1 notitle

Out[11]:
	dummy variable is t for curves, u/v for surfaces
In [12]:
set output "./vec-fig-men.svg"
replot
set output

直線電流による磁場 $\boldsymbol{B}$

plot 版

In [13]:
%gnuplot inline svg size 480,480 fixed enhanced font 'Noto Sans CJK JP,16'
In [14]:
reset
unset xtics
unset ytics
unset border
set zeroaxis
set xrange [-3.5:3.5]
set yrange [-3.5:3.5]
set size ratio 1

# 磁力線用座標ファイル
set print "xy-line.txt"
do for [i=0:6] {
  # 電線に近い方が密となるようにしてみる
  r = 3*(2./3)**i
  do for [j=0: 64] {
    theta = pi/32 * j
    x = r * cos(theta)
    y = r * sin(theta)
    print sprintf("%8.4f  %8.4f", x, y)
  }
  # with lines するので適宜空行をいれる
  print " "
}
set print

# ベクトルの始点用座標ファイル
# 矢印の混雑防止のため適宜間引く
set print "xy.txt"
do for [i=0:5] {
  r = 3*(2./3)**i
  do for [j=0: 3] {
    theta = pi * j + i*pi/2
    x = r * cos(theta)
    y = r * sin(theta)
    print sprintf("%8.4f  %8.4f", x, y)
  }
}
set print

# 磁場 B
Bx(x, y) = - y/(x**2 + y**2)
By(x, y) =   x/(x**2 + y**2)
# 磁場 B の大きさ
B(x, y) = sqrt(Bx(x, y)**2 + By(x, y)**2)
# 規格化された hat B
hBx(x, y) = Bx(x, y)/B(x, y)
hBy(x, y) = By(x, y)/B(x, y)

# 矢印の長さ調整用
scaling = 0.1

# 電線の位置
set label 1 point pt 6 ps 1.5 lc "red" at 0, 0
set label 2 center at first 0,0 "・" tc "red"

set title "直線電流による磁力線"
unset colorbox
# green.pal 改
# https://github.com/Gnuplotting/gnuplot-palettes
set palette defined (\
  0 '#E5F5E0',\
  1 '#C7E9C0',\
  2 '#A1D99B',\
  3 '#74C476',\
  4 '#41AB5D',\
  5 '#238B45',\
  6 '#005A32',\
  7 '#004000' \
)
plot "xy-line.txt" u 1:2:(log(B($1,$2))) w l lw 3 lc palette notitle,\
  "xy.txt" u 1:2:(scaling*hBx($1,$2)):(scaling*hBy($1,$2)):(log(B($1,$2))) \
  w vec notitle lw 7 lc palette filled head

In [15]:
set output "./vec-fig07d.svg"
replot
set output

splot 版

In [16]:
%gnuplot inline svg size 640,640 fixed enhanced font 'Noto Sans CJK JP,16'
In [38]:
reset
unset xtics
unset ytics
unset ztics
#unset border
set zeroaxis
set xrange [-3.5:3.5]
set yrange [-3.5:3.5]
#set size ratio 1

# 磁力線用座標ファイル
# 電流の後ろ側
set print "xyz-linep.txt"
do for [z=0:2] {
  do for [i=0:6] {
    r = 3*(2./3)**i
    do for [j=0: 32] {
    theta = pi/32 * j
      x = r * cos(theta)
      y = r * sin(theta)
      print sprintf("%8.4f  %8.4f  %8.4f", x, y, z)
      # with lines するので適宜空行をいれる
      print ""
    }
    # with lines するので適宜空行をいれる
    print ""
  }
  # with lines するので適宜空行をいれる
  print ""
}
set print

# 磁力線用座標ファイル
# 電流の前側
set print "xyz-linem.txt"
do for [k=0:2] {
  do for [i=0: 6] {
    r = 3*(2./3)**i
    do for [j=32: 64] {
    theta = pi/32 * j
      x = r * cos(theta)
      y = r * sin(theta)
      z = 1*k
      print sprintf("%8.4f  %8.4f  %8.4f", x, y, z)
      # with lines するので適宜空行をいれる
      print ""
    }
    # with lines するので適宜空行をいれる
    print ""
  }
  # with lines するので適宜空行をいれる
  print ""
}
set print

# ベクトルの始点用座標ファイル
# 矢印の混雑防止のため適宜間引く
set print "xyz.txt"
do for [k=0:2] {
  do for [i=0:5:2] {
    r = 3*(2./3)**i
    do for [j=0:3] {
      theta = pi * j + i*pi/2 
      x = r * cos(theta)
      y = r * sin(theta)
      z = 1*k
      print sprintf("%8.4f  %8.4f  %8.4f", x, y, z)
    }
  }
}
set print

# 磁場 B
Bx(x, y) = - y/(x**2 + y**2)
By(x, y) =   x/(x**2 + y**2)
# 磁場 B の大きさ
B(x, y) = sqrt(Bx(x, y)**2 + By(x, y)**2)
# 規格化された hat B
hBx(x, y) = Bx(x, y)/B(x, y)
hBy(x, y) = By(x, y)/B(x, y)

# 矢印の長さ調整用
scaling = 0.1

# 電線
# set arrow 1 from 0,0,0 to 0,0,3.8 lw 3.5 lc "red" lt 4 filled head
set print "densen.txt"
  print 0,0,0,0,0,2.8 
set print

set title "直線電流による磁場"
unset colorbox
# green.pal 改
# https://github.com/Gnuplotting/gnuplot-palettes
set palette defined (\
  0 '#E5F5E0',\
  1 '#C7E9C0',\
  2 '#A1D99B',\
  3 '#74C476',\
  4 '#41AB5D',\
  5 '#238B45',\
  6 '#005A32',\
  7 '#004000' \
)
set key inside sample 1
set view 70,,,
set xyplane 0
set zrange [0:3]
splot \
  "xyz-linep.txt" u 1:2:3:(log(B($1,$2))) w l lw 3 lc palette notitle,\
  "densen.txt" w vec lc "red" lw 4 filled head title "電流", \
  "xyz.txt" \
    u 1:2:3:(scaling*Bx($1,$2)):(scaling*By($1,$2)):(0):(log(B($1,$2))) \
    w vec title "磁場" lw 6 lc palette filled head, \
  "xyz-linem.txt" u 1:2:3:(log(B($1,$2))) w l lw 3 lc palette notitle

In [18]:
set output "./vec-fig08d.svg"
replot
set output