・第6回 もう一度基本に戻る-ガスで押される液体

 第5回で,実際のノズル周りの流れを計算しようとして躓いてしまった.そこで,もう一度,基本に立ち返って,ガスで押される液体の流れを定量的に評価してみよう.


・基本ケースの作成(case00)

 空間的メッシュ密度を増やすために,カーフ近傍のみをモデル化した新しいケースを作った.当面,下ガスはなし.

解析領域 0.005×0.002×0.004[m]
メッシュ 40×40×40,等間隔
時間刻み 1μs
Light fluid(気体) 空気(@20℃)
Heavy fluid(液体) 水(@20℃)
境界条件 手前と上を除く全面流出境界
重力 -z,9.8m/s^2
スラブ(白) 断熱個体
液体 初期条件として与える.0.3×0.15×1.0[mm^3]
上ガス流入量 -20.0[m/s]

q1ファイル


手前側,白く浮いているのが液体.上ガスは直径2mmの円盤から下向きにガスが出てくるという条件で与えた.メッシュは,まだテストなのでこんなもん.


 計算結果は上の様になった.第5回のように,物理量が定義されない謎の空間も生じない.液体は,上側はガスに押されて移動するが,下側はそのままとどまる.意外だったが,物理的にあり得ない動きではない.1msの間に液体上面が動いた距離は約0.5mm.一方,ガスの圧力(200Pa),面積から計算すると上面に掛かる力は10^-5Nで,加速度は液全体を一つの剛体とすると200m/s^2である.1msの間に動く距離は0.1mmとなる.上側が流される現象を加味すれば,観察された動きは理にかなっている.

・ガスを窒素に変える(case01)

 続いて,light fluidを27℃の窒素に変更し,同じ計算を繰り返した.ガスの温度が7K,平均分子量が1程度違うだけで,観察される現象に違いは無いものと思われた.が,下のように,またしても奇妙な空白が観測された.
q1ファイル

これで,何かがおかしいと気がついた.そこで,Case00とCase01のq1ファイルを比較してみた.すると,驚愕の事実を発見!!



上は,Case00とCase01のq1ファイルをファイル比較ソフト「DF」で見たものをキャプチャしたものである.黄色がCase00,右がCase01.SEMモデルは,流体は基本的にHeavy fluidで解き,SURNで液体と気体の境界を決定する,単相流モデルである.ところが,Case01は,流体の密度が窒素ガスの密度に置き換わっている・・・

・Case01の誤りを修正(case02)

 Case01で流体密度が窒素ガスのものに切り替わってしまったのが,Case01の奇妙な空白の原因と考えた.どうやら,SEMモデルでガス種を変更するときに,RHO1やENULなどの変更を抑制する処理を忘れているようだ.つまり,VR Editorのバグ.

 おそらく,このバグは,流体の種類を変更すると無条件に物性値が変わってしまうというものと推測されるので,もう一度「heavy fluidを水に再定義する」手順を踏めば元に戻るだろう.



・・・予想通り.
q1ファイル


 ところが,計算結果を見て二度びっくり.Case01と全く変わっていない.つまり,SEMモデルにおいては,q1ファイルのRHO1やENULなどの情報は使われていないらしい.

マニュアルを検討したところ,SEMモデルでは,heavier fluidの物性,lighter fluidの物性はGroup 19の

IPRPSA = 67 ;IPRPSB = 17

で定義され,流体物性テーブルPROPSから引いてくるのだという.

Phoenicsは,いろいろな流れが計算できるので重宝だが,マニュアルの整備がとっても悪い.たとえばSEMモデルの記述においても,上で問題になったphase 1の物性値はマニュアルでは

RHO1 = GRND10
ENUL = GRND10

でなくてはならないはずなのに,なってないし.ちなみに,手で上のように書き換えても,VR Editorが勝手に元に戻してしまう.というわけで,唯一頼れるドキュメントであるはずのマニュアルも,問題解決の助けにはならなかった.

 

・Case02の時間分割を細かく(case03)

 物理量が定義されない空間が出来てしまう原因は,時間か空間の刻みが粗すぎるからではないかと推測した.空気だとOKで,窒素だとNGというのも考えにくいのだが.ともかく,原因を確かめるために,Δtを1μsから0.5μsに変更した.

解析条件 Case02と同じ(N2+水)
時間刻み 0.5μs

q1ファイル

結果は,Case02と全く同じであった.つまり,時間刻みは充分細かいと言うこと.

 

・その他の試み

 その後も,様々な試みを行った.全部失敗したのでq1ファイルと結果はつけないが,


という試みが全て失敗.

・Eureka!

 あれから,何日がたち,何回の試計算を行っただろうか.とうとう問題解決の糸口を発見した.SEM法には唯一,ユーザーが設定可能なパラメータがある.

これを変更してみたところ,驚くべき効果が!
 変数RLOLIMとRUPLIMはデフォルトではそれぞれ0.4,0.6に設定されている.SURNの値が0.4より小さいセルは全て液体,0.6より大きいセルは全て気体で満たされていると判断するわけだ.上の数値はPhoenicsの推奨値だが,これを変えるとどのような効果が(あるいは不利益が)生じるかはマニュアルには記載されていない.参考文献リストも

Jun L, Spalding DB 1988 "Numerical simulation of flows with moving interfaces." PHOENICS Journal of Computational Fluid Dynamics Vol 10, No 5/6, pp 625-637 1988. Published by Pergamon Press

Hamill IS, Jun L, Waterson N 1991 "A model for the simulation of three-dimensional mould-filling processes with complex geometries." Presented at the International Conference on Mathematical Modelling of Materials Processing, Bristol, 23-25 September 1991

Spalding DB, Jun L 1988 "Numerical simulation of flows with moving interfaces." Published in PCH PhysioChemical Hydrodynamics Vol 10, No. 5/6, pp 625-637 1988

と,入手が難しそうなものばかり.今のところは,仕方がないのでその理由もわからずにオプションを変更することにする.方針としては,Lower limitとHigher limitの間のwindowを小さくしてゆくことだ.どんな不利益が起こるかわからないので,すこしずつ慎重に行うことにした.その結果,

RLOLIM=0.495
RUPLIM=.505

でも計算が無事収束することを確認.計算結果はCase00とほとんど変わることはなかった.しかも,嬉しいことに,空白のセルが全く生じなくなった.



・ところが,次の障害が!

 上の計算結果をよく見てみると,液体の質量が徐々に減っていることがわかる.右上のAverage valueはSURNの平均値だから,水の質量が保存されていれば一定の値になるはず.

 確かに,液体は流れ出すのではなく,少しずつ削られているように見える・・・しばし悩んだ後,その原因がわかった.質量保存とあわせて考えると何のことはない.水は,徐々に系から流れ出していたのである.ただし,霧となって.

 下の図は,上の計算結果の表示範囲を変えてプリントアウトしたものである.SURN~0の領域がバックグラウンドに埋もれないようにした.すると,驚くべきことに,SURN<0.1の領域が,あたかも煙のように液体から離れて,ガスに乗って下流へと飛ばされているのが見える.つまり,ガスに押された液滴が,砂のようにばらばらになって吹き飛んでいるという計算結果が得られていたのである.

 この結果は,直感的に予想される結果と大きく異なるが,「SEMモデルは表面張力を考慮しない」ことを思い出すと,理由は明らかである.液体で定義されたセルに自らを一つの固まりに保とうとする力が働かないため,外力が働くと外力が働いた部分だけが吹き飛ばされるのである.これでは,レーザー加工の物理を研究するためには不都合.

SEMモデルに表面張力を導入しよう!

・まとめ