Go to Qiita Advent Calendar 2024 Top
LoginSignup
2
3

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

More than 5 years have passed since last update.

OpenFOAMAdvent Calendar 2017

Day 22
@masazzz(Masashi Imano)

標準的van Driest減衰関数型LESフィルタ幅の実装

Last updated at Posted at 2017-12-22

標準的van Driest減衰関数型LESフィルタ幅の実装

本記事は，オープンCAE学会V&V委員会の第2回 OpenFOAMコード検証勉強会での「OpenFOAMにおけるvan Driest減衰関数の標準的実装とSmagorinsky LESモデルによるチャネル流れの解析講習」(今野雅)の発表資料(LuaTex形式)の一部をpandocを用いてMarkdown形式に変換し，修正したものです．

はじめに

OpenFOAM における LES フィルタ幅 vanDriestDeltaは，幾何的な LES フィルタ幅に van Driest 減衰関数を乗ずる一般的な van Driest 減衰関数型 LES フィルタ幅の定義と異なる．ここでは，OpenFOAM-v1612+ に対して，一般的な van Driest 減衰関 数型 LES フィルタ幅の実装を行い，チャネル流のチュートリアルで動作テストする．

van Driest減衰関数型LESフィルタ幅のコードのコピー

まず， カスタムライブラリのソースコードを置くディレクリを作成する．なお， 場所やディレクトリ名は任意で良い．

run
cd ..
mkdir src
cd src

元のvanDriestDeltaのコードをコピーし， ファイル名を変更する． なお，#で始まる行は前のコマンドの出力である．

find $FOAM_SRC -name vanDriestDelta
# /FullPathToSrc/TurbulenceModels/turbulenceModels/LES/LESdeltas/vanDriestDelta

cp -a `find $FOAM_SRC -name vanDriestDelta` vanDriestStandardDelta
cd vanDriestStandardDelta
ls
# vanDriestDelta.C  vanDriestDelta.H

mv vanDriestDelta.C vanDriestStandardDelta.C
mv vanDriestDelta.H vanDriestStandardDelta.H

wmakeを用いてソースコードをコンパイルするにはMakeディレクトリが必要となるが，元のvanDriestDeltaのディレクトリの上流のディレクトリのどこかにvanDriestDelta用のMakeディレクトリがあるので， これをコピーする．

cp -a $FOAM_SRC/TurbulenceModels/turbulenceModels/Make .

gitを用いたコード履歴管理

これからソースコードやMake内の設定ファイルを修正していくので，変更点を記録するために， gitを用いて管理する．

git init
# Initialized empty Git repository in /FullPath/src/vanDriestStandardDelta/.git/

git add .
git commit -m "Initial commit."
# master (root-commit) f6455fc] Initial commit.
#  4 files changed, 363 insertions(+), 0 deletions(-)
#  create mode 100644 Make/files
#  create mode 100644 Make/options
#  create mode 100644 vanDriestStandardDelta.C
#  create mode 100644 vanDriestStandardDelta.H

ここで， mオプションで指定するコミットメッセージは任意である． また，OpenFOAMではコンパイルする際に， lnIncludeディレクトリを作成し，その下にソースコードがインクルードするファイルへのシンボリックを作成する．さらに， Makeディレクトリの下に，linux64GccDPInt32Optのようなオブジェクトファイルなどを置くディレクトリを作成する．これらのディレクトリ下のファイルは自動で生成されるものなので，以下のような内容で.gitignoreを作成して，gitが無視するように設定しておくと便利である．

