以下の内容は、まだ未完成です。
(Click here for explanation in English.)
今回の円管内流速分布計算のソフトウェアは,第3.1節および第3.2節で説明した最初のソフトウェアを発展させたもので,これの詳しい説明は,第3.3節で行っています.本節を読む前に,それらの節の内容を十分理解している必要があります.今回提供するソフトウェアでは,円管の直径とほぼ同じ長さの管路長で計算しますが,表面粗さを模擬する鋸刃の数は整数なのでRe数により長さの違いが少し出ます.計算途中で中間結果を出力しながら定常に達した後も,滑らかなデータを得るために計算を続けますが,定常判定は簡単ではないので,一応,510ステップまで計算を続けて終了するようになっています.ここで,1ステップとは,本DSMC計算の分子移動と分子間衝突を分離している最小時間間隔Δtm(DTM=1.1μs)の200倍(220μs)を意味し,途中の1次結果ファイルの出力間隔でもあります.計算に慣れてくれば,独自判断で定常の達成を判定して,途中終了をすることは可能ですが,どちらにしても計算終了までにはかなりの時間が必要です.一例を挙げれば,一つのRe数の計算完了までに要する時間は,CPUにRyzen9_9950X3D(9950Xでも9950X3D2でも計算時間はほぼ同じ)を使用すると約44時間(32 thread使用),Ryzen7_6800H を用いるとその倍の約88時間(16 thread使用)となる(ただし,共に,M.2のSSDを使用した場合).
入力するデータは,① 計算に使用するスレッド数(最大数を推奨),② 最初からの計算か継続計算か,③ Re数 (Re数の入力は,1000から7000に制限される.その理由は,1000以下では統計的変動が激しくなるし,7000以上では圧縮性が出てくる可能性が考えられるから),④ Bird法か U-sys法かの選択,⑤ 鋸刃で近似する表面粗さの3つのパラメータを,Re数から自動決定するか,あるいは独自に値を定めるかの選択(独自に定める場合は続けて3つのデータを入力).⑥ 表面粗さなし(拡散反射面)か,表面粗さ有りかの選択,⑦ 初期速度分布(一様分布か放物線分布)の選択,⑧ 結果を gnuplot用にするか tecplot用にするかの選択.それぞれ,事前に決めておいて入力する.なお,表面粗さを考慮しない壁面(滑らかな面)の場合は,拡散反射面(乱反射面)で計算することになる.
平均流速は Re数に比例しておおよそ10 m/sから70 m/s まで変化するが,計算の基礎となる分子速度は約350m/sであるため,平均流速が低いほど,結果に現れる統計的変動は大きくなる.したがって,滑らかな平均値を求めるために必要なサンプル数,すなわち計算の繰り返し数が多く必要となる(具体的な方法は後述).Re数を変更する方法としては,密度や代表長を変更することも考えられるが,その場合,DSMC計算において物理空間を分割するセルの細かさが計算結果にどのような影響を与えるかを前もって詳しく調べる必要がある.そこで今回は平均流速を変えることだけでRe数の調整を行っている.
10個の出力データの並びは,x 座標,y 座標(r座標),数密度 n,流速 u,u/U(流速u/平均流速U),流速 v,温度 T,x方向温度 Tx,yz方向温度 Tyz,マッハ数 M である.長さは円管の直径(代表長さ)によって無次元化され,速度は分子の最大確率速さあるいは平均流速によって無次元化されていることに注意してください.数密度と温度はそれぞれの初期値で無次元化される.なお,これらの出力データは,第3.2節での計算とは異なっているので注意してください.
●● 計算の実際 ●●
以下に,円管内流速分布を求める計算の手順を述べるが,計算の前提として,第2章「試し計算」や第3章・第2節 (3.2節) に書かれた内容を熟知し,実際にその計算を体験しているものとして話を進めるので,もしそうでない場合は,先にそれらの計算を実行して,問題なく結果が得られるかどうかを確かめてください.
(A) 計算に必要なファイル
計算に必要な4つのファイルは,①計算ソフトウェア本体(AVT5F.exe),②65個の乱数列の初期値ファイル(ransuuM65),③45ステップ毎の結果ファイルを複数個平均する場合のソフトウェア(SUM009.exe),④その他のファイル1つ(libiomp5md.dll).それらは,1つの圧縮ファイル Tube_Turbu (TTurbu) に保存してあるので,以前に述べた「2. 試し計算」等の場合と同様に,それをダウンロードしたのち,解凍作業を行って4つのファイルに戻します.ダウンロードフォルダ内あるいはデスクトップ上に TTurbu というフォルダができ,その中に4つのファイルが入っているので AVT5F.exe を実行する.実行しようとした際,もし,Microsoft Defender SmartScreen が「SmartScreenを使用できません」と言ってきても,SmartScreenは賢くないので,気にせず,「実行」してください.なお,計算実行が行われるコマンドプロンプト画面は,そのウィンドウ幅を170以上にすると見やすい表示になる.
Tube Turbu (TTurbu)
(B) 計算の開始(キーボードから複数のデータを入力)
①計算ソフトウェア本体(AVT5F.exe)を立ち上げると,使用しているPCのスレッド数が表示されるので,その数値以下の数字を,実際に使用するスレッド数として入力します.1から最大値まで以外の数値が入力されると,最大値が選択されます(推奨値は最大値です).
②次に,初めての計算の場合は 0 を,継続計算の場合は 1 を入力します.継続計算の場合は,uvw333やtfile2 等のファイルがすでに作られている必要があります.継続計算の場合,1を入力するとすぐに計算が開始されるが,初めてする計算の場合は,さらに以下の入力が必要になる.
③Re数(1000~7000)を入力します.範囲以外の数値が入力されると,再度,入力要求がされる.
④Bird法で計算するかU-system(U-sys法)で計算するかを選択します.基本的には,U-sys法で計算するものであり,Bird法は比較の計算をしたい人のために付けてあります.なお,本計算でBird法と呼ぶものは,実際には,修正Bird法を指しています.
⑤表面粗さを鋸刃形状で近似するための3つのパラメータを,Re数から自動決定するか,あるいは独自に値を定めるかの選択をする (独自に定める場合は続けて3つのデータを入力)
⑥円筒の内壁面を,「滑らかな面」にするか「粗い面」にするかを指定します.滑らかな面を選択した場合は,直前の入力⑤は無視される.なお,滑らかな面では,Re数に関わらず層流を示す放物線分布の流速が得られるはずです.
⑦計算を開始するにあたり,流速の初期分布の選択をします.平均流速で一様な分布にするか,放物線分布にするかを選択します.なるべく,最終結果とは異なる初期分布を選んだ方が,速度分布が徐々に変化していく様子が観察できます.
⑧結果の出力形式を,gnuplot形式にするか,tecplot形式にするかを選択します.なお,gnuplot用とtecplot用の2つの結果ファイルの違いは,tecplot用では先頭の2行が余計についているだけなので,手作業でそれを削除すれば gnuplot用に変更できる.2つの結果ファイルを見比べて,変更する方法は各自で見つけてください.
以上の入力が全て終わると,計算が開始されますが,まず最初に,次の表示がなされます(下記の数字1は,0かも知れません).
B-or- U, Wall (S or R), Flat-Palabo, GNU_or_TEC
1 1 1 1
数字は順に,今回の計算が,Bird法かU_sys法か,滑らかな面か粗い面か,初期速度が一様分布か放物線分布か,結果の出力形式はgnuplotかtecplotかを示しており,同じ内容が,IN-DATAというファイルの先頭にも書かれています.また,ID50-A40-DTH0-DX2というファイルには,1番目の数値に「その逆数が鋸刃形状のゆらぎ幅」が,2番目の数値に「鋸刃の角度」が,4番目の数値に「鋸刃間隔を決めるための係数」が書かれている.
(C)計算実行時の様子
計算実行中に注意することは,その時点で計算がどこまで進んでいるかや,結果のファイルが順調に出力されているかということです.まず,画面表示の中で,次の表示部分に注目してください.
——————————————————-
@@ Re = ….. ….. ….. …..
——————————————————-
最初のRe数ですが,必ずしも入力したReの値と完全には一致していないかも知れません.これは,分子の速度を乱数で発生させているので,多少の違いが出てくるものです.2つ目からの数値は,まず,次の説明を理解してください.すなわち,
● DSMC計算の最小時間単位である「分子移動と分子間衝突を分離している微小時間」が 200回(これを本計算の1ステップと呼ぶ)繰返される毎に,サブフォルダDATxの中に1次結果ファイルsjetc001, sjetc002, ….. が作られていき,それが15個貯まるとDATx内に2次結果ファイルAAA-001-015が作成され,それまでのsjetcファイルは削除される(上記のファイル番号は順に大きくなる).さらに,AAA-xxx-yyy ファイルが3個できる(45ステップ進む)と,自身のフォルダにOUT-001-045(3次結果ファイル)が作成されて,同時に,フォルダDATx内にある一番古いAAA-***-*** ファイルが削除される.作成された3次結果ファイルは,gnuplotかtecplotで読み取れる形式になっている.15ステップ進むごとに1個のOUT-ppp-qqq ファイル(45ステップ分の平均)が作られるが,連続するOUT-***-*** のファイルは部分的に同じデータを使った平均であることに注意してほしい.
ディスプレイ表示のRe数の次の a/b (=IOUT1/IOUT0) ですが,a は,次に出力される1次結果ファイルが何番目のものか,b は,1次結果ファイルを何回出力したら2次結果ファイルが作成されるかを示すもので,2次結果ファイルが作成されると a は1に戻る.その次の数字は,計算開始から何回目の「最小時間単位」の計算を行っているか,その次は,次回作成される1次結果ファイルは通算で何個目かを示す.数字を3つ飛ばした後の数値(おそらく35.0000)は,鋸刃の角度(degree)を示している(角度を独自に指定した場合はその数値).そして, Bird法かU-sys法かの識別,先に指定した初期速度分布の種類と並ぶ.なお,これとほぼ同じ内容(一部は異なる)が,JB-AKTN2-AKTH というファイルに書き込まれ,1次結果ファイルが作成されるごとに更新される.なお,1次結果ファイルが作成されるたびにuvw333やtfile2等のファイルが作られ,これらのファイルが存在すれば,途中から計算を継続することができる.もしフォルダ内にそれらが存在していないと継続計算はできない.それ以外の画面表示は,プログラムの開発過程で利用したもので,それぞれ意味はあるが,その説明は省略する.
(D)定常状態に達したかどうかの判断とそれ以降の結果の平均計算
計算が定常状態に達したかどうかの判断は,続いて得られる結果が,統計的変動の範囲内で変化しなくなったかどうかによるので,面倒な数式による判断を避けるのなら,とりあえずは得られた結果の図を目で比較して,経験によって収束判定をしてもらうしかありません.なお私の経験では,計算が400ステップを超えれば大抵は定常状態に達していると考えられますが,510ステップで計算終了となっているので,それまで待ってもらえば良いのではないかと考えます. 個々の3次結果ファイルは,45ステップの計算の平均で得られたものなので,大きなRe数における結果ならそのまま使えますが,Re数が小さくなると平均流速が小さくなって統計的変動が大きくなります.そこで,3次結果ファイルを寄せ集めて平均を作り,より滑らかな結果にすることが考えられます.それを行うためのツールが,SUM009です.その使い方は,SUM009を立ち上げ,まず,平均を計算するための一番最後となる3次結果ファイルの最終ステップ番号を入力します(3次結果ファイルがOUT-aaa-bbbなら,bbbの数値).次に,その結果ファイルからいくつ遡ったデータを用いるかを,2~4の数字で指定します.一見すると,遡る数値が小さいようにも思えますが,データが重複する2次結果を避けて遡るので,2を指定すると,3つ前の3次結果ファイル(6個の2次結果平均となる)に,3を指定すれば6つ前の3次結果ファイル(9個の平均),4を指定すると12個の平均と言ったように,かなり幅広い平均を取ることになります.その最終結果のファイル名はO22=xxx-yyy=Z となる.ここで,xxxとyyyは2次結果番号,Zは2~4の数値です.
結果表示にtecplotを使う場合では,結果の図の一枚一枚を連続でつなぎ合わせてアニメーションを作成する機能が備わっているが,私の知る限り,gnuplotにはその機能はない.したがって,gnuplotでアニメーションを作成するには,ネットで公開されているアニメーション作成用の無料ソフトを探してそれを利用してください.なお,gnuplot による結果表示に関しては,ここでの説明は省略するので,第2章や第3章・第2節(3.2節),第4章・第2節(4.2節)の説明を読んで下さい.あるいは,ウェブで公開されているgnuplotの説明(数多く存在する)を探して読んでください.