計算エラーについて

[長濱谷翔伍]

openfoam-v2012の buoyantPimpleFoamを用いて、STLｓにあるような形状内部の温度分布を解析しようとしました。初期条件300K、heatに温度境界条件として370K、innerに速度境界条件として100rad/sの回転を与えています。Allrun.preでSTLsを統合し、mergedSTL.fmsと mergedSTL.stlを出力しました。Foamrun3Dでメッシュの作成、並列分割（processor0 , processor1)、ソルバーの実行を行っています。（Foamrunは2D用のため、使用しておりません。）

計算は途中まで出来ているように見えますが、止まってしまい以下のようなエラーが出ます。また、  buoyantPimpleFoamのチュートリアルを参考に作ったデータもいくつかあるため、条件に誤りがある可能性もあります。

--------------------------------------------------------------------------
Primary job  terminated normally, but 1 process returned
a non-zero exit code. Per user-direction, the job has been aborted.
--------------------------------------------------------------------------
--------------------------------------------------------------------------
mpirun noticed that process rank 0 with PID 0 on node DESKTOP-O7QGVCM exited on signal 8 (Floating point exception).
--------------------------------------------------------------------------

全てのファイルを添付しましたので、どなたかわかる方はご教授頂きたいです。

宜しくお願い致します。

[Masashi Imano]

今野です。

問題が並列計算によるものか、それ以外かを切り離すために、

こちらでは非並列計算をしてみましたが、エラーは以下となりました。

なお、Foamrun3DでのソルバはbuoyantBoussinesqPimpleFoam となっているので、

buoyantPimpleFoamではなく、buoyantBoussinesqPimpleFoam を用いました。

また、p_rghの線形ソルバも、小並列時にはPCGよりも速いであろうGAMGに変更しています。

system/fvSolution

```

    p_rgh
    {

        solver          GAMG;
        smoother        DIC;

```

```bash

buoyantBoussinesqPimpleFoam  | tee log.buoyantBoussinesqPimpleFoam

```

上記の実行方法で、以下のログとなりました。

ただし、OpenFOAM for macOS v2206( https://github.com/gerlero/openfoam-app )による実行結果です。

 log.buoyantBoussinesqPimpleFoam

```

(snip)

Courant Number mean: 0.0128869 max: 0.500519

deltaT = 0.000153058
Time = 0.00173489

PIMPLE: iteration 1
DILUPBiCGStab:  Solving for T, Initial residual = 1, Final residual = 9.89184e-07, No Iterations 47
GAMG:  Solving for p_rgh, Initial residual = 0.999957, Final residual = 0.001586, No Iterations 2
time step continuity errors : sum local = 1.62981e+08, global = 6.07348e-06, cumulative = 4.52581e-05
GAMG:  Solving for p_rgh, Initial residual = 0.944499, Final residual = 0.011611, No Iterations 1000
time step continuity errors : sum local = 2.5687e+06, global = 694349, cumulative = 694349
DILUPBiCGStab:  Solving for epsilon, Initial residual = 0.99835, Final residual = 2.38791e-07, No Iterations 4
bounding epsilon, min: -0.00170518 max: 3.54314e+14 average: 1.40995e+12
DILUPBiCGStab:  Solving for k, Initial residual = 9.89801e-06, Final residual = 1.39436e-06, No Iterations 1000
bounding k, min: -8.05256e+08 max: 1.14758e+12 average: 6.17646e+10
ExecutionTime = 7.63 s  ClockTime = 8 s

Courant Number mean: 1.55957e+08 max: 1.92803e+10
deltaT = 7.93856e-15
Time = 0.00173489

PIMPLE: iteration 1
DILUPBiCGStab:  Solving for T, Initial residual = 1, Final residual = 9.96341e-07, No Iterations 43
GAMG:  Solving for p_rgh, Initial residual = 1, Final residual = nan, No Iterations 1000
time step continuity errors : sum local = nan, global = nan, cumulative = nan
GAMG:  Solving for p_rgh, Initial residual = nan, Final residual = nan, No Iterations 1000


--> FOAM FATAL IO ERROR: (openfoam-2206)
Wrong token type - expected scalar value, found on line 0: word 'nan'

```

見事に発散していますが、kがかなり負になっています。

```

bounding k, min: -8.05256e+08 max: 1.14758e+12 average: 6.17646e+10

```

なお、k-*モデルを用いた乱流計算をしていて、Floating point exceptionのエラーが出る時は、

k、epsilonの値が負になっていることが原因である場合が多いです。

そこで、0/k をみると、以下のようにかなり小さい値となっているので、

これをとりあえず1e-8などとある程度余裕がある値に修正すると、とりあえず発散しなくなりました。

(

0/k

```

internalField   uniform 1.35e-13;

boundaryField
{
    adiab
    {
        type            kqRWallFunction;
        value           uniform 1.35e-13;

```

↓

```

internalField   uniform 1e-8;

boundaryField
{
    adiab
    {
        type            kqRWallFunction;
        value           $internalField;

// 他のパッチも同様に  valueは$internalFieldとしました

```

さらに、以下のfunctionをsystem/controlDictに追加して、

解析している場の範囲をログに出力しましたが、

温度Tがアンダ/オーバーシュートしていました。

(温度)

system/controlDict

```

functions
{
    fieldMinMax1
    {
        type            fieldMinMax;
        libs            (fieldFunctionObjects);
        writeToFile     no;
        log             yes;
        location        yes;
        mode            magnitude;
        fields          ( p_rgh U k epsilon T );
    }
}

```

なお、上記のfieldMinMaxのfunctionについては、以下のドキュメントや
Further informationのTutorialを参考にしてください。

https://www.openfoam.com/documentation/guides/latest/doc/guide-fos-field-fieldMinMax.html 

そこで、0/*の場の設定をみてみると、

backやfrontといったwall型のパッチでの速度Uの境界条件がzeroGradientとなっているので、

本来の境界条件は不明ですが、

とりあえずslipに変更したら、温度のアンダー/オーバーシュートは収まりました。

0/U

```

    back
    {
        type            zeroGradient;
    }
    front
    {
        type            zeroGradient;
    }

```

↓

```

    back
    {
        type            slip;
    }
    front
    {
        type            slip;
    }

```

他にも以下の設定が気になりましたが、とりあえず以上です。

system/controlDict

```

maxCo       1; 

```

(ソルバのログでCourant Numberのmaxが1を超えてしまわないか)

system/fvSolution

```

        momentumPredictor no;
```

system/fvSchemes

```

laplacianSchemes
{

        default         Gauss linear uncorrected;

(snip)

snGradSchemes
{
        default         uncorrected;

```

(checkMeshによる最大非直交性角度は約45度と小さいので、非直交補正をかけても不安定になりにくい)

     Mesh non-orthogonality Max: 45.3762 average: 9.54797

)

constant/transportProperties

```

// Reference temperature
TRef            300;

```

(温度範囲は300-370[K])

[長濱谷翔伍]

今野様

返答が遅くなり大変申し訳ございません。

ご教授頂いた通りに変更したところ、何とか解析を進めることが出来ました。

ご協力いただき本当にありがとうございました。

今後もよろしくお願いいたします。

2023年11月7日火曜日 12:59:40 UTC+9 Masashi Imano: