第 2 章 計算方法 47
2.3 計算手法
2.3.1 解析手順
図41に，本計算のフローチャートを示す。まず，速度場について計算を行った後に，全体の 電磁場を求める。これらを用いて，アークプラズマと電極の温度などを解析する。これらの計 算を収束するまで繰り返し計算し，目的とする時間まで行う。
また，2.2節に示した微分方程式を有限体積法を用いて各グリッドに対して離散化を行い，得 られた代数方程式をTDMA（TriDiagonal-Matrix Algorithm）法により算出する(77)。アーク空 間を解析する最小のグリッド長さは，それぞれ0.1 mm× 0.1 mm× 0.1 mmとした。
3.5―3.6，4.1―4.3節，並びに，5.1節の横磁界印加時の解析では，定常解析を行うため，2.2
節に示した支配方程式の非定常項を省き，計算を行う。3.3―3.4，3.7節，並びに，4.4―5.2節 の縦磁界印加時の解析では，非定常解析を行う。この時，3.3―3.4，3.7節では，電流100 Aで，
定常計算を行ったときの電位や温度などの物理量を初期値として計算を行う。また，4.4節以降 の解析では，電流100 A，縦磁束密度を0，10，20 mTでそれぞれ定常解析を行ったときの物理 量を初期値として計算を行う。ただし，溶融池形状の形成過程を解析するため，母材の温度の
初期値は300 Kにしてから計算を行う。
58
Current continuity equation
Electrical potential, Current density Boundary conditions
Convergence?
Yes
Yes No
No
End time?
+Dt
END Start
Calculation structure
Initial conditions
Momentum conservation equation Mass conservation equation Pressure, Velocity
Maxwell ampere equation
Vector potential, Magnetic flux density
Energy conservation equation Enthalpy
Velocity, Temperature Electromagnetic field
Temperature conversion Temperature
Calculate properties (Gas, Electrode)
Convection and diffusion equation Metal vapor concentration
図 41 本計算のフローチャート Fig.41 Flowchart of the calculation.
(1) 速度場
速度場は式(2)の運動量保存式を各方向に展開した以下の式(15)〜(19)をそれぞれ解くこと により求める。式(15)，(16)は，3.3―3.4節で使用し，式(17)，(18)，(19)は，3.5節以降で 用いる。圧力と速度場の解析には質量保存式と運動量保存式を連成して解析する手法である SIMPLER（Semi-Implicit Method for Pressure-Linked Equations Revised）法を用いた(78, 79)。 なお，SIMPLER法の具体的な計算方法は，付録に載せた。連成して式を解く際に必要な質量 密度は，気体の状態方程式とSahaの式を連立して求めた粒子組成の各粒子数密度を基に算出 している。また，外力である電磁力は電流密度と磁束密度の外積から算出した。
· 運動量保存式（軸方向）
∂
∂t(ρvz) + 1 r
∂
∂r(rρvrvz) + ∂
∂z (ρvz2)
= − ∂pr
∂z + ∂
∂z (
2η∂vz
∂z )
+ 1 r
∂
∂r (
rη∂vz
∂z +rη∂vr
∂r )
− 2 3
∂
∂z (
η (∂vz
∂z +1 r
∂(rvr)
∂r ))
+jrBθ+ρg (15)
· 運動量保存式（径方向）
∂
∂t(ρvr) + 1 r
∂
∂r
(rρv2r) + ∂
∂z (ρvrvz)
= − ∂pr
∂r +1 r
∂
∂r (
2rη∂vr
∂r )
+ ∂
∂z (
η∂vz
∂r +η∂vr
∂z )
− 2 3
∂
∂r (
η (∂vz
∂z +1 r
∂(rvr)
∂r ))
− η2vr
r2 −jzBθ (16)
· 運動量保存式（x方向）
∂(ρvx)
∂t + ∂
∂x(ρvxvx) + ∂
∂y(ρvyvx) + ∂
∂z(ρvzvx)
= −∂pr
∂x + { ∂
∂x(η(2∂vx
∂x −2
3div⃗v)) }
+ { ∂
∂y(η(∂vx
∂y + ∂vy
∂x)) }
+ { ∂
∂z(η(∂vx
∂z +∂vz
∂x)) }
+jyBz−jzBy (17)
· 運動量保存式（y方向）
∂(ρvy)
∂t + ∂
∂x(ρvxvy) + ∂
∂y(ρvyvy) + ∂
∂z(ρvzvy)
= −∂pr
∂y + { ∂
∂y(η(2∂vy
∂y − 2
3div⃗v)) }
+ { ∂
∂z(η(∂vy
∂z + ∂vz
∂y )) }
+ { ∂
∂x(η(∂vx
∂y +∂vy
∂x)) }
+jzBx−jxBz (18)
60
· 運動量保存式（z方向）
∂(ρvz)
∂t + ∂
∂x(ρvxvz) + ∂
∂y(ρvyvz) + ∂
∂z(ρvzvz)
= −∂pr
∂z + { ∂
∂z(η(2∂vz
∂z − 2 3div⃗v))
} +
{ ∂
∂x(η(∂vz
∂x +∂vx
∂z )) }
+ { ∂
∂y(η(∂vy
∂z +∂vz
∂y )) }
+jxBy−jyBx+ (ρ0−ρ)g (19) ここで，x[m]，y[m]，z[m]は，奥行方向，径方向および軸方向の位置であり，r[m]は，円 筒座標系における径方向の位置を示す。また，軸方向の運動量保存式では，重力を，z方向の 運動量保存式では，重力と浮力を考慮しており，ρ0は，300 Kのアルゴンの質量密度である。
本研究では，横風に対する抗力と，横磁束密度と縦磁束密度に起因されるローレンツ力を考 慮することで，横風下における外部磁界印加時のアークの解析を行う。まず，横風に対する抗 力は，式（18）に示すY方向の対流項が横風流速によって変化することで考慮される。また，
横磁束密度に起因されるローレンツ力は，外部磁束密度により式（18） に示されるローレンツ 力の項であるjzBxのBxが増加することで考慮される。一方，縦磁束密度に起因されるローレ ンツ力は，外部磁束密度により式（17），（18） に示されるローレンツ力の項であるjyBz，jxBz のBzが増加することで考慮される。
なお，3.5―3.6節，4.1―4.3節，5.1 節の横風吹付けのみの解析や横磁界を印加した解析にお いて，粘性応力は，質量密度と粘性係数を一定と仮定した条件で導出される簡略化されたもの を使用した(80)。この粘性項を簡略化したことによる実験と計算とのアーク偏向距離の差は小 さいと考えられる。これは，陰極径が大きく陰極先端の角度が鈍角であるため，電流密度が先 端部で集中せず，電流密度の低下によりアークの温度と流速が小さくなるためである。
また，溶融池内の対流を計算する際には，以下の式に示す溶融池内の対流を誘起する4つの 駆動力を考慮した。
· ローレンツ力
FLz2D =jrpoolBθpool (20a)
FLr =jzpoolBθpool (20b)
FLx =jypoolBzpool−jzpoolBypool (20c)
FLy =jzpoolBxpool−jxpoolBzpool (20d)
FLz =jxpoolBypool−jypoolBxpool (20e)
· アークの流れによる引きずり力 FD = ∂
∂z (
η∂vr
∂z )
(21a)
FDx = ∂
∂z (
η∂vx
∂z )
(21b) FDy = ∂
∂z (
η∂vy
∂z )
(21c)
· マランゴニ力 FM = ∂
∂z(∂γ
∂T
∂T
∂r) (22a)
FM x= ∂
∂z(∂γ
∂T
∂T
∂x) (22b)
FM y = ∂
∂z(∂γ
∂T
∂T
∂y) (22c)
· 浮力
FB =ρgβV(T −T0) (23)
ここで，式(20a)は，２次元軸対称円筒座標系の軸方向のローレンツ力，式(22a)―(22c)の γは表面張力(N/m)であり，∂γ/∂Tは表面張力の温度特性に該当し，∂γ/∂T は-4.58× 10−4 と した(81)。浮力の式(23)において，T0= 300 K は常温であり，βV = -4.95× 10−51/K は体積膨 張率(82)とした。3.3―3.4節と，4.4節以降における縦磁界印加時の計算において，これらの駆
動力は，SUS304の融点1,750 K以上の陽極領域に，運動量保存則の生成項として与え（ひきず
り力のみは運動量保存則自身に含まれる），溶融池内の対流を解析した。
(2) 金属蒸気挙動の記述
金属蒸気は，陽極材料SUS304の主成分である鉄のみ陽極から発生すると仮定した。鉄蒸気 の拡散係数に関しては，式(13)の粘性近似式を使用した。また，βAr ，βF eはともにβi =(Dti ρi) /ηi で定義される無次元化された定数で，添え字iはアルゴンまたは鉄蒸気であり，添え字 tは他方の元素を示している。βは，理論上ではアルゴン，ヘリウム，水素，窒素，酸素，炭酸 ガスなど様々なガスに対して1.2〜1.543の範囲の数値となるが(83)，多くの実験データに基づ く平均値としてβAr =βF e =1.385を仮定した。この近似式は，室温から30,000 Kの高温に至る まで比較的に精度が高く(84)，本シミュレーションに適していると考えられる。
62
陽極表面上における鉄の蒸気圧は文献(85)を参考に，陽極表面の温度から決定した。陽極表 面上の鉄蒸気濃度は飽和蒸気圧と周囲の圧力によって決定される。飽和蒸気圧は，溶融池表面 の温度によって決定し，陽極表面上の鉄蒸気濃度を求めた後、鉄蒸気濃度分布は式(12)を解い て求めた(86, 87)。
なお，3.5―3.6節，4.1―4.3節, 5.1節の横風吹付けのみの解析や横磁界を印加した解析では，
物理現象の単純化のために，金属蒸気は考慮せず計算を行った。
(3) 電磁場
次に，電磁場の解析方法について述べる。電位は以下の式(24)，(25)のポアソン方程式より 算出する。
1 r
∂
∂r(rσ∂ϕ
∂r) + ∂
∂z(σ∂ϕ
∂z) = 0 (24)
∂
∂x (
−σ∂ϕ
∂x )
+ ∂
∂y (
−σ∂ϕ
∂y )
+ ∂
∂z (
−σ∂ϕ
∂z )
= 0 (25)
式(24)，(25)で求められた電位に対して，式(9)のオームの法則より各方向の電流密度を算 出することができる。式(9)を展開すると以下の式(26a)〜(26d)となる。
jr =−σ∂ϕ
∂r =σEr (26a)
jx =−σ∂ϕ
∂x =σEx (26b)
jy =−σ∂ϕ
∂y =σEy (26c)
jz =−σ∂ϕ
∂z =σEz (26d)
更に，式(10)，(11)よりアーク内に流れる電流による磁束密度分布を求めることができる。
まず，式(10)よりベクトルポテンシャルを求める。式(10)を展開すると以下の式(27a)〜(27c) となる。
∂Ax2
∂x2 +∂Ax2
∂y2 + ∂Ax2
∂z2 =−µjx (27a)
∂Ay2
∂x2 + ∂Ay2
∂y2 + ∂Ay2
∂z2 =−µjy (27b)
∂Az2
∂x2 +∂Az2
∂y2 +∂Az2
∂z2 =−µjz (27c)
次に，式(11)により磁束密度を求める。式(11)を展開すると以下の式(28a)〜(28c)となる。
Bx = ∂Az
∂y −∂Ay
∂z (28a)
By = ∂Ax
∂z − ∂Az
∂x (28b)
Bz = ∂Ay
∂x − ∂Ax
∂y (28c)
なお，3.3―3.4節で使用する２次元軸対称円筒座標系のシミュレーションでは，式(24)，(26a)，
(26d)で，電流密度と電界の算出を行い，以下の式(29)より磁束密度の算出を行う。
1 r
∂
∂r(rBθ) =µ0jz (29)
また，式(25)，(26b)，(26c)，(26d)，(28a)，(28b)，(28c)は，3.5節以降で用いる。
以上のように求められた電流密度や磁束密度を用いて，運動量保存式の外力として電磁力を，
エネルギー保存式の生成項としてジュール熱を代入することによりアークの電磁界の影響を考 慮した速度場，並びに，温度場を求めることができる。
(4) 温度場
式(6)を展開すると式(30)，(31)となる。電極界面では，電極とアークとのエネルギー授受 を考慮しており，詳細は2.3.3項で説明する。
また，式(30)や，(31)を解いた後に，エンタルピーを温度に変換している。
· ２次元軸対称円筒座標系のエネルギー保存則
∂
∂t(ρh) + 1 r
∂
∂r(rρvrh) + ∂
∂z(ρvzh)
= 1 r
∂
∂r (rκ
Cp
∂h
∂r )
+ ∂
∂z ( κ
Cp
∂h
∂z )
+ (jr2
σ + jz2 σ
)
−U (30)
· ３次元直交座標系のエネルギー保存則
∂(ρh)
∂t + ∂
∂x(ρvxh) + ∂
∂y(ρvyh) + ∂
∂z(ρvzh)
= ∂
∂x ( κ
Cp
∂h
∂x )
+ ∂
∂y ( κ
Cp
∂h
∂y )
+ ∂
∂z ( κ
Cp
∂h
∂z )
+ jx2 σ + jy2
σ + jz2
σ −U (31)
64