長岡工業高等専門学校研究紀要

論第 54 巻（2018）

文

三次元有限要素法による磁束密度分布の解析－ クーロンゲージの仮定の有無に関する解析結果の比較 －鋤柄 あかね1・倉橋 貴彦1・井山 徹郎2
1 長岡技術科学大学（Department of Mechanical Engineering, Nagaoka University of Technology）
2 機械工学科（Department of Mechanical Engineering, National Institute of Technology, Nagaoka College）

Numerical Analysis for Distribution of Magnetic Flux Density
Based on the Finite Element Method in Three Dimensions- Comparison of Numerical Results in Case ofwith and without Assumption of Coulomb Gauge Akane SUKIGARA1, Takahiko KURAHASHI1 and Tetsuro IYAMA2
Abstract
In this paper, we present the numerical analysis of the magnetic density distribution based on the three
dimensional finite element method. In this study, the comparison of the numerical results was carried out
under the conditions with and without assumption of the coulomb gauge. Consequently, we confirmed thatappropriate magnetic density distribution can be obtained by considering the condition of the coulomb
gauge.
Key Words : three dimensional model, finite element method, magnetic field analysis, coulomb gauge

１．はじめにを踏まえ，本論文では，クーロンゲージの仮定の有無による磁束密度分布の解析結果の比較を行う．
1)

近年，モータ設計には磁場解析技術 （図-1）が用いられており，電気自動車の普及により需要が高まることが予想される．磁界解析は，磁気ベクトルポテンシャルを導入した手法が明らかにされているが，磁気ベクトルポテンシャルの解の一意性が保証されない場合，実際の現象と異なる磁束密度分布の乱れが生じるとされている 2)．そこで解の一意性を保証するための条件として，磁気ベクトルポテンシ
ャルの発散を0とするクーロンゲージが用いられる．
一方，解析領域が空気と磁性体で構成されるような透磁率が一定でない場合，クーロンゲージの仮定は不適切であるといった指摘もあり 3)，クーロンゲージの有無について議論がなされている．この様な点

- 6 -

図-1

ANSYS による磁場解析結果の一例

鋤柄 あかね・倉橋 貴彦・井山 徹郎

２．三次元有限要素法による定式化
2.2
磁束密度分布の解析に関する定式化の説明を行うために，式(1)に示すマクスウェルの方程式を導入する．

クーロンゲージを仮定しない場合

式(4)の両辺に対し重み関数ベクトル w を乗じて，要素領域 Ωe に対し積分すると，式(7)に示す重み付
き残差方程式が得られる．

∇×𝑯=𝑱+∇×𝑬=−𝜕𝑫𝜕𝑡𝜕𝑩𝜕𝑡∫Ω 𝒘 ⋅ (∇ × (𝜈∇ × 𝑨))𝑑Ω = ∫Ω 𝒘 ⋅ 𝑱0 𝑑Ω(1)𝑒∇⋅𝑩=0∇⋅𝑫=𝜌(7)𝑒

式(7)左辺に式(8)に示すベクトル公式 3)を用いる．
ここに，Hは磁界の強さ，Eは電界の強さ，Bは磁∇ ⋅ (𝒂 × 𝒃) = 𝒃 ⋅ (∇ × 𝒂) − 𝒂 ⋅ (∇ × 𝒃)

束密度，Dは電束密度，Jは電流密度を示す．静磁界解析においてはマクスウェル方程式の時間微分項を無視する．式(2)に示す構成方程式を導入する．
𝑩 = 𝜇𝑯(8)式(8)を用いて式(7)の左辺を展開すると式(9)が得られる．

(2)∫Ω 𝜈(∇ × 𝒘) ⋅ (∇ × 𝑨)𝑑Ω + ∫Ω 𝜈∇ ⋅ ((∇ × 𝑨) × 𝒘)𝑑Ω𝑒𝑒

