鉄鋼材料のデンドライト樹間液相流動の透過率の定量化
京都工芸繊維大学 機械工学系 教授 高木 知弘（代表研究者）
京都工芸繊維大学 機械工学系 助教 坂根 慎治（共同研究者）
1．緒言
鉄鋼材料の鋳造においてマクロ偏析の予測および制御は極めて重要であり，数値シミュレーション研究が多く行われてきた[1-5]．マクロ偏析の形成過程は複雑であり，その数値シミュレーションでは多くの現象を考慮する必要がある．なかでも固液共存域つまりデンドライト樹間の液相流動は最も重要な現象の一つである．
マクロ偏析シミュレーションにおいて，デンドライト樹間の液相流動の容易さは，透過率(permeability)によって表現される．一般に，透過率は Darcy 則[6]1U
 K P

(1)を用いて求められる．ここで，U は平均流速ベクトル，μ は粘性係数，K は透過率テンソル，P は平均圧力勾配ベクトルである．透過率評価は 1960 年ごろから行われ，当初は実際の合金の凝固過程でデンドライトを成長させ，そこに強制対流を流すことで実測されていた[7-13]．しかしながらこの方法では液相流動時にデンドライト形態が変化し，また具体的なデンドライト形態情報を得ることができない課題があった．その後，デンドライト形態を X 線観察やシリアルセクショニングによって取得し，その形態に対して樹間流動シミュレーションを行う透過率評価が行われた[14-19]．この手法では，実際の凝固組織を用いた安定な透過率評価が可能であるが，特別な装置が必要となり系統的な透過率予測は困難であった．1990 年ごろにフェーズフィールド(phase-field: PF)法が現れ，滑らかな表面を有するデンドライトの再現が可能となり，PF 法と液相流動計算による透過率評価が試みられた．しかしながら，PF シミュレーションは計算コストが高く，デンドライトの一部を用いた透過率評価しかできなかった[20-22] ．研究者のグループでは，graphicsprocessing unit (GPU)を多数搭載したスパコンを用いた大規模 PF シミュレーションによって複数の柱状デンドライト成長の再現に成功し[23]，ここで得られた組織に対して格子ボルツマン(lattice Boltzmann: LB)法によって樹間液相流動を計算することで，高精度な透過率予測法を開発した[24]．これを用いて，柱状デンドライトの垂直方向と平行方向の無次元透過率が，Kozeny-Carman (KC)の式で精度よく表現できることを示した[24, 25]．
KC 式は，

