代数方程式ソルバーの比較

2013年8月30日

はじめに

OpenFOAM の代数方程式ソルバーの比較を行う。

使用バージョン

OpenFOAM 2.2.1

代数方程式ソルバーの比較

OpenFOAM の代数方程式ソルバーの比較を行う。

simpleFoam のチュートリアルケース pitzDaily を使用する。チュートリアルケースの p の relTol を 0.05 にし、k, epsilon の緩和係数を 0.5 に設定、renumberMesh を適用したものを比較対象とし、このケースにおける最適設定を探る。

設定

基本

    p
    {
        solver          PCG;
        preconditioner  DIC;
        tolerance       1e-06;
        relTol          0.05;
    }

    "(U|k|epsilon)"
    {
        solver          PBiCG;
        preconditioner  DILU;
        tolerance       1e-05;
        relTol          0.1;
    }

GAMG

    p
    {
        solver          GAMG;
        smoother        GaussSeidel;
        agglomerator    faceAreaPair;
        nCellsInCoarsestLevel 10;
        mergeLevels     1;
        nPreSweeps      2;
        directSolveCoarsest true;
        cacheAgglomeration on;
        tolerance       1e-06;
        relTol          0.05;
    }

PCG/GAMG preconditioner

    p
    {
        solver          PCG;
	
        preconditioner
        {
            preconditioner  GAMG;
            smoother        GaussSeidel;
            agglomerator    faceAreaPair;
            nCellsInCoarsestLevel 10;
            mergeLevels     1;
            nPreSweeps      2;
            directSolveCoarsest true;
            cacheAgglomeration on;
        }
        
        tolerance       1e-06;
        relTol          0.05;
    }

smoothSolver

    U
    {
        solver          smoothSolver;
        smoother        DILU;
        tolerance       1e-05;
        relTol          0.1;
    }

coupled

    U
    {
        type            coupled;
        solver          PBiCCCG;
        preconditioner  DILU;
        tolerance       (1e-05 1e-05 1e-05);
        relTol          (0.1 0.1 0.1);
    }

比較対象

ケース反復回数Execution Time [s}
p: PCG, U: PBiCG75141.81

計算マシンの CPU は Core i7 870 (4 コア, 定格 2.93 GHz)。1 CPU で計算。

p/GAMG

ケース反復回数Execution Time [s}
GussSeidel, faceAreaPair, nPreSweep 074731.26
algebraicPair74731.43
nPreSweeps 275030.4
directSolveCoarsest true74730.83
cacheAgglomeration on74728.69
nPreSweeps 2, directSolveCoarsest true, cacheAgglomeration on75027.63

p の GAMG の各設定の比較。このケースでは、nPreSweeps 2, directSolveCoarsest true, cacheAgglomeration on がもっとも計算時間が短かったので、これを採用する。

p/GAMG 2

ケース反復回数Execution Time [s}
GaussSeidel75027.63
DIC75027.41
DICGaussSeidel75029.95
FDIC75027.47
nonBlockingGuassSeidel75027.92
symGaussSeidel74731.06

smoother の比較。たいして変わらないが、DIC が多少速いらしい。

p/PCG

ケース反復回数Execution Time [s}
DIC75141.81
FDIC75142.8
GAMG (GaussSeidel)75028.66
diagonal75170.32
none753142.36

PCG の preconditioner の比較。GAMG が速いが、それ以外であれば DIC でよさそうである。

p/smoothSolver

ケース反復回数Execution Time [s}
GaussSeidel756180.23
DILU753245.74

smoothSolver の smoother の比較。かえって時間がかかっているため、中止。

U/GAMG

ケース反復回数Execution Time [s}
GaussSeidel75646.06
DILU76546.04

U の GAMG の smoother の比較。かえって時間がかかっているため、中止。

U/PBiCG

ケース反復回数Execution Time [s}
DILU75141.81
diagonal76436.69
none76342.83

PBiCG の preconditioner の比較。diagonal が速い。

U/smoothSolver

ケース反復回数Execution Time [s}
GaussSeidel78443.23
DILU76541.95
nonBlockingGaussSeidel78443.26
symGaussSeidel78743.42