ここに，μ は透磁率を示す．また，式(3)で定義され= ∫Ω 𝒘 ⋅ 𝑱0 𝑑Ω (9)る磁気ベクトルポテンシャル A を導入する．
𝑩=∇×𝑨𝑒(3)式(2), (3)をマクスウェル方程式（式(1)）の第 1 式に式(9)の左辺第 2 項にガウスの発散定理を適用すると，式(10)が得られる．

代入すると支配方程式（式(4)）が得られる．
𝜈 ∫Ω (∇ × 𝒘) ⋅ (∇ × 𝑨)𝑑Ω + 𝜈 ∫Γ ((∇ × 𝑨) × 𝒘) ⋅ 𝒏𝑑Γ𝑒

∇ × (𝜈∇ × 𝑨) = 𝑱0𝑒(4)= ∫Ω 𝒘 ⋅ 𝑱0 𝑑Ω (10)𝑒

ここに，ν は磁気抵抗率を示し，透磁率μの逆数である．J0 は領域に流入する外部電流を示す．

2.1

ここで，式(11)に示すベクトル公式 3)を用いる．

境界条件𝒂 ⋅ (𝒃 × 𝒄) = 𝒃 ⋅ (𝒄 × 𝒂) = 𝒄 ⋅ (𝒂 × 𝒃)

文献3)より，境界条件をディリクレ境界（第1種境界）およびノイマン境界（第2種境界）に対して，それぞれ式(5)，(6)のように定義する．

(11)

式(11)を用いて式(10)の左辺第二項を展開すると式(12)が得られる．

𝑨 × 𝒏 = 𝒄1 , ∇ ⋅ 𝑨 = 𝑐2 𝑜𝑛 Γ1(5)𝜈 ∫Ω (∇ × 𝒘) ⋅ (∇ × 𝑨)𝑑Ω + 𝜈 ∫Γ 𝒘 ⋅ (𝒏 × (∇ × 𝑨))𝑑Γ(∇ × 𝑨) × 𝒏 = 𝒄3 , 𝑨 ⋅ 𝒏 = 𝑐4 𝑜𝑛 Γ2(6)= ∫Ω 𝒘 ⋅ 𝑱0 𝑑Ω (12)

Γ1 はディレクレ境界，Γ2 はノイマン境界を示す．
c1~c4が0のとき，磁界はディリクレ境界においては境界面に平行であり，ノイマン境界においては垂直𝑒𝑒𝑒

式(12)の左辺第 2 項を移行し，式(6)を代入すると，式(13)が得られる．に入射する．また，nは外向きの単位法線ベクトルを示す．

- 7 -三次元有限要素法による磁束密度分布の解析

式(19)を用いて式(4)左辺を展開すると式(20)となる．

𝜈 ∫Ω (∇ × 𝒘) ⋅ (∇ × 𝑨)𝑑Ω𝑒= ∫Ω 𝒘 ⋅ 𝑱0 𝑑Ω + 𝜈 ∫Γ 𝒘 ⋅ 𝒄3 𝑑Γ𝑒𝜈∇(∇ ⋅ 𝑨) − 𝜈∇2 𝑨 = 𝑱0(13)𝑒(20)式(18)より式(20)の左辺第 1 項は消去でき，クーロンここで，重み関数ベクトルwと磁気ベクトルポテンゲージを仮定した支配方程式（式(21)）が得られる．