1  f S KS 2SkC f S2－ 123 －3(2)によって表される[26]．ここで，fs は固相率，SS は固液界面積 A S L を固相体積 VS で除して得られる比界面積(SS = ASL /VS )，kc は KC 係数である．PF 法と LB 法を用いた研究によって，柱状デンドライトに垂直な流れでは kc = 9 [24]，柱状デンドライトに垂直な流れでは kc = 3 [25]で無次元透過率 KSS 2 を高精度に予測可能であることを示した．一方で，これらの評価は Al-3wt%Cu の特定の条件での評価であり，この結果がどの程度一般性があるのかを確かめる必要がある．そのためには，デンドライト成長と樹間液相流動の計算を高速化し，様々なデンドライト形態に対して系統的な評価を行う必要がある．
本研究では，適合格子細分化法(adaptive mesh refinement: AMR)を複数の GPU を用いた並列計算に適用する parallel-GPU AMR を開発し，透過率評価を高速化した．また，これを用いた Fe-Si 合金に対する系統的な柱状デンドライトに対する透過率評価を行った．
2．透過率評価法の高速化
透過率評価は，
手順 1: PF シミュレーションによるデンドライト組織予測
手順 2: LB シミュレーションによる樹間液相流動計算
手順 3: Darcy 則による透過率計算の流れで行う[24]．手順 1 では，文献[24]と同じ２元合金凝固の定量的 PF モデル[27]を用いた．手順 2 では，文献[24]で用いた単一緩和時間 LB モデルではなく，二緩和時間 LB モデルを用い幅広い粘性係数の適用を可能とした[28]．また，固液界面に interpolated
bounce-back を適用し界面をシャープに表現した[28]．液相流動計算においては，流動方向に周期境界を適用し，その方向に一定圧力勾配を与えることで液相流れを生じさせた[24,25]．最終的に，手順 2 で得られた流動場から平均流速を求め，式(1)を用いることで透過率を計算する．
実験的に得られた組織に対して液相流動計算を行う場合，一般には流入流出側に液相のみのバッファー領域を設ける必要があり，周期境界の適用はできない[29]．一方，本研究で用いる手法では，組織と流動を数値計算によって求めるため，周期境界の適用が可能であり，条件設定の自由度が高いメリットを有する．
デンドライト成長の PF シミュレーションと樹間液相流動の LB シミュレーションを高速化するために，parallel-GPU AMR を開発した．Parallel-GPU AMR は，まず溶質拡散のみを考慮したデンドライト成長問題において開発し[30]，それを液相流動と固体運動を伴う PF 法と LB 法を連成したモデルに拡張し[31]，これを透過率評価に適用した．本研
究で構築した透過率評価に対する parallel-GPU AMR では，ブロック AMR に対してさらにマザーリーフ法[32]を適用し，ブロック間通信量を減らすことで約 1.5 倍の計算高速化を達成した[33]．1.5 倍の高速化は計算規模が小さい場合はあまり有効ではないが，本研究で対象とする大規模計算では大きな成果である．

－ 124 －

3．柱状デンドライト成長シミュレーション
Fe-5.4wt%Si 合金を対象にした柱状デンドライトに対する透過率評価を行った．デンド
ライト成長の PF シミュレーションでは，文献[30]に示す物性値を用いた．温度勾配を一定とする frozen temperature approximation と，デンドライト先端がある座標に到達したあとは moving frame algorism を用いる一方向凝固シミュレーションを行うことで柱状デンドライト構造を得た[23]．温度勾配は G = 3 K/mm，冷却速度は R = 5 K/min と設定した．これらの条件において，まず，準備計算によって格子サイズをΔx = 1.25 µm と決定した．先の研究において，柱状デンドライトは六角形配列が支配的となること[23]，複数の柱状デンドライト組織の樹間液相流動の透過率は周期的な正六角形配列で高精度に再現できることを確認しており，計算コスト低減と柱状デンドライト成長方向への全固相率達成を目的として，柱状デンドライトの六角形周期配列を採用した．図 1(a)に正六角形配列を示し，この周期配列を再現する最小領域を右図に示している．
この領域を用いる場合，x，y 方向いずれもミラー対称境界を設定する．ここで，λは一次枝間隔，Lx は領域の x 方向サイズ，Ly は領域の y 方向サイズである．ブロック分割 AMR を用いる場合，領域サイズが最大ブロックサイズに依存するため，正六角形配列では Lx : Ly = 1 :

3 の領域サイズ比のところを Lx : Ly = 1 : 2 とした．また，成長方向は全固相率を達成するために，液相線と固相線の間隔から決まる凝固長と，一次枝先端前方の拡散長の合計に少し余裕を持たせ，Lz = 8.4 mm (nz = 8.4/Δx = 6,720 grids)とした．

図 1 柱状デンドライトの六角形配列
様々な凝固形態に対する透過率評価を行うために，一次枝間隔λ を変化させることを考えた．図 1(a)右図の条件で２本の柱状晶を成長させる場合，λ が小さすぎると 1 本が淘汰されて 1 本だけの一次枝が成長し，λ が大きすぎると枝分かれによって 2 本よりも多い一次枝が成長する．そこで，表 1 に示すような Case 0 – 4 の 5 つの領域サイズに対して予備計算を行った．結果として，Case 0 では淘汰が生じ，Case 4 では枝分かれが生じた．このため，2 本が安定して成長した Case 1 – 3 を採用した．

