用いたケースファイル
https://github.com/reuscaelum/CFD_cases
これを解説していきます。
実行環境
- ubuntu 18.04
- OpenFOAM 6
- cfMesh-Development
用いたCADファイル
今回の記事用にモデルロケット風なCADデータを用意しました.
機体の進行方向が+Z軸方向です.
解析領域の作成
cfMeshのsurfaceGenerateBoundingBoxコマンドを使用します.(例:surfaceGenerateBoundingBox rocket.stl rocket_box.stl 1.5 1.5 4 4 20 20)
入力するSTLファイルと出力するSTLファイルを指定して,その後に解析領域の-X方向,+X方向,-Y方向,+Y方向,-Z方向,+Z方向(単位はm)の大きさを指定します.物体の後流の流れに影響を及ぼさない程度に大きくしておいたほうがいいらしいです.今回は過剰に取りすぎた感があります.
メッシュ生成
今回使用したmeshdict
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | cfMesh: A library for mesh generation |
| \\ / O peration | |
| \\ / A nd | Author: Franjo Juretic |
| \\/ M anipulation | E-mail: franjo.juretic@c-fields.com |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "system";
object meshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
surfaceFile "rocket_box.stl";
maxCellSize 5.0;
boundaryCellSize 0.1;
minCellSize 3;
localRefinement
{
"modelrocket.*"
{
cellSize 0.0005;
}
}
boundaryLayers
{
//nLayers 15;
//thicknessRatio 1.2;
//maxFirstLayerThickness 0.5;
patchBoundaryLayers
{
"modelrocket.*"
{
nLayers 10;
thicknessRatio 1.02;
maxFirstLayerThickness 0.0001;
allowDiscontinuity 0;
}
}
}
renameBoundary
{
defaultName fixedWalls;
defaultType wall;
newPatchNames
{
"modelrocket.*"
{
newName rocket;
type wall;
}
"zMax"
{
newName inlet;
type patch;
}
"zMin"
{
newName outlet;
type patch;
}
}
}
// ************************************************************************* //
systemの中にあるmeshdictでメッシュの設定をします.まずsurfaceFileで読み込むstlファイルを指定します.この指定するstlは先程解析領域で出力したstlファイルを指定します.
そのあとにメッシュのサイズなどを設定する項目がありますが飛ばします.
localRefinementの項目でstlファイル周辺のメッシュを細分化します.今回は周辺のメッシュを0.5mmに細分化します.
次にboundaryLayersの項目で境界層の設定をします.CFD OnlineのY+ Estimate https://www.cfd-online.com/Tools/yplus.php
などを利用して十分な境界層メッシュが生成されるようにします.今回は乱流モデルとしてK-ω SSTを用いるので10層の境界層メッシュを挿入します.
次に境界層メッシュの成長率と1層目の最大厚さ,不連続なメッシュを許容するかという設定をします.
次に境界の名前と性質を変更します.流入面,流出面を定義しておきます.
meshdictを保存したらcartesianMeshを実行してmeshを生成します.
checkMeshをしてメッシュのアスペクト比や歪度などを見てダメそう(アスペクト比が200超とか)だと設定を変えてメッシュをもう一度生成したほうが良さそうです.
またparaViewで狙い通りのメッシュが生成されているかも見ておいたほうが良さげです.数字上では良くても,壊れているメッシュが生成されていたりします.(今回は大丈夫でした)
今回生成したメッシュ
乱流モデルの選択&初期値の設定
乱流モデルはconstantの中のturbulence Propertiesで選択します.今回はk-ω SSTを用います.詳しい説明はこちらをどうぞhttps://turbmodels.larc.nasa.gov/sst.html
/*--------------------------------*- 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 dictionary;
object turbulenceProperties;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
simulationType RAS;
RAS
{
RASModel kOmegaSST;
turbulence on;
printCoeffs on;
}
// ************************************************************************* //
次に初期値を設定します.今回は-Z方向に100[m/s]の空気が流入,流出面で大気圧開放するようにします.
tutorialのmotorbikeのファイルを持ってきて今回のケースに合うように境界の名前を書き換えたりします.
今回はK-ω SSTを用いるのでUとpの他にk,nut,omegaの初期値を設定します.
/*--------------------------------*- 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 volVectorField;
location "0";
object U;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
dimensions [0 1 -1 0 0 0 0];
internalField uniform (0 0 0);
boundaryField
{
rocket
{
type noSlip;
}
inlet
{
type fixedValue;
value uniform (0 0 -100);
}
outlet
{
type zeroGradient;
}
fixedWalls
{
type noSlip;
}
}
// ************************************************************************* //
/*--------------------------------*- 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 p;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
dimensions [0 2 -2 0 0 0 0];
internalField uniform 0;
boundaryField
{
rocket
{
type zeroGradient;
}
inlet
{
type zeroGradient;
}
outlet
{
type fixedValue;
value uniform 0;
}
fixedWalls
{
type zeroGradient;
}
}
// ************************************************************************* //
次にsystemのcontroldictで解析の設定をします.
applicationでどのソルバーで解析するかを設定します.今回はsimpleFoamを選択します.
次に何回反復計算するかを設定します.simpleFoamを用いるのでendTimeだけ気にすればいいと思います.今回は1000回反復計算させてみます.deltaTは今回は定常計算なので1にします.
次にfunctionsで物体にかかる力のログ出力の設定をします.
またforceCoeffsを読み込んで抗力係数を算出します.
/*--------------------------------*- 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 dictionary;
object controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
application simpleFoam;
startFrom latestTime;
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
{
forces
{
type forces;
functionObjectLibs ("libforces.so");
outputControl timeStep;
timeInterval 1;
log yes;
patches (rocket);
p p;
U U;
rho rhoInf;
log true;
rhoInf 1;
CofR (0 0 0);
}
#include "forceCoeffs"
}
// ************************************************************************* //
次に離散化スキームと代数方程式ソルバーの設定をfvSchemesとfvSolutionで設定します.ここは自分が説明するよりPENGUINITSさんのページを見るのが早いと思うので見てくださいhttp://penguinitis.g1.xrea.com/study/OpenFOAM/index.html
実行させるシェルスクリプトを準備
tutorialからAllrunとAllcleanをコピペして使います.decomposeParを使って並列計算させています.メッシュをみたいだけ用にAllmeshというシェルスクリプトも用意しました.
#!/bin/sh
cd ${0%/*} || exit 1 # Run from this directory
# Source tutorial run functions
. $WM_PROJECT_DIR/bin/tools/RunFunctions
surfaceGenerateBoundingBox rocket.stl rocket_box.stl 1.5 1.5 4 4 20 20
runApplication cartesianMesh
runApplication checkMesh -constant
runApplication decomposePar
runParallel $(getApplication)
runApplication reconstructParMesh -constant
runApplication reconstructPar -latestTime
runApplication foamLog log.simpleFoam
#------------------------------------------------------------------------------
#!/bin/sh
cd ${0%/*} || exit 1 # Run from this directory
# Source tutorial clean functions
. $WM_PROJECT_DIR/bin/tools/CleanFunctions
# Remove surface and features
cleanCase
#------------------------------------------------------------------------------
加えてこの上の階層に次のようなシェルスクリプトを置いて複数のケースを一回のシェルスクリプトの実行で回せるようにします.
#!/bin/sh
cd ${0%/*} || exit 1 # Run from this directory
dir_path="${0%/*}"
dirs=`find $dir_path -maxdepth 1 -type d -not -name .`;
for dir in $dirs;
do
cd $dir
./Allrun
cd ..
done
#------------------------------------------------------------------------------
これを実行すると解析が始まります.
結果
ノーズ後端で剥離してるのがなんとなくわかりますね
感想
技術系記事を書くのが初めてなので至らないところがあったら申し訳ないです.全国のロケッティアが参考になればと思い書いてみました.圧縮性のソルバーもrunApplicationの指定変えて,pの初期値の次元を合わせたらできるので超音速機作る際はそちらを使えばいいと思います.
質問があればtwitter:@8900_blueskyまでどうぞ
次はスペースバルーンの軌道シミュレーションの記事についてor圧縮性ソルバーを用いた解析について書こうと思います.(もっとスペースバルーンやるサークル増えて…)
参考にしたページ
https://damassets.autodesk.net/content/dam/autodesk/www/campaigns/sim-day-2015/Autodesk_CFD_%E8%A7%A3%E6%9E%90%E3%83%86%E3%82%AF%E3%83%8B%E3%83%83%E3%82%AF.pdf
http://penguinitis.g1.xrea.com/study/OpenFOAM/index.html
https://turbmodels.larc.nasa.gov/sst.html