シャルAを四面体一次要素により式(14)のように近−ν∇2 𝑨 = 𝑱0似する．
𝑨 = ∑4𝑖=1 𝑁𝑖 𝑨∗𝑖{𝒘 = ∑4𝑖=1 𝑁𝑖 𝒘𝑖∗

(14)

(21)

式(21)の両辺に対して重み関数ベクトル w を乗じて，要素領域 Ωe に対し積分すると，式(22)に示す重み付き残差方程式が得られる．

式(13)に式(14)に示す補間関数を代入すると式(15)が得られる．

∫Ω 𝒘 ⋅ (−𝜈∇2 𝑨)𝑑Ω = ∫Ω 𝒘 ⋅ 𝑱0 𝑑Ω𝑒𝜈 ∫Ω (∇ × (∑4𝑖=1 𝑁𝑖 𝒘𝑖∗ )) ⋅ (∇ × (∑4𝑖=1 𝑁𝑖 𝑨∗𝑖 ))𝑑Ω𝑒(22)

ここで，式(22)左辺に対してグリーンの定理を適用= ∫Ω (∑4𝑖=1 𝑁𝑖 𝒘𝑖∗ ) ⋅ 𝑱0 𝑑Ω + 𝜈 ∫Γ (∑4𝑖=1 𝑁𝑖 𝒘𝑖∗ ) ⋅ 𝒄3 𝑑Γ (15)𝑒𝑒𝑒

す る と ， 式 (22) は 式 (23) の よ う に 変 形 で き る ．
−𝜈 ∫Γ 𝒘 ⋅ (∇𝑨 ⋅ 𝒏)𝑑Γ + 𝜈 ∫Ω ∇𝒘 ⋅ ∇𝑨𝑑Ω𝑒

重み関数の任意性より両辺任意変数w*i を消去し積𝑒= ∫Ω 𝒘 ⋅ 𝑱0 𝑑Ω(23)𝑒分を実行することで，式(16)に示す有限要素方程式が得られる．

式(23)の左辺第1項を移行し，式(6)を代入すると，
𝐴𝑥1𝐽0𝑥1𝐽0𝑦𝐴𝑦11𝑉𝜈𝑉𝜈𝑉𝑒 [𝐵] 𝐴𝑧1 = 𝑒 𝐽0𝑧 + 𝑒 ∇𝑐4 144⋮⋮⋮{ 1}{ 𝐽0𝑧 }{ 𝐴𝑧4 }

式(24)が得られる．
(16)𝜈 ∫Ω ∇𝒘 ⋅ ∇𝑨𝑑Ω = ∫Ω 𝒘 ⋅ 𝑱0 𝑑Ω + 𝜈 ∫Γ 𝒘 ⋅ (∇𝑐4 )𝑑Γ𝑒𝑒𝑒(24)ただし，行列[B]の成分は式(17)に示す．
ここで，重み関数ベクトルwと磁気ベクトルポテン
𝜕𝑁1 𝜕𝑁1 𝜕𝑁1 𝜕𝑁1+𝜕𝑧 𝜕𝑧𝜕𝑦 𝜕𝑦𝜕𝑁1 𝜕𝑁1−𝜕𝑥 𝜕𝑦
𝜕𝑁1 𝜕𝑁1−𝜕𝑥 𝜕𝑧⋮
𝜕𝑁4 𝜕𝑁1−[𝜕𝑥 𝜕𝑧𝜕𝑁1 𝜕𝑁1−𝜕𝑦 𝜕𝑥
𝜕𝑁1 𝜕𝑁1 𝜕𝑁1 𝜕𝑁1+
𝜕𝑧 𝜕𝑧𝜕𝑥 𝜕𝑥
𝜕𝑁1 𝜕𝑁1−𝜕𝑦 𝜕𝑧⋮
𝜕𝑁4 𝜕𝑁1−𝜕𝑦 𝜕𝑧

𝜕𝑁1 𝜕𝑁1−𝜕𝑧 𝜕𝑥
𝜕𝑁1 𝜕𝑁1−𝜕𝑧 𝜕𝑦𝜕𝑁1 𝜕𝑁1 𝜕𝑁1 𝜕𝑁1+
𝜕𝑦 𝜕𝑦𝜕𝑥 𝜕𝑥⋮
𝜕𝑁4 𝜕𝑁1 𝜕𝑁4 𝜕𝑁1+
𝜕𝑦 𝜕𝑦𝜕𝑥 𝜕𝑥

⋯
⋯
⋯
⋱
⋯

𝜕𝑁1 𝜕𝑁4−𝜕𝑧 𝜕𝑥
𝜕𝑁1 𝜕𝑁4−𝜕𝑧 𝜕𝑦
𝜕𝑁1 𝜕𝑁4 𝜕𝑁1 𝜕𝑁4+
𝜕𝑦 𝜕𝑦𝜕𝑥 𝜕𝑥⋮
𝜕𝑁4 𝜕𝑁4 𝜕𝑁4 𝜕𝑁4+𝜕𝑦 𝜕𝑦
𝜕𝑥 𝜕𝑥 ]

(17)

シャルAを四面体一次要素により近似する．式(24)に式(14)に示す補間関数を代入すると式(25)が得られる．
𝜈 ∫Ω ∇(∑4𝑖=1 𝑁𝑖 𝒘∗𝑖 ) ⋅ ∇(∑4𝑖=1 𝑁𝑖 𝑨∗𝑖 )𝑑Ω𝑒= ∫Ω (∑4𝑖=1 𝑁𝑖 𝒘𝑖∗ ) ⋅ 𝑱0 𝑑Ω + 𝜈 ∫Γ (∑4𝑖=1 𝑁𝑖 𝒘𝑖∗ ) ⋅ (∇𝑐4 )𝑑Γ𝑒

2.3

クーロンゲージを仮定する場合

𝑒(25)
重み関数の任意性より両辺任意変数w*i を消去し積式(18)に示すクーロンゲージを仮定する 1)．

