はじめに

先日OpenFOAM-ｖ1712がリリースされました。
自分は今までversion5.xを使っていたため、plusシリーズの機能の違いを追うことを兼ねてv1712での変更点についていろいろ試しています。
ここではCHT関連での変更点である "sub-looping" についてまとめてみました。

Release Note 直訳

(原文)
New energy sub-looping for CHT solvers
Sub-looping has been added to the Conjugate Heat Transfer (CHT) solvers as a convergence aid for the energy > solution.
The sub-looping is defined, e.g. in the SIMPLE dictionary in the $FOAM_CASE/system/fvSolution file:

CHTソルバーのための新機能 "sub-looping"
CHTのエネルギー計算の収束性向上のためにsub-loopingが実装された。
sub-loopingは$FOAM_CASE/system/fvSolutionに以下のように記述される。

fvSolution

SIMPLE 
{ 
    energyCoupling 
    { 
        iterations      10; 
        // timeStart       0; 
        // timeEnd         100; 
        interval        0; 

        // Convergence criteria to stop looping 
        convergence 
        { 
            h           1e-3; 
        } 

        // Names of function objects to fire with execute(int) when looping 
        onLoop          ( ); 

        // Names of function objects to fire with execute(int) when converged 
        onConverged     ( externalCoupled ); 

        // Names of function objects to fire with execute(int) when loop ends 
        // without convergence 
        onEnd           ( ); 
    } 
}

sub-looping機能の中身

これがどんなものかを知るためにソースコードを見てみる。

src/finiteVolume/cfdTools/general/solutionControl/loopControl/loopControl.C

bool Foam::loopControl::loop()
{
    bool active = (index_ < total_);   // as per status()

    if (active)
    {
        operator++();

        converged_ = checkConverged();

        if (converged_)
        {
            time_.functionObjects().execute(onConverged_, index_);
            stop();
            return false;
        }
        else if
        (
            interval_ && !(index_ % interval_)
         && !onLoop_.empty()
        )
        {
            time_.functionObjects().execute(onLoop_, index_);
        }
    }
    else if (index_)
    {
        // Not active, the loop condition has now exiting on the last subloop

        if (!converged_ && !onEnd_.empty())
        {
            time_.functionObjects().execute(onEnd_, index_);
        }
    }

    return active;
}

このloop()という関数がこれの肝のようだ。

上限回数以内 && 未収束 → functionObject(onLoop)を実行 → 継続
上限回数以内 && 収束 → functionObject(onConverged)を実行 → 終了
最終回 → functionObject(onEnd)を実行 → 終了

となっている。
つまり 3種類の発動条件に分けてfunctionObjectをそれぞれ仕込む ことができる。
これがCHTの収束性にどのよう繋がっていくのかを確認していく。

tutorial

実際にどのように使っているかを用意されたtitorialを実行して確認してみる。
tut/heatTransfer/chtMultiRegionSimpleFoam/externalCoupledHeater
ちなみに実行した時のlogファイルはこちら。

初期条件

ヒーターの底面は500Kで固定。
上と下に空気と水がそれぞれ流れる構造のモデル。

実行結果

とりあえず動かしてみる。
比較のためにsubLoop（externalCoupled）がない条件でもやってみた。

結果にほとんど違いは見られない。

残差

実行中の残差を見てみる。

収束性は確かに良くなっている（気がする）。

実行ログ

subLoopありのログの最後の方は

log.chtMultiRegionSimpleFoam

Time = 497


Solving for fluid region bottomWater
DILUPBiCGStab:  Solving for Ux, Initial residual = 3.2629515e-16, Final residual = 3.2629515e-16, No Iterations 0
DILUPBiCGStab:  Solving for Uy, Initial residual = 5.5632581e-16, Final residual = 5.5632581e-16, No Iterations 0
DILUPBiCGStab:  Solving for Uz, Initial residual = 1.5460601e-16, Final residual = 1.5460601e-16, No Iterations 0
DILUPBiCGStab:  Solving for h, Initial residual = 9.1113794e-06, Final residual = 1.5077409e-07, No Iterations 1
Min/max T:300 486.15836
GAMG:  Solving for p_rgh, Initial residual = 9.0325089e-08, Final residual = 9.0325089e-08, No Iterations 0
time step continuity errors : sum local = 1.1849235e-06, global = -6.7494025e-08, cumulative = -637.7413
Min/max rho:1000 1000

Solving for fluid region topAir
DILUPBiCGStab:  Solving for Ux, Initial residual = 0.0048056485, Final residual = 9.9021778e-06, No Iterations 1
DILUPBiCGStab:  Solving for Uy, Initial residual = 0.0047238216, Final residual = 5.0157246e-06, No Iterations 1
DILUPBiCGStab:  Solving for Uz, Initial residual = 0.017141788, Final residual = 2.8018083e-05, No Iterations 1
DILUPBiCGStab:  Solving for h, Initial residual = 0.0011804489, Final residual = 1.9577222e-05, No Iterations 1
Min/max T:300 403.29228
GAMG:  Solving for p_rgh, Initial residual = 0.014836675, Final residual = 5.2483909e-05, No Iterations 3
time step continuity errors : sum local = 9.8117121e-05, global = 4.7875108e-05, cumulative = -637.74125
Min/max rho:0.86187249 1.1586213

