レビュー
離散的手法に基づく高炉数式モデルの研究開発動向
キーワード:
ironmaking,
blast furnace,
mathematical model,
discrete model,
discrete element method,
particle method
2014 年 100 巻 2 号 p. 198-210
詳細
2014 年 100 巻 2 号 p. 198-210
From the backgrounds of the recent trend towards low reducing agent operation of large blast furnace and application of diversified charging modes of various burdens, an advanced mathematical model of blast furnace is required. Although conventional models based on the continuum model have been widely used, these models are not suitable for the recent demands. The discrete models such as discrete element model (DEM) and particle method are expected to precisely simulate the discontinuous and inhomogeneous phenomena in the recent operation conditions. With the discrete model, the microscopic information on each particle in the packed bed can be obtained besides the overall phenomena in blast furnace. The visual information can be obtained to understand the in-furnace phenomena with high spatial resolution. The liquid dripping and movement of fines in the lower part of blast furnace can be well simulated with DEM and particle method such as Moving Particle Semi-implicit Method (MPS). Moreover, the optimum bed structure for low reducing agent operation is being clarified by application of Eulerian-Lagrangian method. This review summarizes the recent progress on the mathematical model based on the discrete model.
高炉は製鉄所の中核となるプロセスであり，所定の生産量を確保し，製鋼工程に一定品質の鉄源を円滑に供給するために絶えず安定な操業が求められている。一方，高炉は主に固気液の3相から構成され，常温から2000°Cを超える温度領域までも含む不均一な移動層反応器である。安定な操業達成のために，高炉内の複雑な固気液の流動，伝熱，反応などに関わる情報を集め，目的に応じた適切な操業設計，解析を行うことが求められている。しかし，高炉内部を知る計測情報は非常に限定されている。よって，従来から高炉内現象の予測，解析，操業設計のために，高炉数式モデルの研究開発，利用が重視されてきた。また，近年，我が国では5000m3を超える大型高炉が多数出現し，その大型高炉の安定操業を前提にした還元材比低減の達成，原料装入方法の多様化1)，フェロコークスなど新原料の適用検討2)，さらに将来に向けたより抜本的な低炭素化技術の探索などを背景に，高精度で高炉特有の局所的な情報も得られる新しい高炉モデルの研究開発のニーズが高まってきている。
従来，連続体モデルをベースにした高炉モデルの研究開発は活発に行われ，高炉内全体を把握する上で成果を挙げてきた。固気液に加え粉も加えた4流体モデルで一定の完成を見，高炉のマクロ的なプロセス変数の把握はほぼ達成されている3,4)。一方，高炉内は個々の粒子運動から構成される充填層であり，その構造によって気相，液相の運動が決定されるため，モデルの骨格となる高炉内の粒子運動の記述は実現象に最大限，近いことが望ましい。上述のような至近のニーズを考慮すると，個々の粒子情報を直接的に反映できるモデルが期待される。近年，計算機能力の向上もあって，離散的な手法の代表でもある離散要素法DEM(Discrete Element Method)の適用によって，粒子単位の位置情報から高炉内の充填層構造を構築し，炉口部での装入物分布形成から羽口付近の不均一な領域に至るまでの固体流れをビジュアルに表現できるようになった5)。さらに炉下部の液滴の流れについても粒子法の適用によって不規則な運動形態を取る液相移動現象に対応できるモデル解析手法が開発されつつある6)。DEMも含めた高炉内の固気液粉など各相のモデル化に関する基礎的研究については2007年にDong，A.Yuらによって総括されている7)。高炉モデルの離散的手法は発展段階にあり，関連するモデル研究開発は徐々に進歩し，新しい成果が現れつつある8)。本稿では離散的手法をベースに，実高炉への適用，高炉解析を目的に行われてきた至近の高炉モデル研究の動向と今後の展望について述べる。
高炉モデル研究の総括的な流れをFig.1に示す。微分収支に基づく反応工学的モデルは歴史的に1970年頃のYagi and Muchiによる1次元モデルにさかのぼる9)。これは鞭モデルと称されている。以降，2次元，3次元モデル，3次元非定常モデルへと発展してきた10,11)。さらに微粉炭吹き込みの普及と連動して固気液の3相に粉を組み入れた4流体モデルの考案もなされている3,4)。空間次元の拡張，対象とする相の増加によって得られる情報量が豊富になり，多くの成果をもたらした。
Progress of mathematical model of blast furnace.
高炉内現象の原点に戻ると，高炉は鉱石，コークス粒子によって構成される充填層である。装入物は装入後，徐々に降下し，下部では融着帯を経て炉芯を形成し，レースウエイに向かって不均一な運動形態を取り，コークスは消失する。その間の粒子運動，性状変化によって充填構造が定まり，高炉内のガス流れ，液流れはその構造に支配される。よって，緒言で述べたように，充填層の構成粒子の運動の記述が高炉モデルの中で重要な役割を果たす。
Fig.2に高炉モデルに対する近年のニーズを示す。本図に示すように，層状装入に加え，コークス混合装入の採用，炭材内装鉱，フェロコークスなど原料の多様化によって，層の構造設定への自由度，粒子特性の識別への対応力も必要とされている。また，大型高炉前提の還元材比低減指向の中で，高炉下部の炉芯形成，炉芯付近の不均一な固相運動などについて，より詳細な情報が必要とされている。炉下部の液滴下現象についても，その運動は滴の分散，合体などを繰り返す不連続現象であり，集合体としての液相挙動は高炉の生産性，安定操業に関わるため，その特性を表現できるモデルが望ましい。従来のモデルの基本構造は一般的に固体については疑似流体近似，液体，粉体についても固定格子依存の連続体的な定式化によるものであり，高炉固有の3次元的な離散的情報が平滑化されやすい。それらを失わずに現象の数式化ができれば，至近のニーズにも応えることができる。
Change of operating condition and in-furnace phenomena.
粒子運動を離散的に扱うDEM法は非連続固体流れの運動解析などに用いられ，固体粒子全ての運動方程式を解くことで計算領域の粒子挙動を非定常的に追跡できる。高炉内を構成する粒子要素に基づくDEMモデルをFig.3に示す。DEMによる基礎方程式は同図に示すように，個々の粒子間および粒子−壁面間の接触点からの応力とモーメントによる並進運動および回転運動から構成される12)。接触力はバネ，ダッシュポット，摩擦スライダーによって近似したVoigtモデルより算出される。この基本形の拡張により，炉頂からレースウエイに至るまで，粒子固有の力学的挙動に基づいた高炉内全体の固体運動のシミュレーションが可能になる。
Element of discrete element model in blast furnace and fundamental equations.
近年，高炉条件に即した計算手法の進展，計算速度，計算機容量の進歩によって，高炉へのDEMの適用性が拡大している。DEMの構成パラメーターは少なく，これは計算実行上の長所でもあるが，反面，その決定は重要となる。また，高炉適用時には実用的に扱える粒子個数の問題，またDEMでは粒子を真球として扱うために焼結鉱，コークスなど不定形粒子の特性表現に関する課題を解決しなればならない。
DEMの計算速度の改善については，基礎式中の主要パラメーターであるヤング率の値が計算ステップ間隔を通じて関係する。本来，機械的物性値は粒子固有の値に設定すべきであるが，Uedaらはヤング率を0.2~5.0GPaと幅広く変更し，高炉全体のマクロ的解析にはその低減も可能と考え，ヤング率の選択によって計算加速の可能性を示している13)。また，高炉内の粒子総個数は数億個に達し，実用的なDEM計算の範疇を超えている。実高炉対象の計算時には，クラスター化など粗視化近似による個数削減，対称性を考慮した周期境界条件との併用による計算領域の縮小を行い，計算至便性の向上が図られている。不定形粒子の扱いについては，Ademaらは小型スケールの高炉を対象に，Fig.4に示すように，球形の集合体で不定形粒子の形状を表現し，高炉内降下挙動への影響を検討している14)。その結果をFig.5に示す。高炉下部における炉芯高さ，壁近傍の粒子流れに粒子形状の影響が現れることを報告している。Akashiらも同様に，複数粒子の集合体として表現する手法を提示している15)。さらにDEMでは接触摩擦係数，回転摩擦係数の選択によっても，粒子形状の影響考慮は可能である。筆者らは小型冷間実験との対比によって接触摩擦係数，回転摩擦係数の適正値を選択し，近似的に形状因子の影響を表せることを示している16)。計算の目的，計算範囲によって選択肢がある。
Particle shapes and dimensions14).
Pie-slice model simulation results ((a) initial packing, (b) spherical, (c) non-spherical 1, (d) non-spherical 2, (e) non-spherical 3)14).
DEMの高炉への適用の嚆矢としてはKajiwaraらの研究が挙げられる17)。2次元小型スケールモデルを対象に，ベルからの装入物分布形成時の粒子運動解析に適用し，実験との対比を行い，その適用可能性を示した。その後，上述のような計算機能力の進歩もあって，適用範囲，適用条件が拡大し，高炉スケールでの条件，高炉内全体での計算が可能になった。NouchiらはDEM計算による高炉内の粒子間の応力解析に注目し，Fig.6に示すように，高炉内の応力分布は幅広い数値分布を持つ集合として表されることを示している18)。同時に，網目状の応力ネットワーク形成をDEM計算で示している。これらの知見は高炉内の粒子崩壊の解析に有力な情報となる。Hoらはベルレス高炉装入面状での装入時の下層の崩壊について冷間実験との照合を行いながら，DEMの有効性を示している19)。Mioらは粒度に対応した回転摩擦係数の使い分けによって，ホッパー内，ベルレスシュート上の粒径毎の粒子運動と粒度偏析の解析を行っている20)。さらに，その延長としてベルレス高炉の分布形成過程に関する3次元シミュレーションを展開している21)。計算例として高炉の炉口部における装入物分布形成のシミュレーション計算結果をFig.7に示す。Yu and Saxenは炉口付近における装入物降下に伴う上方からのペレット粒子浸透に関して，粒度分布を考慮したDEM計算により解析を実施している22)。冷間実験と対比した計算例をFig.8に示す。炉口部の表層以下の直接観察は困難であるが，DEMモデルは過渡的な分布形成，層内の粒度偏析，粒子浸透などについて3次元性も含めた多くの情報を提供でき，分布制御の改善に役立つ。
Calculated normal stress distribution and stress field by DEM18).
Snapshots of cross section of charged material from chute21).
Comparison of experimental and DEM simulation results for the expansion experiment with a pellet layer on a coke layer (Coke fractions are illustrated in dark gray (small), black (intermediate) and white (large), while pellets are depicted in light gray)22).
DEMは個々の粒子間の力学的相互作用から全体系を表現するために，容器形状などの境界条件の設定，停滞域と流動域の識別などが容易である。よって，空間の異方性に強く，高炉下部の不均一な固体降下運動，非対称的な不均一な固体運動の解析にも柔軟に対応できる。Fig.9，10に小型高炉における装入物分布，並びに羽口半数が閉塞した場合の粒子運動解析を示す5)。Fig.9は標準条件の計算結果である。層状装入による炉内の充填構造が再現でき，特定粒子の追跡により流線，タイムラインが把握できる。羽口閉塞は高炉トラブル時に頻発する事象であり，上方への影響解析は早期立ち上げに重要である。閉塞羽口を設定したFig.10の計算結果では，閉塞羽口直上は停滞部となるが，シャフト中段では徐々に閉塞の影響は薄れ，シャフト上部では閉塞部と反対側の降下速度が速まる結果が得られている。これらは，従来の冷間実験と良く一致する23)。また，DEMは非定常計算であり，不規則な粒子移動にスリップのような荷下がり不安定時のスリップの頻発部位の提示など，不連続な現象についての情報も得られる。Fig.10c)はせん断応力と法線応力の比をスリップ指数Sliとして表したものである5)。レースウエイ奥，炉腹部壁近傍など粒子速度，方向の変化が大きい個所に，スリップ指数Sliの高い降下不安定領域が現れやすいことなどを推測できる。
Descending behavior of burden in blast furnace ((a) Shape of layers, (b) Calculated stream line).
Descending behavior of burden in blast furnace where half of tuyeres are inactive ((a) Shape of layers, (b) Calculated stream lines, (c) Slipping factor).
DEMで得られた粒子配置，空隙率分布などの粒子情報とCFD計算(Computational Fluid Dynamics)によるガス流れ解析の連結によって，充填層内の固気運動の同時解析が可能になる。DEMの基礎運動方程式にガスと粒子の相互作用力を入れ，気相の運動方程式と連立させればよい。Nakano and Yamaokaは小規模スケールの条件でレースウエイ近傍の粒子運動の安定性解析を行っている24)。Zhouらは試験高炉スケールで炉下部の粒子運動とガス流の同時解析を行っている(Fig.11)25)。実高炉スケールでは，Umekageらは高炉下部を対象に9°の扇形の計算領域でコークス運動とガス流れ計算を行い壁付着の影響を解析している26)。この計算例をFig.12に示す。さらに，これに続く研究では1600万個の粒子を含む90°扇形の高炉下部を対象に，粒子運動とガス流れの精緻な解析を行っている27)。
Snapshots of solid flow and corresponding gas flow in the cohesive zone25).
Calculated instantaneous coke particle velocity vector on the actual blast furnace for 7 degree downward tuyere direction (a) without scaffoldings, b) with scaffoldings)26).
高炉全体を対象にした計算の場合，鉱石，コークスの層構造などDEM計算情報(Lagrangian method)をCFD(Eulerian method)の計算メッシュに高精度に転写しなければならない。すなわち，層を構成する各粒子の位置，空隙率などDEM計算による粒子単位の離散的な情報をCFDの格子情報に高精度に置き換えてガス流れ計算を行い，再びガスの抗力を考慮したDEM計算を繰り返すことになる。この計算手順の概念をFig.13に示す28)。CFD計算で扱う運動方程式の粒子情報は実粒子径に基づくものであり，DEM計算で粗視化近似を施している場合，その扱いには工夫を要する。Natsuiらは注目粒子を中心に周囲粒子との立体的な幾何学的配置を考慮した検査空間を設け，個々の粒子周辺の空隙率をセル−粒子中心間距離で重み付けし，3次元的な粒子情報を失わずにCFD計算格子に写す方法を提案している28)。
Combination of discrete element model of solid movement and continuum model of gas flow.
Fig.14に実高炉対象の3次元DEM-CFD計算の実例として高炉炉口部における鉱石装入時のガス流れ計算結果を示す28)。計算領域は半裁である。鉱石装入時に過渡的に炉口中心で流動化が生じることは古くから実験的に確認されているが29)，非定常的な粒子運動とガス流れ計算を連結した本モデル計算によって，この現象の再現ができる28)。Fig.15に5000m3クラスの高炉内全体にわたるガス流存在下での3次元的な装入物分布，ガス流速分布，粒子へのガス抗力分布を示す。Fig.16は粒子移動速度，法線応力分布を示す。粒子へのガス抗力分布からレースウエイ奥と炉頂中心付近に粒子に大きなガス抗力が働いている状況が観察される。粒子移動速度から炉芯に相当する降下遅滞部が現れ，応力分布からはその炉芯領域に大きな荷重が集中することなどがわかる。これらの結果は高炉の安定制御に関わる重要な情報となる。
Transient gas velocity change at the throat during ore particle charging.
Burden distribution, gas velocity and gas drag force calculated by DEM-CFD model.
Particle velocity and normal stress calculated by DEM-CFD model.
この応用として，今後の低炭素化の一手段として注目されるシャフトガス吹き込みの3次元のガス流れ解析も行われている30,31)。大型高炉，小型高炉におけるシャフトガス吹き込み時の3次元半裁のガス流れ計算結果をFig.17に示す。シャフ下部の補助羽口から吹き込まれたガス流の影響範囲が認識でき，高炉内容積の大小によってその影響領域が異なることがわかる。シャフトガス吹き込みは古い概念であり，均一充填層の冷間実験が既に行われているが32)，層内部のガス流の観察は困難である。本シミュレーションでは層状の充填構造，3次元的な補助羽口の設置を考慮し，ガスの浸透の効果が予測できる。その影響範囲は吹き込みガスの相対的量比に依存するなど定量的な結果も得られている30)。
Effect of shaft gas injection on gas flow in blast furnace.
高炉の融着帯は通気性評価で重要な領域である。融着帯において鉱石粒子は還元の進行と共に，荷重の影響を受けながら軟化，収縮が進行する。その鉱石粒子の軟化，収縮など性状変化のモデル化は困難を伴う。DEMでは粒子運動の精緻な移動表現が可能であるが，真球を対象としており，軟化などの性状変化への対応はそのままでは難しい。筆者らは荷重軟化試験における通気性変化の結果を参考に，融着帯部においてDEM計算上の粒子ヤング率を変更し，その領域における粒子のめり込みによって融着帯の空隙率低下を表す方法を考案している33)。その概念をFig.18に示す。
Modeling of softening and shrinkage of ore particles by DEM.
高炉内に上記の鉱石軟化モデルを取り入れることにより，融着帯の特性が評価できる。Fig.19は鉱石軟化を取り入れた標準的な条件の高炉内の空隙率分布，ガス流れベクトルを示す33)。本計算では塊状帯では鉱石，コークスのヤング率を1.0GPa，融着帯では鉱石のヤング率を0.02GPaと変更している。設定された融着帯領域において，融着帯特有のコークススリット流などガス流れの可視化が可能となる。その応用例として，Fig.20に二つのコークス比条件，350kg/t，240kg/tにおける融着帯付近におけるガス流れを示す33)。低コークス比の240kg/tの場合，コークススリットの薄層化によってスリット部の通気抵抗が上昇するため，融着帯上部を避け融着帯根部のスリットを通過して壁に衝突するガス流が観察される。Fig.21は等圧面分布を示す。融着帯部の通気抵抗が大きく上昇している状況がわかり，低コークス比時における層状装入の課題も伺える。
Void fraction distribution and gas velocity vector considering cohesive zone.
Burden layer and gas velocity vector for each coke rate condition (red particle: ore, gray particle: coke).
Isobar planes in blast furnace for each coke rate condition.
なお近年，連続体モデルにおいても，異相界面の表現が可能なVOF(Volume of Fluid)の適用によって，層構造を表現できるモデルがNishiokaらによって開発されている34)。この計算例をFig.22に示す。固相運動に関しては疑似流体近似であるが，本モデルは反応も考慮した総合モデルである。DEMモデルによる融着帯部のコークススリット流と類似した結果も現れており，対比として注目される。
In-furnace states by BF model34).
融着帯の通気性改善を目的に鉱石層へのコークスの混合装入が一部の高炉で既に実施されている1)。将来的な低コークス比操業時には，融着帯のコークススリットは薄層化するため，コークス混合装入はその部位の通気性改善のために有効と推測される。荷重軟化試験装置での効果確認は既に行われているが1)，DEMでは，層状装入，混合装入を問わず分布構造の変更が容易であり，実高炉対象の定量的評価が可能である。著者らは融着帯モデルの発展型として，混合装入時の通気性改善効果を解析している35)。
混合装入の基本的な評価例として，荷重軟化試験装置を対象にDEM計算による収縮後の層構造変化，粒子間応力分布をFig.23に示す35)。コークス混合比率の増大によって層の収縮は緩和され，応力分布から主にコークスが荷重を受ける状況が観察される。Fig.24はヤング率と空隙率，コークス混合比率と相対圧損の関係の関係を示す。コークス混合によって鉱石の軟化収縮に伴う空隙率低下は改善され，体積比で50%のコークス混合は鉱石単味の場合より約20%に圧損を軽減できるなど，コークス混合装入の効果が定量的に把握できる。Fig.25は，240kg/tの低コークス比条件において，層状装入と鉱石とコークスを完全に混合した場合の等圧面分布と融着帯付近のガス流れ変化を示す。混合装入によって融着帯の圧損は減少し，ガス流は融着帯の鉱石層を貫通し，斜め上方に向かうガス流が観察される。このようにDEM計算の充填構造設定への対応力を活かし，混合装入による融着帯の通気性改善効果など，今後の低炭素操業時に最適な層構造設計に有用なデータを得ることができる。
Simulation of shrinkage behavior of ore particles under coke mixed conditions (Vessel diameter:0.12m, Height:0.1m, Particle diameter 0.01m).
Influence of coke mixing rate on void fraction and relative pressure drop in packed bed.
Isobar planes in blast furnace for layered charging and coke mixed charging.
層内の固相反応は粒子単位の反応から構成される。従来は計算格子内で加成的に扱われることが多いが，本来，層内の伝熱と反応は空間的に非線形な事象である。伝熱，反応は互いに連結した事象であり，粒子の局所的な配置によって，その挙動が変化することが予測される。近年，鉱石，炭材の近接配置によって還元機構を制御する考えもある36)。固定格子での計算では，格子間で粒子情報が平均化され，局所的な粒子配置の影響解析には対応しにくい。DEMによって粒子を個別に識別し，注目粒子毎に反応，伝熱を組み込み，CFDによるガス側移動現象を連立させれば，粒子単位で高炉内現象を精緻に捉えることができる。著者らはこのEuler-Lagrange的手法に注目し，第一段階として鉱石，コークス粒子から成る約150m3の固定層を対象に3次元の熱物質移動の直接解析を行っている37)。その計算フローをFig.26に示す。本モデルは伝熱に関して粒子毎の対流，輻射，接触伝熱を考慮し，粒子単位の反応識別により，個々の反応熱の影響を反映できる。
Calculation flow by Eulerian-Lagrangian method.
Fig.27は層状装入(CaseA)，混合装入(CaseB)の二つの条件における充填層内の一定時間後の粒子単位の還元率分布を示す。共に壁流の影響が観察される。層内の平均還元率には両者に差がない。拡大図をFig.28に示す。層状装入ではコークス層への鉱石粒子の浸透が見られ，鉱石層内では下部の鉱石還元率が上方粒子よりわずかに小さい。その下層にあるコークス層のソルーションロス反応による低温化が上方に影響している。このように，反応，伝熱が相互に影響し，充填層内は非常に不均一な状態で反応が進行し，層構造設計がそれに大きく関わるなどの新しい知見が得られつつある。
Distribution of reduction degree in packed bed for layered and mixed structures.
Enlarged view of reduction degree of ore particles in packed bed.
さらに充填構造を任意に設定できる利点を生かし，高反応性コークスの混合装入時の粒子単位の熱移動，反応率解析の計算例も報告している38)。Fig.29は高反応性，あるいは通常コークスの小塊コークスから成る混合層を鉱石層の下部あるいは上部に配置した計算結果である。全体に高反応性コークス使用時には低温化によって鉱石の還元進行は遅くなる。また，鉱石層下部に高反応性コークスを配置すると，コークスガス化の吸熱反応によって上層の鉱石層温度変化を経由し，鉱石還元の進行もその影響を受けるなど，微視的な伝熱，反応の相互作用を捉えることができる。本手法によって高反応性コークスの最適配置，鉱石の被還元性との関係など，新しい機能を持つ装入物使用時の操業諸元決定に活用できる。
Enlarged view of reduction degree of ore particles in packed bed including high reactivity coke and conventional coke.
高炉下部では融着層で発生したスラグ，溶銑がコークス充填層内を滴下し，炉床部に達する。その間の液相の流れは高炉の生産性，操業安定性，さらに溶銑成分に影響するために重要な現象と認識されている。炉下部の液流れに関してはレースウエイ近傍のプローブによる直接観察などを除いて，観察，測定手段がほとんどないために，数学モデルによる理解が重要である。
炉下部の液の流れのモデル化に関して，従来ではDarcyの式を適用した連続体による表現が行われてきた39)。その後，液流れを独立した粒子運動として液分配を確率論的に捉え，液流れを定量化する手法が開発されている40,41)。また，液の分配モデルと連続体モデルをカップリングさせた計算モデルも研究されている42)。液流れ全体については，融着帯直下では融液は「つらら」状に滴下し43)，以降，冷間実験などの観察から，液滴の分裂集合，合体，液滴の変形など複雑な現象から液流れ現象は構成されると考えられる44,45)。さらに，その液流れ挙動には液相とコークス充填層との濡れも影響する。
6・2 粒子法による液流れのシミュレーション以上の観点から，液滴の分散，合体，固相との相互作用を合理的に表現する手法として，Lagrange的手法である粒子法によるシミュレーションが注目されている。粒子法を用いた場合，格子設定が不要であり，界面の複雑な計算をすることなく，自由表面流れの解析が可能である。その一つであるMPS法(Moving Particle Semi-implicit)は46)，非圧縮性の自由表面流れを解析するために開発された手法であり，液相を有限個の疑似的な粒子の集合体で表現し，分散相としての自由表面を含む液相運動を比較的容易に記述することができる。MPS法を用いた表面張力および濡れ性モデルもすでに多く提案されており，近年では粒子間ポテンシャルを用いた表面張力・濡れ性2Dモデルが提案されている47)。このモデルは，分子動力学(MD)法と同様のアルゴリズムを採用して，粒子間に界面張力を与えるため，アルゴリズムが明快であるという長所がある。
NishiokaらはMPS法による高炉融着帯以下の領域における巨視的な液滴下シミュレーションを行っている48)。Konらは特定部位の充填層領域に注目し，MPS法によって，疑似的な粒子から構成される充填粒子の設定によって固相との濡れ，相互作用を考慮した液流れモデルを提案している。そのモデルの構造と主要な基礎式をFig.30に示す6)。充填層内で単一液滴が分散，滴下する挙動に関する計算結果をFig.31に示す6,49)。本図は動粘度をスラグ相当から溶銑まで大きく変えた場合の結果であり，液物性によって液相の動的挙動が大きく変わることがわかる。この手法は液物性の設定に自由度が高く，動的な界面変形を捉えることができ，融着帯直下などの局所的な液滴下の解析には適している。今後の課題は計算領域の拡大と気相との相互作用の組み込みである。
Concept of MPS method and fundamental equations for modeling.
Effect of viscosity on liquid flow in packed bed calculated by MPS method.
高炉下部では固気液相に加えて，粉体の蓄積，運動挙動の解析が高炉内の通気性確保の上で重要である。特に微粉炭吹き込み時には未燃チャーがレースウエイで多量に発生し，通気性に影響する。それ故，粉体の運動モデルが1990年代から提唱され50,51)，マクロ的な粒子群の運動解析，空隙率変化などについては4流体モデルで完成を見ている3,4)。さらに，高炉下部の局所的な粒子の蓄積，閉塞など不連続現象についての知見が得られれば，粉体運動に関係する全貌が明らかになる。
筆者らは充填層内の粉体の微視的な挙動に注目し，Fig.32に示すようなDEMによる充填層と粉体，粉体と粉体間の相互作用から成る運動モデル開発を行っている52,53)。Fig.33はガス流に伴う下方からの抗力によって上方に運動する粒子群の計算結果である。充填粒子の下面に粉体が付着し，その間隙では上方に粉体が飛散する状況が観察できる。これらの結果は冷間実験と一致する54)。Fig.34は充填粒子の間隙に閉塞しつある粉体にかかる応力を示したものである52)。狭い間隙に粉体が集中すると，下方からの粉体が付着し，粉体間に荷重がかかりブリッジが生成されるなど閉塞のメカニズムが把握できる。今後，充填層構造の変化した場合，また羽口など横方向から上方に気流の運動方向が変化する場合への解析適用が期待される。
Element model of fine particle movement by DEM.
Fine particles movements in packed bed for various gas velocity conditions.
Snapshots of calculated normal stress distribution for vertical cross section.
従来，高炉の計測情報が不足し，炉内そのものがブラックボックスに近いことから，まず炉内全体を把握するためにマクロ的な数式モデルの研究に主眼がおかれ，操業設計，解析にその役目を果たしてきた。しかし，冒頭に述べた今後の高炉モデルへのニーズを考えると，高炉全体の把握から局所的なものまで多岐に，また空間解像度もより高いレベルのものが要求されている。離散的な手法は素現象から組み立てるものであり，従来手法と異なり，粒子単位の情報に基づく装入物分布形成，降下挙動，粒子運動に連動したガス流れ，融着帯の性状変化，液相および粉体の運動などのモデル的記述が原理的に可能である。徐々に適用範囲も拡張されつつある。上記のニーズに原理的に対応できる。離散的手法の特徴を生かし，粒子単位のミクロスケールからマクロ規模まで，個々の粒子情報に基づいた高炉固有の特性が明らかになってきた。その全体像をFig.35に示す。高炉内の各種の挙動について，ビジュアルな情報と数値情報が同時に得られ，今まで見えていなかった高炉の姿が把握できるようになった。ただし，その全体系への拡張に対しては発展途上であり，反応と伝熱を加えたモデルは実高炉スケールでの流動，反応，伝熱などを統合した総合モデルにまで至っていない。総合モデルに関して，現時点では連続体モデルが優位性を持つが，複雑な層構造，混合装入などへの対応力には課題も有すると思われ，離散的手法の展開が期待される。
Various mathematical models based on discrete method for blast furnace.
Table 1には対象スケールと計算モデル構造の関係を示す。各モデルには空間解像度(Spatial resolution)，計算負荷(Computational workability)など，それぞれに差異，特徴がある。解像度と計算負荷は相反する。最終目的は統合モデルの構築であるが，当面は使用目的に応じ，必要とする空間解像度，計算負荷の条件を考慮してこれらのモデルを階層的に併用していくことが現実的であると思われる。
Characteristics of various models.
近年，計算機性能は着実に進歩しつつあり，今後，要素モデルの統合化に向けては，不連続的な粒子情報を保存して解像度も優れ，同時に連続体的手法の長所も保有するEuler-Lagrange的手法の更なる展開によって，徐々に目標に近づけると思われる。同時に，数値計算技術だけでなく，局所的な情報と全体とのバランスに基づいた合理的なモデルの考案など，現象の数理的記述方法にも工夫が必要とされる。それら両面からの取り組みによって，反応特性を含む粒子特性，充填層構造から組み上げた高炉の総合的解析，高炉の制御因子に基づく様々な思考実験が可能となり，机上でのVirtual高炉の実現も遠くないと思われる。
高炉モデルに関する至近のニーズと，それに対応しうるモデルとして離散的手法によるモデル開発動向を総括した。以下，内容を記す。
高炉内の構造を決める粒子運動の記述が高炉モデルで非常に重要であるが，離散要素法は高炉内の固体運動を精緻に記述することができ，高炉全体の固体運動を表現できるところまで至っている。その固体運動をガス流れモデルを連成させることにより3次元的な高炉内のガス流れだけでなく，粒子への荷重，荷下がりの安定性など様々な情報が入手できるようになった。溶融帯の鉱石軟化にも適用できる。高炉下部の液滴の運動，粉体の運動解析にも粒子法など離散的手法の応用によって，新しい解析手法として研究が進みつつある。最後に各モデルの特徴比較と今後の展望について述べた。
