【OpenFOAMのTips】覚えておくと便利だけど忘れがちなもの
こんにちは（@t_kun_kamakiri)
メモに残しておきます。 たまに追加していくので困ったときに知っておくと役に立つかもしれません。
気になる部分は目次からジャンプしてください_(._.)_
OpenFOAMv2212
目次- 流体解析で計算がうまくいかないとき
- bashのプロンプトを変更する
- OpenFOAMグループ
- snappyHexMeshでメッシュが思った通りに生成されていないとき
- 流量を知りたい
- 残差と連続式の誤差を出力
- 実験データと比較したい
- 境界条件を楽にする
- 並列計算
- 特定のフォルダの容量調べる
- 境界条件の候補を調べる
- ソースコードを探す
- 検索したい文字列のファイルを探す
- メッシュ品質の基準を検索
- メッシュ品質の悪い部分を確認
- viコマンド
- 圧縮ファイルの解凍
- tar.gz
- zip
- treeコマンドが使えないとき
- OpenFOAMの環境変数の確認
- ケースファイルの初期化
- postProcess
- sampling
- 壁面摩擦応力(wallShearStress）の出力
- yPlusの出力
- 特定値の流速のみを表示(ParaView)
- WSLでメモリ不足になったとき
- 環境変数が展開されない
- 境界条件の参考
- 非圧縮流れの場合(simpleFoma、pimpleFoamなど）
- 熱流体の場合
- 熱流体の場合（マルチリージョン）
- モデルの形状線を表示
- モデルのスケール変換
- クーラン数の確認
- 参考書
流体解析で計算がうまくいかないとき
OpenFOAMをメインに流体解析でうまくいかないときの原因や対策などのメモを残しています。
ご参考ください_(._.)_
あわせて読みたい【流体解析で発散】解析がうまくいかないとき
bashのプロンプトを変更する
ターミナルのパスが長すぎたりすると改行したり見づらかったりしますよね。 今回は以下のように時刻と改行を入れた設定に変えるようにします。
以下をターミナルで打って、
vi ~/.bashrc 1 vi ~/.bashrc一番最後の行に下記を追加します。
PS1='\ek\e\\\][\e[32m\u@\h \e[35m\t \e[34m\w\e[0m]\n\$ ' 1 PS1='\ek\e\\\][\e[32m\u@\h \e[35m\t \e[34m\w\e[0m]\n\$ '「:wq」と打てば保存してvimから抜けることができます。
↓こちらにするとパスの文字が黄色になります。
PS1='\ek\e\\\][\e[32m\u@\h \e[36m\t \e[33m\w\e[0m]\n\$ ' 1 PS1='\ek\e\\\][\e[32m\u@\h \e[36m\t \e[33m\w\e[0m]\n\$ '好きなカラーで設定してください。
OpenFOAMグループ
こちらに質問と回答が蓄積されているので同じ悩みの方に出会えるかもしれません。
ただ良い回答を得るためには良い質問の仕方があります。
あわせて読みたいOpenFOAM Google Groupに投稿する際に抑えるポイント
snappyHexMeshでメッシュが思った通りに生成されていないとき
- 「locationInMesh (2.1 0.09 -1.1);」の座標がおかしい
- メッシュサイズが大きくて形状変化の解像度をとらえられていない
- stlファイルをよく見ると欠けている
だいたいどれか。
並列計算しているなら結合しないとParaViewでは見れない。 並列分割してフォルダ内にprocessorがあるなら、ParaviewのCase TypeをDecomposed Caseに変えると並列分割した状態でも読み込めます。
流量を知りたい
あわせて読みたい- 【OpenFOAM】境界面での流量計算(volFieldValue)
- 【OpenFOAM】任意の断面での流量計算(cellZone)
- 【OpenFOAM】任意の断面での流量計算(codedFunctionObject)
①に残差の例もあります。
残差と連続式の誤差を出力
system/controlDictには残差や連続式の誤差を出力する設定を記述します。
system/controlDict
application simpleFoam; startFrom startTime; startTime 0; stopAt endTime; endTime 500; deltaT 1.0; adjustTimeStep no; maxCo 0.9; writeControl timeStep; writeInterval 50; purgeWrite 0; writeFormat ascii; writePrecision 6; writeCompression off; timeFormat general; timePrecision 6; runTimeModifiable yes; functions { continuityError1 { type continuityError; libs ( fieldFunctionObjects ); // Mandatory entries (unmodifiable)// Optional entries (runtime modifiable) phi phi; // Optional (inherited) entries writePrecision 8; writeToFile yes; useUserTime yes; region region0; enabled yes; log yes; timeStart 0; timeEnd 6000; executeControl timeStep; executeInterval 1; writeControl timeStep; writeInterval 50; } residuals { type solverInfo; libs ("libutilityFunctionObjects.so"); fields (U p k omega); } //#includeFunc residuals(U,p,k,epsilon); #includeFunc yPlus; #includeFunc CourantNo; } 12345678910111213141516171819202122232425262728293031323334353637383940414243444546474849505152 application simpleFoam;startFrom startTime;startTime 0;stopAt endTime;endTime 500;deltaT 1.0;adjustTimeStep no;maxCo 0.9;writeControl timeStep;writeInterval 50;purgeWrite 0;writeFormat ascii;writePrecision 6;writeCompression off;timeFormat general;timePrecision 6;runTimeModifiable yes; functions{ continuityError1 { type continuityError; libs ( fieldFunctionObjects ); // Mandatory entries (unmodifiable)// Optional entries (runtime modifiable) phi phi; // Optional (inherited) entries writePrecision 8; writeToFile yes; useUserTime yes; region region0; enabled yes; log yes; timeStart 0; timeEnd 6000; executeControl timeStep; executeInterval 1; writeControl timeStep; writeInterval 50; } residuals { type solverInfo; libs ("libutilityFunctionObjects.so"); fields (U p k omega); } //#includeFunc residuals(U,p,k,epsilon); #includeFunc yPlus; #includeFunc CourantNo; }残差と連続の式による誤差は最低確認するようにしましょう。
連続の式の誤差
- sum local : 誤差絶対値の格子体積重み付け平均
- global : 誤差(符号あり)の格子体積重み付け平均
- cumulative : globalの累積
グラフを出力するPythonスクリプトはこちらです。
連続の式の誤差
import numpy as np import matplotlib.pyplot as plt import pandas as pd import os def graph_func(df_continuityError, imin, imax, iposi,icol): plt.subplot(imin,imax,iposi) plt.plot(df_continuityError['time'],df_continuityError[icol],linestyle='solid',color='b',marker="o",label=icol) plt.xlabel('time') plt.ylabel(icol) plt.legend() plt.grid() def df_continuityError_func(dir_): data_continuityError = np.loadtxt(f'./postProcessing/continuityError1/{dir_}/continuityError.dat') df_continuityError = pd.DataFrame(data_continuityError, columns=['time','Local','Global','Cumulative']) return df_continuityError def df_concat(dir_List): df = pd.DataFrame() for dir_ in dir_List: try: df_ = df_continuityError_func(dir_) df = pd.concat([df, df_]) except Exception: print("Empty") return df dir_List = os.listdir('./postProcessing/continuityError1') df_continuityError = df_concat(dir_List) plt.figure(figsize=(20,3)) graph_func(df_continuityError, 1,3,1,'Local') graph_func(df_continuityError, 1,3,2,'Global') graph_func(df_continuityError, 1,3,3,'Cumulative') print('time : ',df_continuityError['time'].iloc[-1], 'continuityError : ', df_continuityError['Local'].iloc[-1]) df_continuityError.head() 1234567891011121314151617181920212223242526272829303132333435363738394041 import numpy as npimport matplotlib.pyplot as pltimport pandas as pdimport os def graph_func(df_continuityError, imin, imax, iposi,icol): plt.subplot(imin,imax,iposi) plt.plot(df_continuityError['time'],df_continuityError[icol],linestyle='solid',color='b',marker="o",label=icol) plt.xlabel('time') plt.ylabel(icol) plt.legend() plt.grid() def df_continuityError_func(dir_): data_continuityError = np.loadtxt(f'./postProcessing/continuityError1/{dir_}/continuityError.dat') df_continuityError = pd.DataFrame(data_continuityError, columns=['time','Local','Global','Cumulative']) return df_continuityError def df_concat(dir_List): df = pd.DataFrame() for dir_ in dir_List: try: df_ = df_continuityError_func(dir_) df = pd.concat([df, df_]) except Exception: print("Empty") return df dir_List = os.listdir('./postProcessing/continuityError1')df_continuityError = df_concat(dir_List) plt.figure(figsize=(20,3))graph_func(df_continuityError, 1,3,1,'Local')graph_func(df_continuityError, 1,3,2,'Global')graph_func(df_continuityError, 1,3,3,'Cumulative') print('time : ',df_continuityError['time'].iloc[-1], 'continuityError : ', df_continuityError['Local'].iloc[-1]) df_continuityError.head()残差
import pandas as pd import matplotlib.pyplot as plt import os def graph_layout(): plt.grid() plt.legend(loc='best',bbox_to_anchor=(1, 1)) plt.yscale('log') plt.xlabel('time', fontsize=16) plt.ylabel('residual', fontsize=16) def df_residual_func(dir_): df_residual = pd.read_table(f'./postProcessing/residuals/{dir_}/solverInfo.dat',skiprows=1) df_residual = pd.DataFrame(df_residual) return df_residual def df_concat(dir_List): df = pd.DataFrame() for dir_ in dir_List: try: df_ = df_residual_func(dir_) df = pd.concat([df, df_]) except Exception: print("Empty") return df dir_List = os.listdir('./postProcessing/residuals') df_residual = df_concat(dir_List) initial_residial = [data for data in df_residual.columns if "initial" in data] final_residial = [data for data in df_residual.columns if "final" in data] df_residual.plot(x="# Time ",y=initial_residial, figsize=(6, 4)) graph_layout() #df_residual.plot(x="# Time ",y=final_residial) 123456789101112131415161718192021222324252627282930313233343536 import pandas as pdimport matplotlib.pyplot as pltimport os def graph_layout(): plt.grid() plt.legend(loc='best',bbox_to_anchor=(1, 1)) plt.yscale('log') plt.xlabel('time', fontsize=16) plt.ylabel('residual', fontsize=16) def df_residual_func(dir_): df_residual = pd.read_table(f'./postProcessing/residuals/{dir_}/solverInfo.dat',skiprows=1) df_residual = pd.DataFrame(df_residual) return df_residual def df_concat(dir_List): df = pd.DataFrame() for dir_ in dir_List: try: df_ = df_residual_func(dir_) df = pd.concat([df, df_]) except Exception: print("Empty") return df dir_List = os.listdir('./postProcessing/residuals')df_residual = df_concat(dir_List) initial_residial = [data for data in df_residual.columns if "initial" in data]final_residial = [data for data in df_residual.columns if "final" in data]df_residual.plot(x="# Time ",y=initial_residial, figsize=(6, 4))graph_layout()#df_residual.plot(x="# Time ",y=final_residial)残差の分布を出力す量にカスタマイズした有益な議論があります。
あわせて読みたい- 収束性悪化の原因箇所の推定方法
- inabower/OpenFOAM4.1-utilities
実験データと比較したい
だいたいここを見ます。
ベンチマーク Classic Collection Database http://cfd.mace.manchester.ac.uk/ercoftac/doku.php wolf http://www.wolfdynamics.com/tutorials.html?id=94 都市部 https://www.aij.or.jp/jpn/publish/cfdguide/ https://thtlab.jp/PTV/ptv_database.html あわせて読みたい- 【OpenFOAM】流体解析用のベンチマークのサイト
- 【流体の実験データ紹介】CAE解析と実験を比較したい。
境界条件を楽にする
共通の設定の場合
"(inlet|outlet)" { type fixedValue; value uniform (0 0 0); } 12345 "(inlet|outlet)"{ type fixedValue; value uniform (0 0 0);}デフォルト(全て壁条件)
".*" { type noSlip; } 1234 ".*"{ type noSlip;}設定の使いまわし
inlet1 { type fixedValue; value uniform (1 0 0); } inlet2 { $inlet1; } 12345678910 inlet1 { type fixedValue; value uniform (1 0 0);} inlet2{ $inlet1;}
並列計算
並列計算
decomposePar reconstructParMesh -constant -mergeTol 1e-6 mpirun -np 4 renumberMesh -overwrite -parallel 123 decomposeParreconstructParMesh -constant -mergeTol 1e-6mpirun -np 4 renumberMesh -overwrite -parallel特定のフォルダの容量調べる
duコマンド
du -sh /mnt/d/OpenFOAM/ #出力結果 50G /mnt/d/OpenFOAM/ 1234 du -sh /mnt/d/OpenFOAM/ #出力結果50G /mnt/d/OpenFOAM/境界条件の候補を調べる
foamHelp boundary -field U 1 foamHelp boundary -field U以下が出力される。
Create mesh for time = 0 Available boundary conditions for vector field: U SRFFreestreamVelocity SRFVelocity SRFWallVelocity acousticWaveTransmissive activeBaffleVelocity activePressureForceBaffleVelocity advective calculated codedFixedValue codedMixed cyclic cyclicACMI cyclicAMI cyclicSlip cylindricalInletVelocity directionMixed empty exprFixedValue exprMixed extrapolatedCalculated fixedGradient fixedInternalValue fixedJump fixedJumpAMI fixedMean fixedMeanOutletInlet fixedNormalInletOutletVelocity fixedNormalSlip fixedProfile fixedShearStress fixedValue flowRateInletVelocity flowRateOutletVelocity fluxCorrectedVelocity freestream freestreamVelocity inletOutlet interstitialInletVelocity kqRWallFunction mapped mappedField mappedFixedInternalValue mappedFixedPushedInternalValue mappedFlowRate mappedMixed mappedMixedField mappedVelocityFlux matchedFlowRateOutletVelocity mixed movingWallVelocity noSlip nonuniformTransformCyclic outletInlet outletMappedUniformInlet outletPhaseMeanVelocity partialSlip pressureDirectedInletOutletVelocity pressureDirectedInletVelocity pressureInletOutletParSlipVelocity pressureInletOutletVelocity pressureInletUniformVelocity pressureInletVelocity pressureNormalInletOutletVelocity pressurePIDControlInletVelocity processor processorCyclic rotatingPressureInletOutletVelocity rotatingWallVelocity scaledFixedValue sliced slip supersonicFreestream surfaceNormalFixedValue swirlFanVelocity swirlFlowRateInletVelocity swirlInletVelocity symmetry symmetryPlane timeVaryingMappedFixedValue translatingWallVelocity turbulentDFSEMInlet turbulentDigitalFilterInlet turbulentInlet uniformFixedGradient uniformFixedValue uniformInletOutlet uniformJump uniformJumpAMI uniformNormalFixedValue variableHeightFlowRateInletVelocity waveTransmissive wedge zeroGradient End 1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162636465666768697071727374757677787980818283848586878889909192939495969798 Create mesh for time = 0 Available boundary conditions for vector field: USRFFreestreamVelocitySRFVelocitySRFWallVelocityacousticWaveTransmissiveactiveBaffleVelocityactivePressureForceBaffleVelocityadvectivecalculatedcodedFixedValuecodedMixedcycliccyclicACMIcyclicAMIcyclicSlipcylindricalInletVelocitydirectionMixedemptyexprFixedValueexprMixedextrapolatedCalculatedfixedGradientfixedInternalValuefixedJumpfixedJumpAMIfixedMeanfixedMeanOutletInletfixedNormalInletOutletVelocityfixedNormalSlipfixedProfilefixedShearStressfixedValueflowRateInletVelocityflowRateOutletVelocityfluxCorrectedVelocityfreestreamfreestreamVelocityinletOutletinterstitialInletVelocitykqRWallFunctionmappedmappedFieldmappedFixedInternalValuemappedFixedPushedInternalValuemappedFlowRatemappedMixedmappedMixedFieldmappedVelocityFluxmatchedFlowRateOutletVelocitymixedmovingWallVelocitynoSlipnonuniformTransformCyclicoutletInletoutletMappedUniformInletoutletPhaseMeanVelocitypartialSlippressureDirectedInletOutletVelocitypressureDirectedInletVelocitypressureInletOutletParSlipVelocitypressureInletOutletVelocitypressureInletUniformVelocitypressureInletVelocitypressureNormalInletOutletVelocitypressurePIDControlInletVelocityprocessorprocessorCyclicrotatingPressureInletOutletVelocityrotatingWallVelocityscaledFixedValueslicedslipsupersonicFreestreamsurfaceNormalFixedValueswirlFanVelocityswirlFlowRateInletVelocityswirlInletVelocitysymmetrysymmetryPlanetimeVaryingMappedFixedValuetranslatingWallVelocityturbulentDFSEMInletturbulentDigitalFilterInletturbulentInletuniformFixedGradientuniformFixedValueuniformInletOutletuniformJumpuniformJumpAMIuniformNormalFixedValuevariableHeightFlowRateInletVelocitywaveTransmissivewedgezeroGradient Endソースコードを探す
全ての$FOAM_APPを出力しておいて、laplacianFoam.Cの検索に引っかかったものだけを抽出
find $FOAM_APP |grep "laplacianFoam.C" 1 find $FOAM_APP |grep "laplacianFoam.C"$FOAM_APPの中で拡張子.Cだけを出力しておいて、laplacianFoam検索に引っかかったものだけを抽出
find $FOAM_APP -name *.C | grep "laplacianFoam" 1 find $FOAM_APP -name *.C | grep "laplacianFoam"検索したい文字列のファイルを探す
# find ファイルパス -type f | xargs grep -n '検索したい文字列' find $FOAM_TUTORIALS -type f| xargs grep -ln "noSlip" 12 # find ファイルパス -type f | xargs grep -n '検索したい文字列'find $FOAM_TUTORIALS -type f| xargs grep -ln "noSlip"- l：ファイル名のみ出力
- n：行番号を出力
メッシュ品質の基準を検索
checkMeshコマンドでメッシュ情報を確認すると、メッシュ数の他にメッシュ品質に関わる内容も出力されます。 しかし、「メッシュ品質って何なの？」となることもあるので、以下に基準値が記載してあるファイルのありかを調べる方法を記します。
find $FOAM_SRC -name "primitiveMeshCheck" 1 find $FOAM_SRC -name "primitiveMeshCheck"以下にあることがわかる。
/opt/OpenFOAM/OpenFOAM-v2012/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshCheck 1 /opt/OpenFOAM/OpenFOAM-v2012/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshCheckviコマンドで中身を確認。
vi /opt/OpenFOAM/OpenFOAM-v2012/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshCheck/primitiveMeshCheck.C 1 vi /opt/OpenFOAM/OpenFOAM-v2012/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshCheck/primitiveMeshCheck.C品質基準は以下を参照。
Foam::scalar Foam::primitiveMesh::closedThreshold_ = 1.0e-6; Foam::scalar Foam::primitiveMesh::aspectThreshold_ = 1000; Foam::scalar Foam::primitiveMesh::nonOrthThreshold_ = 70; // deg Foam::scalar Foam::primitiveMesh::skewThreshold_ = 4; Foam::scalar Foam::primitiveMesh::planarCosAngle_ = 1.0e-6; 12345 Foam::scalar Foam::primitiveMesh::closedThreshold_ = 1.0e-6;Foam::scalar Foam::primitiveMesh::aspectThreshold_ = 1000;Foam::scalar Foam::primitiveMesh::nonOrthThreshold_ = 70; // degFoam::scalar Foam::primitiveMesh::skewThreshold_ = 4;Foam::scalar Foam::primitiveMesh::planarCosAngle_ = 1.0e-6;メッシュ品質は過去にオープンCAEのサマースクールで発表した際の資料にも書いてあります。
メッシュ品質の悪い部分を確認
メッシュ品質が悪い箇所のメッシュを修正したい場合、ParaViewでその場所を特定することができます。
まず、以下でメッシュ品質項目のvtkファイルを作成します。
foamToVTK -faceSet skewFaces -time 0 foamToVTK -faceSet nonOrthoFaces -time 0 12 foamToVTK -faceSet skewFaces -time 0foamToVTK -faceSet nonOrthoFaces -time 0ls constant/polyMesh/setsで品質の悪い項目のファイルを確認。 constant/polyMesh/sets/ファイル名 ファイル名 = skewFaces
次に、paraviewでvtkを読み込むとメッシュ品質の悪い箇所を確認できます。
- 高アスペクト比（Aspect Ratio） 非常に微細な境界層で現れる。 ソルバーの安定性にとって致命的なものではありませんが、収束速度を著しく低下させる可能性があり
- 非直交性（Orthogonality） 3つの値の範囲を定義する ・n0 < 70 – 安全な値 70 < n0 < 90 – fvSolutionのnonOrthoCorrectorやfvSchemesの数値スキームに対して特別な処理を行う必要あり ・n0 > 90 – シミュレーションに使用できない悪いメッシュ
- 歪度(skewness) 値が大きいと結果の品質（正確さ）が損なわれる可能性がありますが、適度な大きさの歪度であれば、シミュレーションに使用することができる
- Smoothness
- メッシュ間隔の最大変化は20%未満が理想
viコマンド
例えばメッシュ品質基準を調べたい場合。
vi /opt/OpenFOAM/OpenFOAM-v2012/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshCheck/primitiveMeshCheck.C 1 vi /opt/OpenFOAM/OpenFOAM-v2012/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshCheck/primitiveMeshCheck.C#行数を表示 :set number
/primitiveMesh n：下へ検索 N(Shift + n)：上へ検索
- :q （上書きせずに閉じる）
- :q! （編集した場合でも上書きせずに閉じる）
- :wq （上書き保存）
- iでその場所から編集
- Aで最終列に移動して編集
- uで戻す
圧縮ファイルの解凍
tar.gz 圧縮 tar -zcvf xxxx.tar.gz directory 解凍 tar -zxvf xxxx.tar.gz 1234 圧縮tar -zcvf xxxx.tar.gz directory解凍tar -zxvf xxxx.tar.gz zip 圧縮 zip -r xxxx.zip directory 解凍 unzip xxxx.zip 1234 圧縮zip -r xxxx.zip directory解凍unzip xxxx.ziptreeコマンドが使えないとき
pwd;find . | sort | sed '1d;s/^\.//;s/\/\([^/]*\)$/|--\1/;s/\/[^/|]*/| /g' 1 pwd;find . | sort | sed '1d;s/^\.//;s/\/\([^/]*\)$/|--\1/;s/\/[^/|]*/| /g'OpenFOAMの環境変数の確認
こちら
echo $FOAM_TUTORIALS # 出力結果 /opt/OpenFOAM/OpenFOAM-v2012/tutorials 123 echo $FOAM_TUTORIALS# 出力結果/opt/OpenFOAM/OpenFOAM-v2012/tutorialsケースファイルの初期化
初期化したい内容に応じてコマンドを使い分ける。
#設定の初期化 foamCleanTutorials #メッシュを初期化 foamCleanPolyMesh #計算結果の初期化 foamListTimes -rm 12345678 #設定の初期化foamCleanTutorials #メッシュを初期化foamCleanPolyMesh #計算結果の初期化foamListTimes -rmpostProcess
controlDictで設定。
Q1 { // Mandatory entries (unmodifiable) type Q; libs (fieldFunctionObjects); // Optional (inherited) entries field <inpField>; result <fieldResult>; region region0; enabled true; log true; timeStart 0; timeEnd 1000; executeControl timeStep; executeInterval 1; writeControl timeStep; writeInterval 1; } 12345678910111213141516171819 Q1{// Mandatory entries (unmodifiable)type Q;libs (fieldFunctionObjects); // Optional (inherited) entriesfield <inpField>;result <fieldResult>;region region0;enabled true;log true;timeStart 0;timeEnd 1000;executeControl timeStep;executeInterval 1;writeControl timeStep;writeInterval 1;}コマンド実行
postProcess -func Q 1 postProcess -func Qとすると時刻のフォルダにQが出力
sampling
samplingによりフィールドデータをサンプリングするためのサンプリング関数オブジェクトが用意されています。グラフにプロットするための1Dライン、画像として表示するための2D平面や3Dサーフェスのいずれかを使用できます。各サンプリングツールは、controlDictファイルのメイン関数辞書、またはcaseシステムディレクトリの別ファイルのいずれかで辞書に指定されています。データは、Grace/xmgr、gnuplot、jPlotのような有名なグラフ作成パッケージを含む様々なフォーマットで書き出すことができます。
以下はsystem/sampleDict内に記述する例を示しています。
sampleDict { type sets; libs (sampling); setFormat raw; interpolationScheme cell; fields (U p wallShearStress); writeControl writeTime; sets ( z_0.001m { type face; axis xyz; start (-0.01 0.01 0.001); end ( 0.01 0.01 0.001); } z_0.05m { type face; axis xyz; start (-0.01 0.01 0.05); end ( 0.01 0.01 0.05); } z_0.1m { type face; axis xyz; start (-0.01 0.01 0.1); end ( 0.01 0.01 0.1); } z_0.15m { type face; axis xyz; start (-0.01 0.01 0.15); end ( 0.01 0.01 0.15); } ); } 1234567891011121314151617181920212223242526272829303132333435363738394041 sampleDict{ type sets; libs (sampling); setFormat raw; interpolationScheme cell; fields (U p wallShearStress); writeControl writeTime; sets ( z_0.001m { type face; axis xyz; start (-0.01 0.01 0.001); end ( 0.01 0.01 0.001); } z_0.05m { type face; axis xyz; start (-0.01 0.01 0.05); end ( 0.01 0.01 0.05); } z_0.1m { type face; axis xyz; start (-0.01 0.01 0.1); end ( 0.01 0.01 0.1); } z_0.15m { type face; axis xyz; start (-0.01 0.01 0.15); end ( 0.01 0.01 0.15); } );}以下のコマンドで結果をpostProcessingに出力。
postProcess -latestTime -func sampleDict | tee log.sample 1 postProcess -latestTime -func sampleDict | tee log.sample壁面摩擦応力(wallShearStress）の出力
simpleFoam simpleFoam -postProcess -func wallShearStress -noZero -noFunctionObjects postProcess -func sampleDict -noZero -latestTime 123 simpleFoamsimpleFoam -postProcess -func wallShearStress -noZero -noFunctionObjectspostProcess -func sampleDict -noZero -latestTime- simpleFoam 計算実行
- simpleFoam -postProcess -func wallShearStress -noZero -noFunctionObjects 計算実行後に実行する。各時刻にwallShearStressを出力される。 -noZero：0ディレクトリは使わないし出力もしない -noFunctionObjects：controlDictのfunctionObjectは使わない。このオプションを付けないとfunctionObjectの記述に従って出力されてしまう。 各時刻のディレクトリ内にwallShearStressがあればParaViewでも表示できる。
- postProcess -func sampleDict -noZero -latestTime postProcessingにsampleDictができて上記のsampleDictのfieldsを出力してくれる。wallShearStressを書いておくと出力される。
※ただしwallShearStressはinterpolationSchemeをcellPointにしないと出力してくれませんでした。
interpolationScheme cellPoint;//cell; 1 interpolationScheme cellPoint;//cell;interpolationScheme
- cell：Cell-centre value assumed constant over cell
- cellPoint：Linear weighted interpolation using cell values
上記同様のことをyPlusで行えばよい。
- 各時間でyPlusを出力 system/controlDict
もしくは
simpleFoam -postProcess -func yPlus -noZero -latestTime 1 simpleFoam -postProcess -func yPlus -noZero -latestTime-latestTime：最後の結果のみ出力するオプション
これでyPlusをParaViewで確認できます。
- system/sampleDict
「fields (U p wallShearStress yPlus);」にyPlusを追加。
postProcess -func sampleDict -noZero -latestTime 1 postProcess -func sampleDict -noZero -latestTimeこれでpostProcessing/postProcessing/最終時刻/z_0.0m_yPlus.xyに結果が出力されます。
0.0450069 23.1191 0.0452398 0 0.0454727 0 0.0457056 0 ... （以下省略） 123456 0.0450069 23.11910.0452398 00.0454727 00.0457056 0...（以下省略）特定の場所のy+の具体的な値をを知りたい場合は便利です。
特定値の流速のみを表示(ParaView)
ParaViewで設定を変えます。 「Edit」>[Settings]から以下の項目にちゅえっくを入れます。
流速のコンターも描くことができます。特定の流速の範囲に絞って表示することができます。
WSLでメモリ不足になったとき
WSLは何も設定しない場合、ホスト(Windows)メモリの 50%または8GBのどちらか少ない方になるように調整します。 Windows側が16GBのメモリの場合は以下のように、
free -h 1 free -hで確認するとトータル7.7BGB使えるようになっており、あふれた分はスワップで2GBまで保持できるようになっています。
シミュレーションをしているとためにメモリが足らなくなるのでもう少し使えるように設定を変えます。 使えるメモリを80%までに引き上げて、残り足らなくなった分は少々遅くなっても良いのでスワップで保持するようにします。
.wslconfig は以下のパスに作成することで、WSL2 起動時に設定が読み込まれる仕組みになっています。
C:\Users\[ユーザ名]\.wslconfig 1 C:\Users\[ユーザ名]\.wslconfig.wslconfigを以下とします。
.wslconfig
[wsl2] memory=13GB swap=10GB 123 [wsl2]memory=13GBswap=10GBPowerShellを管理者権限でシャットダウンします。
wsl --shutdown 1 wsl --shutdown環境変数が展開されない
例えば以下のように$FOAM_TUTORIALSの後にTABを押して、補完機能を使いながらパスをたどっていきたい場合があるとします。
vi $FOAM_TUTORIALS/[TAB] 1 vi $FOAM_TUTORIALS/[TAB]このようにすると、バックスラッシュが先頭に入って面倒だったりします。
「shopt」で、シェルオプションの設定状況（on／off）を一覧表示します。「shopt シェルオプション名」で表示対象を指定できます。
shopt -s direxpand 1 shopt -s direxpandこれで環境変数が自動で展開されます。
vi $FOAM_TUTORIALS/[TAB] ↓以下に展開される vi /opt/OpenFOAM/OpenFOAM-v2012/tutorials/ 12345 vi $FOAM_TUTORIALS/[TAB] ↓以下に展開される vi /opt/OpenFOAM/OpenFOAM-v2012/tutorials/境界条件の参考
非圧縮流れの場合(simpleFoma、pimpleFoamなど） patch type U p k epsilon nut omega 流速指定 patch Region0 { type fixedValue; value uniform (1 0 0); } Region0 { type zeroGradient; } Region0 { type fixedValue; value uniform 0.001; } Region0 { type fixedValue; value uniform 0.01; } Region0 { type calculated; value uniform 0.001; } Region1 { type fixedValue; value $internalField; } 質量流 patch Region { type flowRateInletVelocity; massFlowRate constant 1; //[kg/s] rhoInlet 1000.0; //[kg/m3] } Region { type zeroGradient; } Region0 { type fixedValue; value uniform 0.001; } Region { type fixedValue; value uniform 0.01; } Region { type calculated; value uniform 0.001; } Region1 { type fixedValue; value $internalField; } 体積流量指定 patch Region0 { type flowRateInletVelocity; volumetricFlowRate constant 1; //[m3/s] } Region0 { type zeroGradient; } Region0 { type fixedValue; value uniform 0.001; } Region0 { type fixedValue; value uniform 0.01; } Region0 { type calculated; value uniform 0.001; } Region1 { type fixedValue; value $internalField; } 法線方向流速 patch Region { type surfaceNormalFixedValue; refValue uniform -1.0; } Region { type zeroGradient; } Region { type fixedValue; value uniform 0.001; } Region { type fixedValue; value uniform 0.01; } Region { type zeroGradient; } Region1 { type zeroGradient; } 静止圧指定 patch Region1 { type zeroGradient; } Region1 { type fixedValue; value uniform 0.0; } Region1 { type zeroGradient; } Region1 { type zeroGradient; } Region1 { type calculated; value uniform 0.001; } Region1 { type fixedValue; value $internalField; } 全圧指定 patch Region1 { type zeroGradient; } Region1 { type totalPressure; p0 uniform 0.0; value uniform 0.0; } Region1 { type zeroGradient; } Region1 { type zeroGradient; } Region1 { type calculated; value uniform 0.001; } Region1 { type fixedValue; value $internalField; } 自然流出流入 patch Region2 { type pressureInletOutletVelocity; value uniform (0 0 0); } Region2 { type totalPressure; p0 uniform 0; value uniform 0; } Region2 { type inletOutlet; inletValue uniform 0.001; value uniform 0.001; } Region2 { type inletOutlet; inletValue uniform 0.001; value uniform 0.001; } Region2 { type calculated; value uniform 0.001; } Region1 { type zeroGradient; } 静止壁 wall Region2 { type fixedValue; value uniform (0 0 0); } Region2 { type zeroGradient; } Region2 { type kqRWallFunction; value uniform 0.0; } Region2 { type epsilonWallFunction; value uniform 0.0; } Region2 { type nutkWallFunction; value uniform 0.0; } Region1 { type omegaWallFunction; } 回転壁 type wall; inGroups 1(wall); Region { type rotatingWallVelocity; origin (0 0 0); axis (0 0 1); omega 0.0;//[rad/s] } Region { type zeroGradient; } Region { type kqRWallFunction; value uniform 0.0; } Region { type epsilonWallFunction; value uniform 0.0; } Region { type nutkWallFunction; value uniform 0.0; } 対称面 type symmetry; inGroups 1(symmetry); Region1 { type symmetry; } Region1 { type symmetry; } Region1 { type symmetry; } Region1 { type symmetry; } Region1 { type symmetry; } 滑り条件 patch 熱流体の場合 U T p_rgh p k epsilon nut alphat 流速 type fixedValue;value uniform (1 0 0);
type fixedValue;value uniform 273.15;
type fixedFluxPressure;value uniform 101325;
type calculated;value uniform 101325.0;
type fixedValue;value uniform 0.001;
type fixedValue;value uniform 0.01;
type calculated;value uniform 0.001;
type calculated;value uniform 0.0;
質量流量 type flowRateInletVelocity;massFlowRate constant 1; //[kg/s]
rhoInlet -1.0; //[kg/m3]
type fixedValue;value uniform 273.15;
type fixedFluxPressure;value uniform 101325;
type calculated;value uniform 101325.0;
type fixedValue;value uniform 0.001;
type fixedValue;value uniform 0.01;
type calculated;value uniform 0.001;
type calculated;value uniform 0.0;
体積流量 type flowRateInletVelocity;volumetricFlowRate constant 1; //[m3/s]
type fixedValue;value uniform 273.15;
type fixedFluxPressure;value uniform 101325;
type calculated;value uniform 101325.0;
type fixedValue;value uniform 0.001;
type fixedValue;value uniform 0.01;
type calculated;value uniform 0.001;
type calculated;value uniform 0.0;
法線方向流速 type surfaceNormalFixedValue;refValue uniform -1.0;
type fixedValue;value uniform 273.15;
type fixedFluxPressure;value uniform 101325;
type calculated;value uniform 101325.0;
type fixedValue;value uniform 0.001;
type fixedValue;value uniform 0.01;
type calculated;value uniform 0.001;
type calculated;value uniform 0.0;
静圧 type zeroGradient; type zeroGradient; type fixedValue;value uniform 101325.0;
type calculated;value uniform 101325.0;
type zeroGradient; type zeroGradient; type calculated;value uniform 0.001;
type compressible::alphatWallFunction;value uniform 0.0;
自然流入出 type pressureInletOutletVelocity;value uniform (0 0 0);
type inletOutlet; inletValue uniform 273.15; value uniform 273.15; type totalPressure;p0 uniform 101325;
value uniform 101325;
type calculated;value uniform 101325.0;
type inletOutlet;inletValue uniform 0.001;
value uniform 0.001;
type inletOutlet;inletValue uniform 0.001;
value uniform 0.001;
type calculated;value uniform 0.001;
type calculated;value uniform 0.0;
すべり壁（温度） type slip; type fixedValue;value uniform 273.15;
type fixedFluxPressure;value uniform 101325;
type calculated;value uniform 101325.0;
type zeroGradient; type zeroGradient; type zeroGradient; type compressible::alphatWallFunction;value uniform 0.0;
すべり壁（断熱） type slip; type zeroGradient; type fixedFluxPressure;value uniform 101325;
type calculated;value uniform 101325.0;
type zeroGradient; type zeroGradient; type zeroGradient; type compressible::alphatWallFunction;value uniform 0.0;
回転壁（温度） type rotatingWallVelocity;origin (0 0 0);
axis (0 0 1);
omega 0.0;//[rad/s]
type fixedValue;value uniform 273.15;
type fixedFluxPressure;value uniform 101325;
type calculated;value uniform 101325.0;
type kqRWallFunction;value uniform 0.0;
type epsilonWallFunction;value uniform 0.0;
type calculated;value uniform 0.001;
type compressible::alphatWallFunction;value uniform 0.0;
全圧 type zeroGradient; type zeroGradient; type totalPressure;p0 uniform 101325.0;
value uniform 101325.0;
type calculated;value uniform 101325.0;
type zeroGradient; type zeroGradient; type calculated;value uniform 0.001;
type compressible::alphatWallFunction;value uniform 0.0;
熱流体の場合（マルチリージョン）流体と固体、もしくは流体（気体）と流体（液体）などを行えるマルチリージョンソルバの場合は、各領域の境界面が共通されており以下のように設定します。
k-omegaSSTを使用した場合の例を示します。
U T p_rgh p k omega nut alphat 液体の場合 type fixedValue; value uniform (0 0 0);thicknessLayers (2.3); kappaLayers (80); type compressible::turbulentTemperatureCoupledBaffleMixed; value uniform 295.15; Tnbr T; kappaMethod fluidThermo; kappa none; type fixedFluxPressure; value uniform 101325; type calculated;
value uniform 101325.0;
type kqRWallFunction; value uniform 1.35e-13; type omegaWallFunction; blending binomial2; value uniform 6.7082e-05; type nutkWallFunction; value uniform 0; type compressible::alphatWallFunction; Prt 0.85; value uniform 0;固体の場合はTとpのみで良い。
T p 流速 type externalWallHeatFluxTemperature; mode flux; q 2000000; kappaMethod solidThermo; value uniform 298.15; type calculated; value uniform 101325; externalWallHeatFluxTemperature Property Description Required Default mode ‘power’, ‘flux’ or ‘coefficient’ yes Q Power [W] for mode ‘power’ q Heat flux [W/m^2] for mode ‘flux’ h Heat transfer coefficient [W/m^2/K] for mode ‘coefficient’ Ta Ambient temperature [K] for mode ‘coefficient’ thicknessLayers Layer thicknesses [m] no kappaLayers Layer thermal conductivities [W/m/K] no relaxation Relaxation for the wall temperature no 1 emissivity Surface emissivity for radiative flux to ambient no 0 qr Name of the radiative field no none qrRelaxation Relaxation factor for radiative field no 1 kappaMethod Inherited from temperatureCoupledBase inherited kappa Inherited from temperatureCoupledBase inheritedモデルの形状線を表示
例えばFreeCADなどでmm単位でモデル作成したとします。 それをmodelというフォルダにmodel.stlとして保存した場合を考えます。
FreeCADは数値しか拾っていないので、SI単位系で流体解析を行う場合には0.001倍にしてm単位として計算する必要があります。
- model/model.stl：mm単位
- model/model_m.stl：m単位
このようにしたいとします。
surfaceConvertでスケール変換し、それを特徴線として出力(eMesh)、さらにParaViewで表示するためにobjファイルにするコマンドは以下です。
surfaceConvert -scale 0.001 model/model.stl model/model_m.stl cp -r model/model_m.stl constant/triSurface/ surfaceFeatureExtract surfaceFeatureConvert constant/triSurface/model_m.eMesh constant/triSurface/model_m.obj 1234 surfaceConvert -scale 0.001 model/model.stl model/model_m.stlcp -r model/model_m.stl constant/triSurface/surfaceFeatureExtractsurfaceFeatureConvert constant/triSurface/model_m.eMesh constant/triSurface/model_m.obj特徴線は以下のように設定しておきます。
system/surfaceFeatureExtractDict
model_m.stl { // How to obtain raw features (extractFromFile || extractFromSurface) extractionMethod extractFromSurface; // Mark edges whose adjacent surface normals are at an angle less // than includedAngle as features // - 0 : selects no edges // - 180: selects all edges includedAngle 150; subsetFeatures { // Keep nonManifold edges (edges with >2 connected faces) nonManifoldEdges no; // Keep open edges (edges with 1 connected face) openEdges yes; } // Write options // Write features to obj format for postprocessing writeObj yes; } 1234567891011121314151617181920212223242526 model_m.stl{ // How to obtain raw features (extractFromFile || extractFromSurface) extractionMethod extractFromSurface; // Mark edges whose adjacent surface normals are at an angle less // than includedAngle as features // - 0 : selects no edges // - 180: selects all edges includedAngle 150; subsetFeatures { // Keep nonManifold edges (edges with >2 connected faces) nonManifoldEdges no; // Keep open edges (edges with 1 connected face) openEdges yes; } // Write options // Write features to obj format for postprocessing writeObj yes;}model_m.objはParaViewで読み込めば特徴線も読み込めます。
モデルのスケール変換
OpenFOAMのメッシュを作成した後に節点座標を「移動・回転・スケール」させることができます。 以下の例は、スケールを1/1000倍にした例です。
# 節点の場合 transformPoints -scale "(0.001 0.0001 0.001)" 12 # 節点の場合transformPoints -scale "(0.001 0.0001 0.001)"また、stlを変換することもできます。
# stl変換の場合 surfaceConvert -scale 0.001 constant/triSurface/model.stl constant/triSurface/model_m.stl 12 # stl変換の場合surfaceConvert -scale 0.001 constant/triSurface/model.stl constant/triSurface/model_m.stlクーラン数の確認
クーラン数は$Co = \frac{u\Delta t}{\Delta x}$で定義されます。
陰解法の場合は、安定性の観点からクーラン数を必ずしも1以下にする必要はなく、$\Delta t$があまりにも小さくなりすぎて計算が終わらない場合は、精度を多少犠牲にしてクーラン数を大きくとり$Delta t$が小さくなりすぎないようにします。 しかし、クーラン数の最大値を指定できますが、最小値は指定できません（たぶん） なので、1度過去に計算したもののクーラン数を確認したい場合があります。
以下のコマンドで各時刻（ステップ）ごとのクーラン数を出力します。
postProcess -func CourantNo 1 postProcess -func CourantNoこれでParaViewでのクーラン数の分布を見ることもできます。
クーラン数の最大、最小を見たい場合は以下のようにsystem/controlDictのfunctionsに記述して、fieldの最大最小を出力しておきます。
fieldMinMax { type fieldMinMax; libs ( "libfieldFunctionObjects.so" ); region region0; writeToFile yes; log yes; location yes; mode magnitude; fields (U p Co ); } 123456789101112131415 fieldMinMax { type fieldMinMax; libs ( "libfieldFunctionObjects.so" ); region region0; writeToFile yes; log yes; location yes; mode magnitude; fields (U p Co ); }system/controlDictのfunctionsを書きます。
CourantNo1 { // Mandatory entries (unmodifiable) type CourantNo; libs (fieldFunctionObjects); // Optional entries (runtime modifiable) rho rho; // Optional (inherited) entries field phi; result Co;//CourantNumberField; region region0; enabled true; log true; timeStart 0; timeEnd 1000; executeControl timeStep; executeInterval 1; writeControl writeTime;;//timeStep; writeInterval 1; } 12345678910111213141516171819202122 CourantNo1{ // Mandatory entries (unmodifiable) type CourantNo; libs (fieldFunctionObjects); // Optional entries (runtime modifiable) rho rho; // Optional (inherited) entries field phi; result Co;//CourantNumberField; region region0; enabled true; log true; timeStart 0; timeEnd 1000; executeControl timeStep; executeInterval 1; writeControl writeTime;;//timeStep; writeInterval 1;}記述後にpostProcess のオプション付きで計算を実行します（出力計算のみ）。
interFoam -postProcess 1 interFoam -postProcesspostProcessingに出力されるので、グラフにするなりすればOK。
例えば、以下のPythonプログラムなど。
import pandas as pd import matplotlib.pyplot as plt dir_ = 20 df_Co = pd.read_table(f'./postProcessing/fieldMinMax/{dir_}/fieldMinMax.dat',skiprows=6, sep="\s+", header=None) df_MaxCo = df_Co[[0,1, 6]] df_MaxCo.columns = ["time", "filed", "MaxCo"] df_MaxCo = df_MaxCo[df_MaxCo["filed"] == "Co"] plt.plot(df_MaxCo["time"], df_MaxCo["MaxCo"]) plt.ylim([0,12]) plt.xlabel("time[sec]") plt.ylabel("MaxCo") plt.legend(loc='best',bbox_to_anchor=(1, 1)) plt.title(f'MaxCo') plt.grid() 1234567891011121314151617 import pandas as pdimport matplotlib.pyplot as plt dir_ = 20df_Co = pd.read_table(f'./postProcessing/fieldMinMax/{dir_}/fieldMinMax.dat',skiprows=6, sep="\s+", header=None)df_MaxCo = df_Co[[0,1, 6]]df_MaxCo.columns = ["time", "filed", "MaxCo"]df_MaxCo = df_MaxCo[df_MaxCo["filed"] == "Co"] plt.plot(df_MaxCo["time"], df_MaxCo["MaxCo"]) plt.ylim([0,12])plt.xlabel("time[sec]")plt.ylabel("MaxCo")plt.legend(loc='best',bbox_to_anchor=(1, 1))plt.title(f'MaxCo')plt.grid()参考書
PENGUINさんサイトを体系的に学べる書籍となっています。ネット記事でも十分勉強できるのですが、OpenFOAMの初学者でOpenFOAMをインストール済みであれば一冊持って置き、体系的に学ぶのが良いでしょう。
OpenFOAMによる熱移動と流れの数値解析(第2版)
Amazon楽天市場Yahoo
Amazonの情報を掲載していますあとは初心者向けに丁寧に解説がされているこちらの書籍もお勧めです。最後の章にはオーバーセットメッシュ（重合メッシュ）の機能を使った解析を最後まで丁寧に解説しているので挫折することはないでしょう。
OpenFOAMの歩き方 (技術の泉シリーズ（NextPublishing）)
川畑 真一 発売日: 2021/02/26Amazon楽天市場Yahoo
Amazonの情報を掲載しています辞書的な参考書としては、こちらが良いです。 Foundation版ですが、大変参考になります。
OpenFOAMライブラリリファレンス
株式会社テラバイト, 人見 大輔Amazon楽天市場Yahoo
Amazonの情報を掲載しています ツイート