こんにちは（@t_kun_kamakiri）(‘◇’)ゞ
前回までで、Salomeを使ってモデル作成とメッシュ作成を行いました。
今回は、salomeで作成したメッシュをOpenFOAMを使って流体解析をしたいと思います。
【やる事】
こういった形状での、段差部分にできる渦生成を流体解析で再現したいと思います。
※今回は、OpenFOAMを使ってみた系なので、実現象と比較して精度的にどうかという議論などは行いません。
勉強になるかと思い、層流モデルと乱流モデルの2つのケースで解析をしてみます。
自分自身のまだ理解していない部分があり、かなり説明が抜けてしまっていますのでご了承ください。
それでもOpenFOAMをはじめて使う人にとっての何かのヒントになればと考えています_(._.)_
なので、この記事は以下の人向けに書いています。
- OpenFOAMを使い始めたけど、自分でモデルを作成して流体解析してみたいOpenFOAM初心者
- お金もかけずに自宅のPCで流体解析をやってみたい方
ひとつひとつの設定方法や設定の意味についての解説は行うことができないので、逐一調べるようにしてみてください。
使用環境
- Windows10
- OpenFOAM:ver6.x(WSLでのOpenFOAMのインストール方法はこちら)
ちなみに僕は流体解析の実務経験もなければ、大学で流体を専門にしていたわけではなく完全独学です(笑)
正直全然わかっていないことだらけです・・・・_(._.)_
解くべき方程式
- 運動方程式（2次元）：ナビエストークス方程式
- 質量保存則：非圧縮条件
今回は、定常解析をするため「simpleFOAM」を使います。
simpleFOAMってもともとどんな用途のためのソルバなのか？って思った場合は以下の表を参考にすれば良いでしょう。
OpenFOAMの「非圧縮性流体ソルバ」を以下に示しておきます。
| adjointShapeOptimizationFoam | 定常/非ニュートン流体の乱流/随伴方程式を使用することで障害物による圧損を考慮してダクト形状を最適化 |
| boundaryFoam | 定常/1次元乱流/通常は流入口での境界層条件の生成に使用 |
| icoFoam | 非定常/ニュートン流体の層流 |
| nonNewtonianIcoFoam | 非定常/非ニュートン流体の層流 |
| pimpleFoam | 非定常/時間刻み幅大/（PISOとSIMPLEを組み合わせた）PIMPLEアルゴリズム |
| pimpleDyMFoam | 非定常/移動メッシュ/（PISOとSIMPLEを組み合わせた）PIMPLEアルゴリズム |
| SRFPimpleFoam | 非定常/時間刻み幅大/単一の回転領域 |
| pisoFoam | 非定常/PISOアルゴリズム |
| shallowWaterFoam | 非定常/回転を伴う非粘性の浅水方程式 |
| simpleFoam | 定常 |
| porousSimpleFoam | 定常/多孔質体（陰解法・陽解法）有りの乱流/MRF（Multiple Reference Frame）機能を使用可能 |
| SRFSimpleFoam | 定常/単一の回転領域 |
この中で今回は「simpleFOAM」という、非定常・非圧縮性の流れを扱う際に使うソルバをベースにします。
非圧縮性の流体は基本的には温度場は解きません。
なぜなら、知りたい物理量は「流速\(\boldsymbol{u}\)」と「圧力\(p\)」の２つであり、解くべき方程式「運動量保存則」と「質量保存則」の２つです。
なので、方程式自体が閉じいるので温度を解く必要がないですし、密度変化が小さい非圧縮性の仮定が成り立つような状況において温度変化による密度への影響というのはとても小さいので温度を考える必要がないというのもあります。
温度変化することによる密度変化を考慮しないといけない場合は、有名なのはブシネスク近似と呼ばれるものですが、ここでは深い話はしません。
今回は、乱流モデルを使わない場合と使う場合の2パターンで計算を実行してみます。
正直、実機データと比較しているわけではないのでどちらが実現象に近い解析設定かはわかりませんが、OpenFOAMでの乱流モデルの設定を知るために設定方法を簡単にまとめておきました。
ケースファイルをコピーする
OpenFOAMを使う場合は各ソルバに用意された豊富なチュートリアルを自身の作業ディレクトリにコピーして使うのが基本です。
【作業内容】
- 「/opt/openfoam6/tutorials/incompressible/simpleFoam/motorBike」を自身の作業ディレクトリにコピーしてください。
- 作業ディレクトリの名前は、「motorBike_test10m-s」（写真と名前が違いますが・・・・）としておきます。
では、次に「motorBike_test1」のファイルを編集していきます。
基本的な設定方法はネットに情報があります。
公式ドキュメントとしてはユーザーガイド（日本語訳）が役に立つので、わからないことがあればこちらを参考に理解をしていけば良いでしょう。
境界のタイプの変更
計算領域の境界にあたる「境界面」には、パッチのタイプというのを指定する必要があります。境界面の設定の意味については「ユーザーガイド（日本語訳）」に書いています。
ここで使用する境界面のタイプについては以下のように黄色マーカーをしておきました。
| patch | 一般的なパッチ |
| symmetryPlane | 対称面 |
| empty | 2 次元形状の前後の面 |
| wedge | 軸対称形状のための，くさび型の前後 |
| cyclic | 周期境界面 |
| wall | 壁面（乱流の壁関数に使用） |
| processor | 並列計算時のプロセッサ間の境界 |
一般的な境界面に対しては「patch」と指定するのですが、壁関数や2次元境界面で使用したい場合は、「patch」以外の設定をする必要があります。
【作業内容】
- 「wall」と「FrontAndBack」の境界面の種類を変更
パッチタイプの設定は「polyMesh/boundary」ファイルで指定することができます。
では、次に乱流モデルの設定について見ていきます。
まずは、乱流モデルを使わない場合の設定から見ていきます。
乱流モデル（層流モデル）を使わない場合
乱流モデルを使わない場合というのは、つまりナビエストークス方程式を直接解くということですね。
余計な設定がないため、とても設定が楽です(^^)/
境界条件の設定
境界条件は以下のように設定をします。
OpenFOAMでは以下のように設定するという事を頭に置きつつ、ファイルの編集を行っていきます。
【作業内容】
- 「0/U」と「0/p」の境界条件を設定する。
まず、「0/U」の設定を以下のようにします。
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 | /*--------------------------------*- C++ -*----------------------------------*\ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org \\ / A nd | Version: 6 \\/ M anipulation | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class volVectorField; location "0"; object U; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // #include "include/initialConditions" dimensions [0 1 -1 0 0 0 0]; internalField uniform $flowVelocity; boundaryField { //- Set patchGroups for constraint patches #includeEtc "caseDicts/setConstraintTypes" #include "include/fixedInlet" inlet { type fixedValue; inletValue uniform (10 0 0); } outlet { type zeroGradient; } wall { type noSlip; } #include "include/FrontAndBack" } // ************************************************************************* // |
そして「0/p」の設定を次のようにします。
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 | /*--------------------------------*- C++ -*----------------------------------*\ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org \\ / A nd | Version: 6 \\/ M anipulation | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class volScalarField; object p; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // #include "include/initialConditions" dimensions [0 2 -2 0 0 0 0]; internalField uniform $pressure; boundaryField { //- Set patchGroups for constraint patches #includeEtc "caseDicts/setConstraintTypes" inlet { type zeroGradient; } outlet { type totalPressure; value $internalField; p0 uniform 0; } wall { type zeroGradient; } #include "include/FrontAndBack" } // ************************************************************************* // |
最後の行の「#include “include/FrontAndBack”」は、少々わかりにくいですが、「0/U」でも「0/p」でも同じ記述を書くため、「include/FrontAndBack」というファイルにまとめています☟。
このように同じ記述や初期設定のファイルなどは、別ファイルにまとめて「#include」でファイルを読み込む方がまとまりが良くなるためにこのようにしています。
（チュートリアルの設定がはじめからこのようになっていました）
乱流モデルの設定（層流モデルとして扱う）
今回は、乱流モデルは扱わないため層流の設定にします。
【作業内容】
- 「simulationType」を「laminar」とします。
laminarは層流という意味です。
「constant/turbulenceProperties」
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 | /*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: v1912 | | \\ / A nd | Website: www.openfoam.com | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class dictionary; object turbulenceProperties; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // simulationType laminar; //以下コメントアウト //RAS //{ // RASModel kOmegaSST; // turbulence on; // printCoeffs on; //} // ************************************************************************* // |
このようにすると層流モデルとして扱われますので、物理量に「U」と「p」以外に設定するものはありません。
計算の制御
【作業内容】
- endTime100;
- functionsに「#includeFunc “residuals”」を設定する
「system/controlDict」
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 | /*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: v1912 | | \\ / A nd | Website: www.openfoam.com | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class dictionary; object controlDict; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // application simpleFoam; startFrom startTime; startTime 0; stopAt endTime; endTime 1000; deltaT 1; writeControl timeStep; writeInterval 100; purgeWrite 0; writeFormat binary; writePrecision 6; writeCompression off; timeFormat general; timePrecision 6; runTimeModifiable true; functions { #includeFunc "residuals" // #include "streamLines" // #include "wallBoundedStreamLines" // #include "cuttingPlane" // #include "forceCoeffs" // #include "ensightWrite" } // ************************************************************************* // |
- 今回は定常解析になるので、「endTime1000;」とすることは、最大で1000サイクル計算させるという意味になります。
- 「#includeFunc “residuals”」は残差を出力するための設定です。
残差を確認するためのファイルを用意
流体解析においては残差は必ず確認した方が良いでしょうね。
全然収束していないのにサイクル数が少なすぎて解析が終了してしまうこともあります。
【作業内容】
- 「/opt/openfoam6/etc/caseDicts/postprocessing/numerical」から「residuals」ファイルを、自身の作業ディレクトリの「system」の中にコピーする。
※「system/controlDict」ファイルの「#includeFunc “residuals”」でこのファイルを読み込んでいます。 - 「residuals」ファイルを編集。
Uとpの残差を出力するようにする。
「system/residuals」
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 | /*--------------------------------*- C++ -*----------------------------------*\ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Version: v1912 \\ / A nd | Website: www.openfoam.com \\/ M anipulation | ------------------------------------------------------------------------------- Description For specified fields, writes out the initial residuals for the first solution of each time step; for non-scalar fields (e.g. vectors), writes the largest of the residuals for each component (e.g. x, y, z). \*---------------------------------------------------------------------------*/ #includeEtc "caseDicts/postProcessing/numerical/residuals.cfg" fields (p U); // ************************************************************************* // |
以上で、解析条件の設定が終わりました。
では、計算実行してみしょう(^^)/
計算実行
ターミナル上で、自身作業ディレクトリ上で以下のコマンドにより「simpleFoam」の解析が実行されます。
1 | $simpleFoam |
計算が実行されると以下のように、なにやら記述が出て・・・・「Time=1000（1000サイクル）」で終了します。
計算が終了したら「End」で返ってきます。
するとフォルダの中身は次のようになっています。
「0」「100」・・・・・「1000」というフォルダができています。
これが100サイクルごとの結果ファイルとなります。
結果の可視化
結果の可視化にはparaviewを使います。
ターミナルで「paraFoam」と打って、paraviewを起動します！
1 | $ paraFoam |
☟ターミナルで打つとこんな感じ・・・
残差を確認
残差の確認には、「postProcessing/residuals/0/residuals.dat」にあるdatファイルをgnuplotで出力してあげればよいです。
「datファイル」は普通の空欄区切りのデータなので、gnuplotでグラフ化すれば良いのですが、以下のようなコマンドを打つことですべてのデータをグラフ化してくれます。
※ちなみに、オプション「-l」は縦軸をlogスケールにするという意味です。
1 | $foamMonitor -l postProcessing/residuals/0/residuals.dat |
以下のように残差の時刻歴が出てきます。
やけに振動していますね（笑）
あんまり、残差も小さくない値で収束しながら振動しています・・・・いいんでしょうかこれ（笑）
乱流モデル（k-ωSST）
では、次に乱流モデルを使って流体解析を行います。
乱流モデルは、流れに応じて様々なモデルが提案されていますが、今回使ってみる乱流モデルは「k-ωSST」です。
k-ωSSTの特徴
「自由流れでの乱流モデル化に強いk-εモデル」と「圧力の逆勾配や乖離を伴う境界層流れのモデル化に強いk-ωモデル」とのハイブリッドな乱流モデルです。
要するに、k-εモデルとk-ωモデルの良いところだけを上手く使うモデルってことです。
☟乱流モデルに関しては、以下の標準的な参考書で勉強しておくと良いでしょう。
境界条件
先ほどの層流モデルでの境界条件は「U」と「p」だけでしたが、今回は乱流モデルとしてk-ωSSTを使ってみますので、
- k：乱流エネルギー
\(k=\frac{3}{2}\big(|\boldsymbol{U_in}|I\big)^2\)
\(I\)：乱流強度
※自由乱流の場合はI=0.05(5%の強さの乱れ）が用いられます。 - omega：比散逸率
\(\omega=\frac{\sqrt{k}}{C^{1/4}_{\nu}L}\)
\(C^{1/4}_{\nu}=0.09\)
\(L\):混合長 - nut：乱流粘性係数
\(\nu_{T}=\frac{k}{\omega}\)
の設定が追加で必要になります。
全体の設定をまとめておくと以下のようになります。
乱流モデルの設定に関しては、公式のページや日本語でまとまっているサイトを参考にするとよいでしょう。
OpenFOAM: User Guide v1912 The open source CFD toolbox
乱流モデルの設定 – PENGUINITIS
「0/k」の設定を次のようにします。
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 | /*--------------------------------*- C++ -*----------------------------------*\ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org \\ / A nd | Version: 6 \\/ M anipulation | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format binary; class volScalarField; location "0"; object k; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions [0 2 -2 0 0 0 0]; internalField uniform 0.1; boundaryField { inlet { type turbulentIntensityKineticEnergyInlet; intensity 0.05; value uniform 1; } outlet { type zeroGradient; } wall { type kqRWallFunction; value uniform 0.1; } FrontAndBack { type empty; } } // ************************************************************************* // |
「0/omega」の設定を次のようにします。
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 | /*--------------------------------*- C++ -*----------------------------------*\ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org \\ / A nd | Version: 6 \\/ M anipulation | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format binary; class volScalarField; location "0"; object omega; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions [0 0 -1 0 0 0 0]; internalField uniform 1; boundaryField { inlet { type turbulentMixingLengthFrequencyInlet; mixingLength 0.005; value uniform 200; } outlet { type zeroGradient; } wall { type omegaWallFunction; value uniform 0; } FrontAndBack { type empty; } } // ************************************************************************* // |
「0/nut」の設定を次のようにします。
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 | /*--------------------------------*- C++ -*----------------------------------*\ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org \\ / A nd | Version: 6 \\/ M anipulation | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format binary; class volScalarField; location "0"; object nut; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions [0 2 -1 0 0 0 0]; internalField uniform 0; boundaryField { inlet { type calculated; value uniform 0; } outlet { type calculated; value uniform 0; } wall { type nutkWallFunction; value uniform 0; } FrontAndBack { type empty; } } // ************************************************************************* // |
流入境界条件に付いて簡単な解説を加えておきます。
- k：乱流エネルギー
\(k=\frac{3}{2}\big(|\boldsymbol{U_in}|I\big)^2\)
\(I\)：乱流強度
※自由乱流の場合はI=0.05(5%の強さの乱れ）が用いられます。 - omega：比散逸率
\(\omega=\frac{\sqrt{k}}{C^{1/4}_{\nu}L}\)
\(C^{1/4}_{\nu}=0.09\)
\(L\):混合長 - nut：乱流粘性係数
\(\nu_{T}=\frac{k}{\omega}\)
なので、
☟「0/k」：乱流強度\(I=0.05\)
1 2 3 4 5 6 | inlet { type turbulentIntensityKineticEnergyInlet; intensity 0.05; value uniform 1; } |
☟「0/omega」：混合長\(L=0.005\)
1 2 3 4 5 6 | inlet { type turbulentMixingLengthFrequencyInlet; mixingLength 0.005; value uniform 200; } |
☟「0/nut」：\(\nu_t=k/\omega\)の計算から求まる
1 2 3 4 5 | inlet { type calculated; value uniform 0; } |
※ちなみに「wall」という名前の境界面は壁関数を使っています。
乱流モデルの設定
今回は、乱流モデルは扱わないため層流の設定にします。
【作業内容】
- 「simulationType」を「laminar」とします。
laminarは層流という意味です。
「constant/turbulenceProperties」
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 | /*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: v1912 | | \\ / A nd | Website: www.openfoam.com | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class dictionary; object turbulenceProperties; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // simulationType RAS; RAS { RASModel kOmegaSST; turbulence on; printCoeffs on; } // ************************************************************************* // |
このようにするとk-ωSSTの乱流モデルモデルとして扱われますので、物理量の設定として「0」というフォルダの中には、
- k：乱流エネルギー
- omega：比散逸率
- nut：乱流粘性係数
が必要になるというわけです。
以上で、解析条件の設定が終わりました。
では、計算実行してみしょう(^^)/
計算実行
ターミナル上で、自身作業ディレクトリ上で以下のコマンドにより「simpleFoam」の解析が実行されます。
1 | $simpleFoam |
計算が実行されると以下のように、なにやら記述が出て・・・・「Time=1000（1000サイクル）」で終了します。
計算が終了したら「End」で返ってきます。
するとフォルダの中身は次のようになっています。
「0」「100」・・・・・「1000」というフォルダができています。
これが100サイクルごとの結果ファイルとなります。
結果の可視化
結果の可視化にはparaviewを使います。
ターミナルで「paraFoam」と打って、paraviewを起動します！
1 | $ paraFoam |
☟ターミナルで打つとこんな感じ・・・
残差を確認
残差の確認には、「postProcessing/residuals/0/residuals.dat」にあるdatファイルをgnuplotで出力してあげればよいです。
「datファイル」は普通の空欄区切りのデータなので、gnuplotでグラフ化すれば良いのですが、以下のようなコマンドを打つことですべてのデータをグラフ化してくれます。
※ちなみに、オプション「-l」は縦軸をlogスケールにするという意味です。
1 | $foamMonitor -l postProcessing/residuals/0/residuals.dat |
以下のように残差の時刻歴が出てきます。
乱流モデルを使わなかった場合と比較すると残差が小さな値で収束していますね。
流体解析の実務経験がないため、残差としての考え方があまりわかっていません_(._.)_勉強します。。。
☟解析結果のファイルはこちらからどうぞ！
motorBike_test10m-s_k-omegaSST
あまり良いモデルとはいいがたいです(笑)
まとめ
ひとまず、エラーなくお金もかけずに流体解析を自宅のPCのみで実現できました(^^)/
前回からのシリーズものになっているので、是非☟こちらの記事もご参考ください。
おすすめ参考図書
☟こちらは、OpenFOAMの日本語書籍が無い中唯一わかりやすい参考書だと思います。
☟こちらもOpenFOMの古いバージョンでの和訳になります。さすがにこちらはバージョンに対する日本語でのケアはしていないので、OpenFOAMに慣れている方は購入しても良いかと思います。僕は「日本語でまとまっている内容」なので少し重宝しています。
☟本書で出てきた乱流モデルに関しては、こちらの参考書を推薦図書として挙げておきます。
乱流モデルについて理論面からの解説もわかりやすく、はじめて読む人も十分理解できる標準的な参考書だと思います。
☟以下に、もっと初心者向けの同人誌を紹介しておきます。
初心者は「ってか、まずどうやってOpenFOAMをインストールするの？」というところからつまずきがちです。
そんな時は、以下の書籍をおすすめします。
改訂新版 OpenFOAMの歩き方 (技術の泉シリーズ（NextPublishing）)
インストール方法とチュートリアルで流体解析を体験・・・ちょっと高度な解析まで解説があります。著者曰くOpenFOAMのバージョンの追跡を行いながら、書籍をアップデートするようなので安心ですね。
もっと易しい内容でOpenFOAMを体験したい方向け
OpenFOAMにはこれといったGUIが用意されているわけではないので、基本的にテキストベースで編集する必要があります。
しかし、無料で使えるOpenFOAM用のGUIを作ってくれている人もいるので、それを使わない手はないと思います。
そこで、以下の記事に初心者でも使いやすいXsimという無料のツールを作っている方がいるので使わせてもらい、記事にまとめてみました。
Follow @t_kun_kamakiri