表皮効果とベッセル関数

はじめに

一般的に送電においては円筒状の銅線が用いられるが、高周波電流を流す場合は、渦電流により表面に電流が流れやすくなる表皮効果の影響により、銅線の抵抗値が増加してしまう。今回は、天下り的ではあるが、ベッセル関数を用いて表皮効果の影響についてPythonを用いてグラフ化することを試みる。(以下の動画は、私のgithubのページに置いてあるプログラムを動かすと出てきます。この様に、周波数を上げるほど、表面付近の電流密度が上昇していくことが分かります)

表皮効果の定量化

まず、モデルとしては、円筒座標を用いて解析するわけだが、ここでは以下の記事の解のみを引用する。

円筒導体における電流分布と表皮効果

電流密度Jz(r)はベッセル関数を使い厳密に解くと、以下の様に表せる。

$J_z(r)=\frac{KI}{2\pi r_0}\frac{J_0(kr)}{J_1(kr_0)}(\because k=\sqrt{-j\sigma \mu \omega })$

ただし、σ:電気伝導率,μ:透磁率,ω:角速度,I:総電流,r:半径,r0:総半径

これをpythonプログラムとして書くと以下の様になる。

Python

import matplotlib.pyplot as plt
import math 
import numpy as np
#ベッセル関数の読み込み
from scipy.special import jv

#円筒銅線の材料パラメーター

#電気伝導率σ
sigma = 1.0
#透磁率μ
mu = 1.0

#角周波数
omega = 10
#電流値
I=1.0
#銅線の半径
r_0=100
k= (-1j*omega*mu*sigma)**0.5
r = np.linspace(-1*r_0,r_0,100)
#電流密度（rの関数）
y = (k*I/(2*math.pi))*jv(0,k*r)/jv(1,k*r_0)

plt.plot(r,y)
plt.show()

結果のグラフ表示

ωを1から100まで動かしたときの結果をグラフ表示すると以下の様になる。

ω=100のとき

表面付近に電流が集中していることが分かる。

ω=10のとき

ω=1のとき

表面付近に電流は集中はしているが影響は小さい。

まとめ

周波数が大きい場合は、電流が銅線の表面付近に流れやすくなる表皮効果の影響が無視できなくなってしまうということが、ベッセル関数を用いた解析モデルによって分かった。ゆえに、流れる電流路の断面積が小さくなるので抵抗値が増加してしまうので、交流機器を作る時はそのような影響も考慮する必要性がある。

参考:円筒表示

計算量は格段に大きくなるのだが、円筒座標にして対称性を考慮せずにコンター図を作成するプログラムを書くと以下の通りである。

import matplotlib.pyplot as plt
import math 
import numpy as np
#ベッセル関数の読み込み
from scipy.special import jv
#3dデータを上手く加工する
from scipy.interpolate import griddata



#円筒銅線の材料パラメーター

#電気伝導率σ
sigma = 1.0
#透磁率μ
mu = 1.0

#角周波数
omega = 10
#電流値
I=1.0
#銅線の半径
r_0=10

k= (-1j*omega*mu*sigma)**0.5

#プロット

NUM=10000                          # サンプリング数
x=0                                     # 点のx値
y=0                                     # 点のy値                                  
#zを求めるために使用
i=0   



#a<r<b つまり -r_0<r<r_0の範囲で乱数を生成する
a=-1*r_0
b=r_0

x_a=[]
y_a=[]
z_a=[]
                                  # ループカウンタ
while i<NUM:                            # NUM回ループ
    x=(b - a) * np.random.rand() + a                  # ランダムな点の位置x a<x<bにおけるランダムな点を作成
    y=(b - a) * np.random.rand() + a                  # ランダムな点の位置y
    if x**2+y**2<r_0**2:                      # 円筒内部ならば
 
        z= (k*I/(2*math.pi))*jv(0,k*(x**2+y**2)**0.5)/jv(1,k*r_0)
        x_a.append(x)
        y_a.append(y)
        z_a.append(z)
    else:                               # 円の外
        #plt.scatter(x,y,color='blue')
        z=0   
    i+=1                                # ループカウンタ増


#データの整理
x_new, y_new = np.meshgrid(np.unique(x_a), np.unique(y_a))
z_new = griddata((x_a, y_a), z_a, (x_new, y_new))


#コンター図で表示する
plt.contourf(x_new, y_new, z_new)
#カラーバー
plt.colorbar()
plt.show()

これにより作成したコンター図は以下の様になる。