－ 125 －

表 1 一次枝間隔と領域サイズ

λ[mm]

Lx [mm]Ly [mm]nx [grids]

ny [grids]

Case 01798016064128

Case 1358160320128256

Case 2537240480192384

Case 3716320640256512

Case 4894400800320640

図 2 に PF シミュレーションの結果得られた柱状組織形態を示す．図 2(a)は全体を x 方向から見た結果，図 2(b)は Case 1 に示す位置の拡大図である．一次枝間隔λ を変化させることで，セル形態(Case 1)からデンドライト形態(Case 2, 3)を表現することができた．また，図 2(b1)のデンドライト先端部のみブロック分割を示している．固液界面にのみ細かいブロックが配置されていることがわかる．ブロックのレベルは 4 段階とし，1 ブロックに83 格子を割り当てた．このため，図 2(b1)の一つの正方形の中に 8 × 8 の格子が入っている．また，4 ブロックを用いた場合の最小格子と最大格子の 1 辺のサイズ比は 1 : 8 となる．図 2(b1)のブロック分割において，ブロックの中に一つ粗いブロックが表示されるのがマザーリーフ法の特徴である．なお，本シミュレーションでは最大で GPU (NVIDIA TeslaA100)を 4 基並列して計算した．GPU の性能は異なるが，先の研究では最大で 128 GPUを並列しており[24]，parallel-GPU AMR の適用によって研究室の GPU サーバでの計算が可能となった．

図 2 一方向凝固 PF シミュレーションにより得られた t = 7 × 106Δt 時の柱状組織（時間増分はΔt = 5.07 × 10-5 s）－ 126 －図 3 は (a)固相率 fs と(b)界面密度 SV の図 2(a)に示すデンドライト先端から下向きの変化である．界面密度 SV は，ある領域内の固液界面積 ASL を領域体積 V で除した量として定義される(SV = ASL /V )．ここでは z 方向に，最大ブロックサイズに対応する 64Δx 毎にfs と SV を算出した．固相率 fs は単調に増加し，ほぼ全固相率を網羅できていることがわかる．界面密度 SV は，いずれもはじめ急激に増加するが，Case 1 では 0.5 mm あたりで変化が小さくなりほぼ直線的に増加し，6 mm あたりで減少した．Case 2 では，1 mm あたりから減少に転じている．Case 3 では，3 mm あたりで最大値をとり，その後減少した．
このように，固液界面の面積は柱状晶の形態と二次枝の発展に依存してその変化の傾向が異なることがわかる．

図 3 図 2(a)に示す柱状組織の(a)固相率 fs と(b)界面密度 SV のデンドライト先端位置から z 方向下向きの変化
4．透過率評価シミュレーション
4. 1. 垂直方向流れ
図 2(a)の柱状晶に垂直方向（x 方向と y 方向）の液相流動を与えるために，図 4 に示す手順に従って LB シミュレーションで用いるモデルを作成した．圧力勾配を与えることで流速を生じさせるため，まずデンドライト先端上部の液相のみの領域を削除する．

－ 127 －

図 4 垂直流れに対する透過率評価モデルの作成方法
次いで，流れ方向に領域の複製と反転を行い，２つの領域を結合した．このように作成したモデルに対して，流れ方向は周期境界，それに垂直方向はミラー対称境界を設定した．
樹間液相流動を計算する LB シミュレーションにおいても AMR を適用しているが，格子サイズを PF シミュレーションの倍Δx LB = 2Δx とすることで，さらに計算コストを低減した．なお，Δx LB = Δx とする計算と結果は同じであることを確認している．LB シミュ
レーションは液相流動が定常状態に達するまで行われた．
図 5 に，図 2(a)の柱状晶に x 方向と y 方向の流動を与えた際の，無次元透過率(a) KxSS 2，(b) KySS 2 と固相率 fS の関係を示す．また，黒の破線で kc = 9 とした式(2)の KC 式を示している．図 5(a)の x 方向流れに関しては，Case 1–3 のいずれも KC 式と比較的よく一致している．詳細にみると，fS < 0.6 の領域において，Case 1 は KC 式に一致しているが，Case 2，Case 3 と二次枝が発達するほど少しずつ KC 式から上側に外れる傾向にあることが確認できる．図 5(b)の y 方向流れに関しては，Case 1–3 はよく一致しているが，fS >0.5 の領域で KC 式から外れる傾向にある．特に fS > 0.7 では急激に低下した．これは，図2(b3)に示すように y 方向の流路が閉ざされたことによる．