.gitignore
lnInclude/
Make/*/

.gitignore についてもgitで管理するように設定する．

git add .gitignore
git commit -m "Add .gitignore."

コンパイルするファイルのリストと生成ライブラリ名の定義

Make/files にはコンパイルするファイル名リストと， 生成ライブラリ名を定義するので，例えば以下のように変更する． ここで， $(FOAM_USER_LIBBIN)は，ユーザーのカスタムライブラリを置くディレクトリを示す環境変数である．なお， ライブラリ名は任意である．

Make/files
vanDriestStandardDelta.C

//LIB = $(FOAM_LIBBIN)/libturbulenceModels
LIB = $(FOAM_USER_LIBBIN)/libvanDriestStandardDelta

クラス名の変更

ソースコード中のクラス名をvanDriestDeltaからvanDriestStandardDeltaに変更する．また，ランタイム型名をvanDriestStandardDeltaからvanDriestStandardに変更する．

sed -i s/vanDriest/vanDriestStandard/g vanDriestStandardDelta.H
sed -i s/vanDriest/vanDriestStandard/g vanDriestStandardDelta.C

gitで変更箇所を確認しておく．

git diff

インクルードパスの設定

本格的に実装を変更する前に， コンパイル可能か調べる．

wmake
# wmakeLnInclude: linking include files to ./lnInclude
# Making dependency list for source file vanDriestStandardDelta.C
# could not open file LESdelta.H for source file vanDriestStandardDelta.C due to No such file or directory

LESdelta.Hが無いというエラーになるので， このファイルを捜す．

find $FOAM_SRC -name LESdelta.H
# /FullPath/src/TurbulenceModels/turbulenceModels/LES/LESdeltas/LESdelta/LESdelta.H
# /FullPath/src/TurbulenceModels/turbulenceModels/lnInclude/LESdelta.H

2つのファイルが見付かっているが， 1つめがファイルの実体であり，2つめのlnInclude以下のファイルは，実体ファイルにシンボリックリンクしたものである．コンパイル時にLESdelta.Hが見付かるよう， Make/filesのEXE_INCに，LESdelta.Hのインクルードパスを追加すれば良いのであるが， 以下に示すように，通常はlnIncludeのディレクトリを指定する． これは， lnIncludeディレクトリには，他のディレクトリに散在するインクルードファイルを，シンボリックリンクという形で集約しているので，実体ファイルのディレクトリを指定するより，インクルードパスの指定が少なくて済むためである． また，デバッグコードを有効にするために， -DDEBUG を追加して， DEBUGを定義にしておく．

Make/files
EXE_INC = \
    -DDEBUG \
    -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
    -I$(LIB_SRC)/finiteVolume/lnInclude \
    -I$(LIB_SRC)/meshTools/lnInclude

LIB_LIBS = \
    -lfiniteVolume \
    -lmeshTools

コード変更前のコンパイル・実行テスト

コンパイル可能かテストする．

wmake

実装を修正する前に， vanDriestStandardDeltaがソルバーで使用可能かテストするため，
channel395のチュートリアル・ケースでテストする． まず，channel395チュートリアルをコピーする．

cp -a $FOAM_TUTORIALS/incompressible/pimpleFoam/LES/channel395 .
pushd channel395

次に， system/controlDictのlibでカスタムライブリlibvanDriestStandardDelta.soを指定する． なお， .soは共有ライブラリ(shared object)の拡張子である． また，1時間ステップで終了して結果を出力するように， endTimeとwriteInterval を変更する．

system/controlDict
libs ("libvanDriestStandardDelta.so");

application     pimpleFoam;

startFrom       startTime;

startTime       0;

stopAt          endTime;

//endTime         1000;                                                                       
endTime         0.2;

deltaT          0.2;

writeControl    timeStep;

//writeInterval   200;                                                                        
writeInterval   1;

さらに， constant/turbulencePropertiesにおいて， LESModelにSmagorinskyを指定し， deltaにvanDriestStandardを指定する．

constant/turbulenceProperties
LES
{
    //    LESModel        WALE;                                                               
    LESModel        Smagorinsky;

    turbulence      on;

    printCoeffs     on;

    //    delta           cubeRootVol;                                                        
    delta           vanDriestStandard;

また， vanDriestStandardの係数設定用のvanDriestStandardCoeffsを用意する．

constant/turbulenceProperties
    //    vanDriestCoeffs                                                                     
    vanDriestStandardCoeffs
    {

ソルバーpimpleFoamが実行できれば， vanDriestStandardのカスタムライブラリが有効になっている事が確認できる．

blockMesh
pimpleFoam

LESフィルタ幅の実装変更

vanDriestStandard のソースコードのディレクトリに復帰する．

popd

先の資料で示した通り，一般的にはvan Driestの減衰関数型フィルター幅Δvは以下で定義される．

Δv = [ 1 − exp (−y+/A+)) ] Δgeo

ここで， y+は無次元壁座標，A+は係数でOpenFOAMのデフォルト値は26，*Δgeo*は幾何的グリッドフィルタ幅であり，通常は格子体積の1/3乗である． 上式を実装するには，vanDriestDeltaのソースを例えば以下のように修正すれば良い． なお，#ifdef DEBUGと#endif間はデバッグ用コードであるため， 必須でないが，このコードはグリッドフィルタ幅deltaの最小値と最大値を表示し，場のファイルを出力すると共に， 格子毎の各項を表示する．

vanDriestStandardDelta.C
    scalar cutOff = wallPointYPlus::yPlusCutOff;
//    wallPointYPlus::yPlusCutOff = 500;
    wallPointYPlus::yPlusCutOff = HUGE;
    wallDistData<wallPointYPlus> y(mesh, ystar);
    wallPointYPlus::yPlusCutOff = cutOff;

/*
    delta_ = min
    (
        static_cast<const volScalarField&>(geometricDelta_()),
        (kappa_/Cdelta_)*((scalar(1) + SMALL) - exp(-y/ystar/Aplus_))*y
    );
*/

    const volScalarField& geometricDelta = geometricDelta_();
    delta_ = (scalar(1) - exp(-y/ystar/Aplus_)) * geometricDelta;