smoothSolver の smoother の比較。DILU が多少速いが、PBiCG よりは遅い。

U/coupled

ケース反復回数Execution Time [s}
PBiCCCG53629.8
PBiCICG52629.34
SmoothSolver (smoother: GaussSeidel)55830.35
symGaussSeidel78743.42

coupled の solver の比較。PBiCICG が多少速い。

U/coupled 2

ケース反復回数Execution Time [s}
PBiCICG, DILU52629.34
diagonal53929.62

preconditioner の比較。DILU でよさそう。

k-epsilon/PCG

ケース反復回数Execution Time [s}
DILU75141.81
diagonal74539.9

k-epsilon の PCG の preconditioner の比較、U と傾向が合うかどうかの確認。やはり diagonal が多少速い。

k-epsilon/smoothSolver

ケース反復回数Execution Time [s}
GaussSeidel75141.81
DILU74540.45

smoothSolver の smoother の比較、U と傾向が合うかどうかの確認。やはり DILU が多少速い。

最適設定

以上の結果から、以下の設定がよさそうである (左から速そうな順番)。

  • p: GAMG (DIC), PCG (DIC)
  • U: coupled (PBiCICG, DILU), PCG (diagonal), smoothSolver(DILU)
  • k-epsilon: PCG (diagonal), smoothSolver (DILU)

ただし、組み合わせで傾向が変わる可能性もあるので、最後に組み合わせでいくつかのケースを比較してみる。

pUk-epsilon反復回数Execution Time [s}
PCG (GaussSeidel)PBiCG (DILU)PBiCG (DILU)75141.81
PCG (GaussSeidel)PBiCG (diagonal)PBiCG (diagonal)76239.33
PCG (GaussSeidel)smoothSolver (DILU)smoothSolver (DILU)76141.21
GAMG (DIC)PBiCG (diagonal)PBiCG (diagonal)76626.92
GAMG (DIC)PBiCG (diagonal)smoothSolver (DILU)75726.96
GAMG (DIC)smoothSolver (DILU)PBiCG (diagonal)76627.27
GAMG (DIC)smoothSolver (DILU)smoothSolver (DILU)76127.04
GAMG (DIC)coupled (PBiCICG, DILU)PBiCG (diagonal)53219.52
GAMG (DIC)coupled (PBiCICG, DILU)smoothSolver (DILU)52819.36
GAMG (GaussSeidel)PBiCG (diagonal)PBiCG (diagonal)76726.84
GAMG (GaussSeidel)PBiCG (diagonal)smoothSolver (DILU)75826.42
GAMG (GaussSeidel)smoothSolver (DILU)PBiCG (diagonal)76627.25
GAMG (GaussSeidel)smoothSolver (DILU)smoothSolver (DILU)76127.14
GAMG (GaussSeidel)coupled (PBiCICG, DILU)PBiCG (diagonal)53219.38
GAMG (GaussSeidel)coupled (PBiCICG, DILU)smoothSolver (DILU)52819.17

U, k-epsilon は smoothSolver よりも PCG (diagnoal) のほうがよさそうであったが、GAMG と組み合わせると差がなくなってしまっている。coupled との併用では、smoothSolver のほうがほんの少し速そうであるが、差は微妙である。GAMG の preconditioner について、GaussSeidel と DIC の差はほとんどなかったので、一応 GussSeidel も試してみると、ほんの少し DIC よりも速いときがあるが、言うほどの差があるわけでもないので、どちらでもよいのだろう。

GAMG 再検討

最適設定のめどがついてきたので、ここで改めて GAMG のパラメタを細かく検討してみる。