－ 128 －図 5 垂直流れに対する無次元透過率 KSS 2 と固相率 fS の関係(a) x 方向流れ (b) y 方向流れ
4. 2. 平行方向流れ
4. 1.節で示した垂直流れの場合は，図 2(a)の全領域を対象とすることができたが，平行方向流れにおいては，流れ方向に固相率が変化するため，領域を取り出して計算を行う必要がある[25]．図 6 に平行方向流れの LB シミュレーションに用いるモデル作成方法を示す．まず，対象とする一部の領域を抽出し，それを複製し z 方向に反転し，最後に２つを結合した．この方法で，固相率を変化させて 11 個のモデルを作成した．これらのモデルに対して z 方向の周期境界のもとで圧力勾配を与えた．他の方向はミラー対称境界としている．

図 6 平行流れに対する透過率評価モデルの作成方法－ 129 －

図 7 に平行方向流れに対する無次元透過率 KzSS 2 と固相率 fS の関係を示す．黒の破線
は，式(2)の KC 式において kc = 3 としたものである．全体的には，いずれも KC 式とよく一致している．詳細にみると，fS < 0.6 の領域において，Case 1 は KC 式より少し小さく，二次枝が発達するほど値が大きくなり，Case 3 は KC 式より高い値となっている．この傾向は図 5(a)の x 方向流れと似ている．