分を実行することで，式(26)に示す辺要素有限要素∇⋅𝑨=0(18)方程式が得られる．

ここで，式(19)に示すベクトル公式 3)を用いる．
∇ × (∇ × 𝒂) = ∇(∇ ⋅ 𝒂) − ∇2 𝒂(19)- 8 -鋤柄 あかね・倉橋 貴彦・井山 徹郎∇𝑁1 ⋅ ∇𝑁100𝜈𝑉𝑒⋮[00∇𝑁1 ⋅ ∇𝑁10⋮000∇𝑁1 ⋅ ∇𝑁1⋮∇𝑁4 ⋅ ∇𝑁1𝐴𝑥1⋯
0𝐴𝑦1⋯0⋯ ∇𝑁1 ⋅ ∇𝑁4 𝐴𝑧1⋱⋮⋮⋯ ∇𝑁4 ⋅ ∇𝑁4 ] { 𝐴𝑧4 }𝐽0𝑥
1𝐽0𝑦1𝑉𝜈𝑉= 𝑒 𝐽0𝑧 + 𝑒 ∇𝑐4 1 (26)
44⋮
⋮
{1 }
{ 𝐽0𝑧 }
ここで，式(26)は x, y および z 成分が互いに干渉しないため，成分ごとに分けて計算することができる．

３．数値解析例

図-2

計算モデル図（1/8 モデル）

文献 2)を参照し，電気学会にて制定された標準比較モデルを用いて，クーロンゲージの条件の有無による計算結果の比較を行う．図-2に計算モデルを示
す ． 100mm × 100mm × 高 さ 200mm の 鉄 心 に 内 幅150mm×外幅200mm×高さ100mmの四角形コイルによる電流を与えるとし，1/8モデルにより解析を行う．x=0, y=0面に対しディリクレ境界条件，その他の面に対しては磁界が垂直に入射する面としてノイマン境界条件を与える．ここで，クーロンゲージを仮定しない場合，c3=0であってもc4=0を満たさないため，x=0, y=0面以外に対してA・n=0の条件をノイマン条件として与える必要がある．
表-1に計算条件を示す．比透磁率は物体の透磁率と真空の透磁率の比であり，磁気抵抗率νは物体の透磁率の逆数（ν=1/μ0μr）である．また，電流密度J0はコイルの断面積1m2あたり3000[AT]の電流が流れることを想定している．ここで，[AT]は[A]にコイルの巻き数を表す無次元数をかけたものである．図-3に計算メッシュを示す．
表-1 計算条件