#ifdef DEBUG
    Info << "vanDriestDelta min/max = " << min(delta_).value() << ", "
         << max(delta_).value() << endl;
    delta_.write();
    forAll(delta_, deltaI)
    {
        Info << "C: " << mesh.C()[deltaI]
            << " y: " << y[deltaI]
            << " ystar: " << ystar[deltaI]
            << " fv: " << scalar(1)-exp(-y[deltaI]/ystar[deltaI]/Aplus_)
            << " deltaGeo: " << geometricDelta[deltaI]
            << " delta: " << delta_[deltaI]
            << endl;
    }
#endif

また， 変数Cdelta_は不要となったので， この変数に関連する部分を修正する．なお， 以下では修正がわかりやすいようにコメントアウトしているが，最終的にはコメントアウトしている行は消去して構わない．

vanDriestStandardDelta.H
    //        scalar Cdelta_;
vanDriestStandardDelta.H
/*                                                                                                                      
    Cdelta_                                                                                                             
    (                                                                                                                   
        dict.subDict(type() + "Coeffs").lookupOrDefault<scalar>("Cdelta", 0.158)                                        
    ),                                                                                                                  
*/
(snip)
   //    coeffsDict.readIfPresent<scalar>("Cdelta", Cdelta_);

実行テスト

channel395チュートリアルでソルバーを実行する．

pushd channel395
pimpleFoam > log.pimpleFoam
more log.pimpleFoam
# vanDriestDelta min/max = 1.04361e-29, 0.088154
# C: (0.05 0.00480001 0.0333333) y: 0.00480001 ystar: 0.00120188 fv: 0.14239 deltaGeo: 0.04 delta: 0.0056956
# C: (0.15 0.00480001 0.0333333) y: 0.00480001 ystar: 0.00120847 f: 0.141671 deltaGeo: 0.04 delta: 0.00566684

格子中心の座標ベクトルCのy成分と， 壁座標のyが一致することを確認する．
また， van Driestの減衰関数fv = 1 − exp (−y/(y*A+))，およびLESグリッドフィルタ幅Δ = fvΔgeoが正しく計算されているか確認する．なお，channel395 チュートリアルでのA+は26である．

##LESグリッドフィルタ幅の分布確認

paraviewを用いて， LESグリッドフィルタ幅delta の分布を確認したのが， 以下の図である．

グリッドフィルタ幅の分布
図：グリッドフィルタ幅の分布

デバッグコードを除いたコンパイル

LESグリッドフィルタ幅deltaの算出結果が正しい事を確認したら，デバッグコードを除いて再コンパイルするため，ソースコードのディレクトリに復帰する．

popd

Make/optionsのEXE_INCを以下のように変更する．

Make/options
EXE_INC = \
    /* -DDEBUG */ \
    -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
    -I$(LIB_SRC)/finiteVolume/lnInclude \
    -I$(LIB_SRC)/meshTools/lnInclude

再コンパイルする．

wmake

最後に，gitを用いてコミットする．

git add Make/ vanDriestStandardDelta.*
git commit -m "Implement standard van Driest damping function type."

まとめ

van Driest 減衰関数型 LES フィルタ幅として一般的な実装を行なった．OpenFOAMにおいて，乱流モデル関係ライブラリのカスタマイズは容易な部類に入るので，カスタマイズ初心者の方は是非試してみてください．

2
3
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
Sign upLogin
2
3

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?