図 7 平行流れに対する無次元透過率 KzSS 2 と固相率 fS の関係
5．考察
図 5 の垂直流れでは，x 方向流れは kc = 9 の KC 式より少し高く，y 方向流れは KC 式より低くなる傾向が得られた．実際の柱状晶は図 1 に示すような完全周期配列ではないため，方向依存性はないと考えられる[24]．このため，図 5 の x 方向と y 方向の平均的な値が実際の柱状晶に対する無次元透過率と考えることができ，その場合，kc = 9 の KC 式はよい一致を示すといえる．高固相率側で y 方向の流路が完全に閉じ，図 5(b)の高固相率側で結果が大きく KC 式から外れているが，実際には流路が完全に閉じることは考えにくく，kc = 9 の KC 式で表現可能であろうと考えている．また，図 5 と図 7 から，デンドライト二次枝が発達すると KC 式よりも高くなる傾向があり，これをどのように考えるかは，より系統的な評価を行う必要があるが，今回の結果の範囲では kc = 9 と kc = 3 で表現可能と考える．
6．結論
デンドライト樹間液相流動の透過率を効率的に評価可能な数値計算法を開発した．適合
格子分割法(AMR)を複数 GPU 並列化し(parallel-GPU AMR)，さらにブロック分割に対してマザーリーフ法を適用し，デンドライト成長を再現するフェーズフィールド(PF)シミュ
レーションと，樹間液相流動を再現する格子ボルツマン(LB)シミュレーションを高速化した．LB シミュレーションでは PF シミュレーションで用いた格子の倍の大きさを用いることで，更なる計算高速化を可能とした．この高性能計算手法を用いることで，Fe4.5wt%Si 合金の一方向凝固シミュレーションを一次枝間隔を変えて行い，セル構造とデ－ 130 －ンドライト構造を作成し，全固相率を対象とする透過率計算を行った．その結果，柱状晶に対して垂直流れは Kozeny-Carmal (KC)の式において kc = 9，平行流れは kc = 3 とすることでよく無次元透過率を表現できるという結果を得た．これによって，マクロ偏析予測計算において透過率テンソルを高精度に導入できる[25, 34]．なお，本研究では等軸晶の透過率評価も計画していたが，高性能計算法の構築に時間を要し，現時点では達成できていない．今後，等軸晶に対する透過率評価も行いたいと考えている．
謝辞
本研究は「公益財団法人 JFE21 世紀財団」の研究助成のもと実施されました．ここに記して感謝の意を表します．
参考文献[1] C. Beckermann, Modeling of macrosegregation: Past, present and future,
Proceedings of the Merton C. Flemings Symposium on Solidification and Materials
Processing, 2000, pp. 297-310.
[2] M.C. Flemings, Our Understanding of Macrosegregation: Past and Present, ISIJ
Int. 40(9) (2000) 833-841.
[3] C. Beckermann, Modelling of macrosegregation: Applications and future needs,
Int. Mater. Rev. 47(5) (2002) 243-262.
[4] E.J. Pickering, Macrosegregation in Steel Ingots: The Applicability of Modelling
and Characterisation Techniques, ISIJ Int. 53(6) (2013) 935-949.
[5] M. Wu, A. Ludwig, A. Kharicha, Volume-Averaged Modeling of Multiphase Flow
Phenomena during Alloy Solidification, Metals 9(2) (2019) 229.
[6] S. Whitaker, Flow in porous media I: A theoretical derivation of Darcy's law,
Transport in Porous Media 1(1) (1986) 3-25.
[7] T.S. Piwonka, M.C. Flemings, Pore formation in solidification, Trans. Metall. AIME
236 (1966) 1157-1165.
[8] D. Apelian, M.C. Flemings, R. Mehrabian, Specific permeability of partiallysolidified dendritic networks of Al-Si alloys, MT 5(12) (1974) 2533-2537.
[9] D.R. Poirier, Permeability for flow of interdendritic liquid in columnar-dendritic
alloys, MTB 18(1) (1987) 245-255.
[10] A.S. Sangani, C. Yao, Transport processes in random arrays of cylinders. II. Viscous
flow, The Physics of Fluids 31(9) (1988) 2435-2444.
[11] C.Y. Liu, K. Murakami, T. Okamoto, Permeability of dendrite network of cubicalloys, Materials Science and Technology 5(11) (1989) 1148-1152.
[12] S. Ganesan, C.L. Chan, D.R. Poirier, Permeability for flow parallel to primary－ 131 －dendrite arms, Mater. Sci. Eng., A 151(1) (1992) 97-105.
[13] D.R. Poirier, P. Ocansey, Permeability for flow of liquid through equiaxial mushy
zones, Materials Science and Engineering: A 171(1–2) (1993) 231-240.
[14] D. Bernard, Ø. Nielsen, L. Salvo, P. Cloetens, Permeability assessment by 3D
interdendritic flow simulations on microtomography mappings of Al–Cu alloys,Materials Science and Engineering: A 392(1–2) (2005) 112-120.
[15] D. Fuloria, P.D. Lee, D. Bernard, Microtomographic characterization of columnar
Al-Cu dendrites for fluid flow and flow stress determination, Mater. Sci. Eng., A494(1-2) (2008) 3-9.
[16] E. Khajeh, D.M. Maijer, Physical and numerical characterization of the neareutectic permeability of aluminum–copper alloys, Acta Mater. 58(19) (2010) 63346344.
[17] J. Madison, J. Spowart, D. Rowenhorst, L.K. Aagesen, K. Thornton, T.M. Pollock,
Modeling fluid flow in three-dimensional single crystal dendritic structures, Acta
Mater. 58(8) (2010) 2864-2875.
[18] J. Madison, J.E. Spowart, D.J. Rowenhorst, L.K. Aagesen, K. Thornton, T.M.
Pollock, Fluid Flow and Defect Formation in the Three-Dimensional Dendritic
Structure of Nickel-Based Single Crystals, Metall. Mater. Trans. A 43(1) (2011) 369380.
[19] C. Puncreobutr, A.B. Phillion, J.L. Fife, P.D. Lee, Coupling in situ synchrotron Xray tomographic microscopy and numerical simulation to quantify the influence ofintermetallic formation on permeability in aluminium–silicon–copper alloys, Acta
Mater. 64(0) (2014) 316-325.
[20] 棗千修, 高橋大喜, 河嶋佳純, 谷川英司, 大笹憲一, 3 次元数値流体解析による柱状デンドライト組織の透過率評価, 鉄と鋼 99(2) (2013) 117-125.
[21] A. Ludwig, A. Kharicha, C. Hölzl, J. Domitner, M. Wu, T. Pusztai, 3D Lattice
Boltzmann flow simulations through dendritic mushy zones, Eng. Anal. Bound.
Elem. 45(0) (2014) 29-35.
[22] B. Böttger, C. Haberstroh, N. Giesselmann, Cross-Permeability of the Semisolid
Region in Directional Solidification: A Combined Phase-Field and LatticeBoltzmann Simulation Approach, JOM 68(1) (2016) 27-36.
[23] T. Takaki, S. Sakane, M. Ohno, Y. Shibuta, T. Shimokawabe, T. Aoki, Primary arm
array during directional solidification of a single-crystal binary alloy: Large-scale
phase-field study, Acta Mater. 118 (2016) 230-243.
[24] T. Takaki, S. Sakane, M. Ohno, Y. Shibuta, T. Aoki, Permeability prediction for flow
normal to columnar solidification structures by large–scale simulations of phase–－ 132 －field and lattice Boltzmann methods, Acta Mater. 164 (2019) 237-249.
[25] Y. Mitsuyama, T. Takaki, S. Sakane, Y. Shibuta, M. Ohno, Permeability tensor forcolumnar dendritic structures: Phase-field and lattice Boltzmann study, Acta
Mater. 188 (2020) 282-287.
[26] P.C. Carman, Flow of gases through porous media, Academic Press, New York,1956.
[27] M. Ohno, K. Matsuura, Quantitative phase-field modeling for dilute alloy
solidification involving diffusion in the solid, Phys. Rev. E 79(3) (2009) 031603.
[28] S. Sakane, T. Takaki, Phase-field lattice Boltzmann method with two-relaxationtime model for dendrite growth of a binary alloy with melt convection, Comput.
Mater. Sci. 186 (2021) 110070.
[29] E. Khajeh, D.M. Maijer, Permeability evolution during equiaxed dendritic
solidification of Al-4.5 wt%Cu, Modell. Simul. Mater. Sci. Eng. 20(3) (2012) 035004.
[30] S. Sakane, T. Takaki, T. Aoki, Parallel-GPU-accelerated adaptive mesh refinement
for three-dimensional phase-field simulation of dendritic growth duringsolidification of binary alloy, Materials Theory 6(1) (2022) 3.
[31] S. Sakane, T. Aoki, T. Takaki, Parallel-GPU AMR implementation for phase-fieldlattice Boltzmann simulation of a settling dendrite, Comput. Mater. Sci. 211 (2022)
111542.
[32] Y. Hasegawa, T. Aoki, High-performance GPU-kernel Generation Using C++
Template for a Computation for Lattice Boltzmann Method with Octree-based
Adaptive Mesh Refinement, IPSJ Transactions on Advanced Computing Systems
9(2) (2016) 34-45.
[33] R. Suzuki, S.A. Sakane, Takauki, T. Takaki, (2023) submitting.
[34] Y. Mitsuyama, T. Takaki, S. Sakane, Y. Shibuta, M. Ohno, Permeability tensor forvarious columnar dendrite structures, IOP Conf. Ser.: Mater. Sci. Eng. 861 (2020)
012029.

－ 133 －