総節点数
総要素数
真空の透磁率μ0 [H/m]比透磁率μr [-]電流密度J0 [A/m2]

2969
122884π×10-71000
1.2×106

図-3

有限要素メッシュ図

図-4～図-6 に磁束密度ベクトルの解析結果を示す．図-4 は磁束密度ベクトルの大きさ，図-5 はy=25mm における内部の磁束密度ベクトルの大きさを示す．また，図-6 は磁束密度ベクトルの z 方向成分を示す．クーロンゲージを仮定しない場合，磁束密度ベクトルの大きさは乱れなく分布しているように見られるが，成分ごとに確認すると磁束密度の分布に乱れが発生している．一方，クーロンゲージを仮定した場合，磁束密度分布に乱れは見られないが，鉄心における磁束密度は文献 2)の値に比べ小さい結果となっている．これは，クーロンゲージの仮定は磁気ベクトルポテンシャルの発散が 0 であり，鉄心と空気の境界において透磁率が変化するにも関わらず磁気ベクトルポテンシャルが連続になることが起因していると考えられる．透磁率の異なる材料の境界を表現するには文献 3)で示される辺要素が有効であるが，境界条件の設定方法が明確ではないという課題もある．

- 9 -三次元有限要素法による磁束密度分布の解析

図-4(a)

磁束密度ベクトルの分布図

図-4(b)（クーロンゲージの仮定をしない場合）

図-5(a)

磁束密度ベクトルの分布図（クーロンゲージの仮定をした場合）

y=25mm 断面における磁束密度ベクトルの図-5(b)

y=25mm 断面における磁束密度ベクトルの分布図（クーロンゲージの仮定をしない場合）

分布図（クーロンゲージの仮定をした場合）

図-6(a) 磁束密度ベクトルの z 方向成分の分布図（クーロンゲージの仮定をしない場合）

図-6(b) 磁束密度ベクトルの z 方向成分の分布図（クーロンゲージの仮定をした場合）- 10 -

鋤柄 あかね・倉橋 貴彦・井山 徹郎

４．おわりに
本研究では，磁気ベクトルポテンシャルを用いた有限要素法による磁界解析において，クーロンゲージの仮定の有無に関する比較を行った．クーロンゲージを仮定しない場合，磁気ベクトルポテンシャルの解の一意性が満たされず，磁束密度分布に乱れが生じることが確認できた．一方，クーロンゲージを仮定した場合，磁束密度分布の乱れもなく解析を行うことができた．また，クーロンゲージを仮定した場合には，x, y および z 方向成分を独立に計算することができるため，計算コストの削減が可能となる．しかし，クーロンゲージを仮定した場合は磁気ベクトルポテンシャルの連続性により透磁率の異なる材料の界面を表現できないという問題が明らかになった．今後，この問題を解決するとされる辺要素有限要素法による磁界解析を検討する．
謝辞：本研究を進めるにあたり，長岡技術科学大学学長戦略的経費（⑤材料科学領域：未来を牽引する高信頼性材料分野）の援助を受けた．また，本論文の解析結果は，九州大学情報基盤研究開発センターの高性能アプリケーションサーバを使用し計算を行ったものである．ここに謝意を表す．
参考文献
1)

CAD CAE 研 究 会 編 ， “ 有 限 要 素 法 解 析 ソ フ トANSYS 工 学 解 析 入 門 ” ， 理 工 学 社 ，pp.190-208，
2005．

2)

金井靖，阿部武雄，仙石正和，飯島泰蔵，飯塚雅博，
武笠幸一，“磁気ベクトルポテンシャルのゲージ条件の検討と三次元静磁界解析のための新しい有限要素定式化”，電気学会論文誌 D，110 巻 3 号，pp.265271，1990．

3)

本間利久，五十嵐一，川口秀樹，“数値電磁力学”，
森北出版，pp.122-125, 237，2002．
(2018. 8. 1 受付)- 11 -