Solving for solid region heater
DICPCG:  Solving for h, Initial residual = 2.219836e-05, Final residual = 1.7029064e-06, No Iterations 1
Min/max T:366.37296 500

Solving for solid region leftSolid
DICPCG:  Solving for h, Initial residual = 0.0088343281, Final residual = 0.00017768469, No Iterations 1
Min/max T:300.45412 300.95541

Solving for solid region rightSolid
DICPCG:  Solving for h, Initial residual = 0.00025888777, Final residual = 8.2377129e-06, No Iterations 1
Min/max T:348.64113 388.62024
ExecutionTime = 4.55 s  ClockTime = 1242 s

externalCoupled: writing data to "/home/inabower/OpenFOAM/tutorials-v1712/heatTransfer/chtMultiRegionSimpleFoam/externalCoupledHeater/comms/heater_topAir/coupleGroup/T.out"
externalCoupled: removing lock file
externalCoupled: waiting for lock file to appear "/home/inabower/OpenFOAM/tutorials-v1712/heatTransfer/chtMultiRegionSimpleFoam/externalCoupledHeater/comms/OpenFOAM.lock"
externalCoupled: wait time = 1
externalCoupled: wait time = 2
externalCoupled: wait time = 3
externalCoupled: wait time = 4
externalCoupled: wait time = 5
externalCoupled: wait time = 6
externalCoupled: wait time = 7
externalCoupled: wait time = 8
externalCoupled: wait time = 9
externalCoupled: wait time = 10
externalCoupled: found lock file "/home/inabower/OpenFOAM/tutorials-v1712/heatTransfer/chtMultiRegionSimpleFoam/externalCoupledHeater/comms/OpenFOAM.lock"
externalCoupled: removing data files written by master
externalCoupled: reading data from "/home/inabower/OpenFOAM/tutorials-v1712/heatTransfer/chtMultiRegionSimpleFoam/externalCoupledHeater/comms/heater_topAir/coupleGroup/T.in"
Time = 498


Solving for fluid region bottomWater

(中略)

Time = 500


Solving for fluid region bottomWater
DILUPBiCGStab:  Solving for Ux, Initial residual = 3.3122644e-16, Final residual = 3.3122644e-16, No Iterations 0
DILUPBiCGStab:  Solving for Uy, Initial residual = 5.4676824e-16, Final residual = 5.4676824e-16, No Iterations 0
DILUPBiCGStab:  Solving for Uz, Initial residual = 1.6142154e-16, Final residual = 1.6142154e-16, No Iterations 0
DILUPBiCGStab:  Solving for h, Initial residual = 8.8495017e-06, Final residual = 1.4645634e-07, No Iterations 1
Min/max T:300 486.15898
GAMG:  Solving for p_rgh, Initial residual = 9.0325088e-08, Final residual = 9.0325088e-08, No Iterations 0
time step continuity errors : sum local = 1.1849235e-06, global = -6.7494026e-08, cumulative = -637.74119
Min/max rho:1000 1000

Solving for fluid region topAir
DILUPBiCGStab:  Solving for Ux, Initial residual = 0.0018100326, Final residual = 2.9745702e-06, No Iterations 1
DILUPBiCGStab:  Solving for Uy, Initial residual = 0.0027191246, Final residual = 7.4261645e-06, No Iterations 1
DILUPBiCGStab:  Solving for Uz, Initial residual = 0.006608857, Final residual = 1.2510442e-05, No Iterations 1
DILUPBiCGStab:  Solving for h, Initial residual = 0.00063135443, Final residual = 1.4960139e-05, No Iterations 1
Min/max T:300 403.29857
GAMG:  Solving for p_rgh, Initial residual = 0.0052644461, Final residual = 4.1412996e-05, No Iterations 3
time step continuity errors : sum local = 7.6304057e-05, global = -8.6578293e-06, cumulative = -637.7412
Min/max rho:0.86185904 1.1586213

Solving for solid region heater
DICPCG:  Solving for h, Initial residual = 2.1549816e-05, Final residual = 1.6529656e-06, No Iterations 1
Min/max T:366.38235 500

Solving for solid region leftSolid
DICPCG:  Solving for h, Initial residual = 0.0086419156, Final residual = 0.00017381682, No Iterations 1
Min/max T:300.45052 300.95162

Solving for solid region rightSolid
DICPCG:  Solving for h, Initial residual = 0.00025156745, Final residual = 8.1000426e-06, No Iterations 1
Min/max T:348.65406 388.62803
ExecutionTime = 4.59 s  ClockTime = 1252 s

externalCoupled: removing data files written by master
externalCoupled: removing data files written by slave
externalCoupled: lock file status=done
End

のようになっている。
今回のfunctionObjectのexternalCoupledでは待機時間が存在するためClockTimeがかかってしまっているが今回は本筋ではないためこの辺には触れない。

sub-loopingの特徴として、このexternalCoupledの発動が各timeStepではなくおよそ5〜10step毎に行われている。
これはhの計算のresidualがfunctionObjectの発動条件として使われていることが理由であり、おそらく普通に毎timeStepでexternalCoupledを実行させるよりも効率的なのだろう。

最後に

他のfunctionObjectでも色々使い方が出てきそうに感じた。
一方で今回使われたexternalCoupledはplusでは以前から入っていたもののようであり、有用そうなのでもう少し掘り下げてみたい。