ケース反復回数Execution Time [s}
base (DIC)52819.36
nCellsInCoarsestLevel 10 -> 552819.32
nCellsInCoarsestLevel 10 -> 2052819.3
mergeLevels 1 -> 252819.19
nCellsInCoarsestLevel 20, mergeLevels 252819.02
nPreSweeps 2 -> 052817.3
nPostSweeps 2 -> 152821.63
nPreSweeps 0, nPostSweeps 152816.99
maxPre/PostSweeps 4 (default) -> 1052819.51
maxPre/PostSweeps 252818.59
maxPre/PostSweeps 10, pre/postSweepsLevelMultiplier 252820
nFinestSweep 2 (default) -> 152820.66
nFinestSweep 352819.32
directSolveCoarsest true -> false52819.35
scaleCorrection true (default) -> false52819.24
interpolateCorrection off (default) -> on52819.79
nPreSweeps 0, nPostSweeps 1, maxPre/PostSweeps 252816.99

nPreSweeps, nPostSweeps, maxPreSweeps, nPostSweeps を減らすのは効果がありそうである。nCellsInCoarsestLevel, mergeLeves は nPreSweeps などと組み合わせるとよい効果はみられない。directSolveCoarsest は今回はオフにしてもたいして変わらないようである。

smoother も再検討する。

ケース反復回数Execution Time [s}
DIC52816.99
GaussSeidel52821.73
DICGaussSeidel52817.77
FDIC52817.25
nonBlockingGuassSeidel52821.82
symGaussSeidel52822.33

GaussSeidel はあまり速くない。

上記の設定を使い、PCG の smoother に GAMG を用いた場合の nVcycles の違いをみる。

ケース反復回数Execution Time [s}
nVcycles 2 (default)52817.93
nVcycles 152818.81
nVcycles 352818.31

変えてもよくはならなそうである。

結論

今回のケースでは、GAMG や coupled が使えない状況も考慮して、つぎの 3 つを用意しておくのがよさそうである。

候補pUk-epsilon
1GAMG (DIC)coupled (PBiCICG, DILU)smoothSolver (DILU)
2GAMG (DIC)PBiCG (diagonal)smoothSolver (DILU)
3PCG (DIC)PBiCG (diagonal)PBiCG (diagonal)

差が微妙なものは、反復回数が少ないほうを選んだ。

候補 1

    p
    {
        solver          GAMG;
        smoother        DIC;
        agglomerator    faceAreaPair;
        nCellsInCoarsestLevel 10;
        mergeLevels     1;
        nPreSweeps      0;
        nPostSweeps     1;
        maxPreSweeps    2;
        maxPostSweeps   2;
        nFinestSweeps   2;
        preSweepsLevelMultiplier 1;
        postSweepsLevelMultiplier 1;
        scaleCorrection true;
        directSolveCoarsest false;
        cacheAgglomeration on;
        interpolateCorrection off;
        tolerance       1e-06;
        relTol          0.05;
    }

    U
    {
        type            coupled;
        solver          PBiCICG;
        preconditioner  DILU;
        tolerance       (1e-05 1e-05 1e-05);
        relTol          (0.1 0.1 0.1);
    }

    "(k|epsilon)"
    {
        solver          smoothSolver;
        smoother        DILU;
        tolerance       1e-05;
        relTol          0.1;
    }

候補 2

    p
    {
        solver           GAMG;
        smoother         DIC;
        agglomerator     faceAreaPair;
        nCellsInCoarsestLevel 10;
        mergeLevels      1;
        nPreSweeps       0;
        nPostSweeps      1;
        maxPreSweeps     2;
        maxPostSweeps    2;
        nFinestSweeps    2;
        preSweepsLevelMultiplier 1;
        postSweepsLevelMultiplier 1;
        scaleCorrection  true;
        directSolveCoarsest false;
        cacheAgglomeration on;
        interpolateCorrection off;
        tolerance        1e-06;
        relTol           0.05;
    }

    U
    {
        solver          PBiCG;
        preconditioner  diagonal;
        tolerance       1e-05;
        relTol          0.1;
    }

    "(k|epsilon)"
    {
        solver          smoothSolver;
        smoother        DILU;
        tolerance       1e-05;
        relTol          0.1;
    }

候補 3

    p
    {
        solver          PCG;
        preconditioner  DIC;
        tolerance       1e-06;
        relTol          0.05;
    }

    "(U|k|epsilon)"
    {
        solver          PBiCG;
        preconditioner  diagonal;
        tolerance       1e-05;
        relTol          0.1;
    }