PBiCGStab: New preconditioned bi-conjugate gradient stabilized solver for asymmetric lduMatrices

using a run-time selectable preconditioner References: Van der Vorst, H. A. (1992). Bi-CGSTAB: A fast and smoothly converging variant of Bi-CG for the solution of nonsymmetric linear systems. SIAM Journal on scientific and Statistical Computing, 13(2), 631-644. Barrett, R., Berry, M. W., Chan, T. F., Demmel, J., Donato, J., Dongarra, J., Eijkhout, V., Pozo, R., Romine, C. & Van der Vorst, H. (1994). Templates for the solution of linear systems: building blocks for iterative methods (Vol. 43). Siam. See also: https://en.wikipedia.org/wiki/Biconjugate_gradient_stabilized_method Tests have shown that PBiCGStab with the DILU preconditioner is more robust, reliable and shows faster convergence (~2x) than PBiCG with DILU, in particular in parallel where PBiCG occasionally diverges. This remarkable improvement over PBiCG prompted the update of all tutorial cases currently using PBiCG to use PBiCGStab instead. If any issues arise with this update please report on Mantis: http://bugs.openfoam.org
2025-12-08 06:57:46 +00:00 · 2016-09-05 11:46:42 +01:00
parent fec5c063f5
commit c339d3018c
58 changed files with 481 additions and 102 deletions
--- a/src/OpenFOAM/Make/files
+++ b/src/OpenFOAM/Make/files
@ -271,6 +271,7 @@ $(lduMatrix)/solvers/diagonalSolver/diagonalSolver.C
 $(lduMatrix)/solvers/smoothSolver/smoothSolver.C
 $(lduMatrix)/solvers/PCG/PCG.C
 $(lduMatrix)/solvers/PBiCG/PBiCG.C
+$(lduMatrix)/solvers/PBiCGStab/PBiCGStab.C

 $(lduMatrix)/smoothers/GaussSeidel/GaussSeidelSmoother.C
 $(lduMatrix)/smoothers/symGaussSeidel/symGaussSeidelSmoother.C
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCG/PBiCG.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCG/PBiCG.C
@ -76,37 +76,25 @@ Foam::solverPerformance Foam::PBiCG::solve
        fieldName_
    );

-    label nCells = psi.size();
+    const label nCells = psi.size();

    scalar* __restrict__ psiPtr = psi.begin();

    scalarField pA(nCells);
    scalar* __restrict__ pAPtr = pA.begin();

-    scalarField pT(nCells, 0.0);
-    scalar* __restrict__ pTPtr = pT.begin();
-
    scalarField wA(nCells);
    scalar* __restrict__ wAPtr = wA.begin();

-    scalarField wT(nCells);
-    scalar* __restrict__ wTPtr = wT.begin();
-
-    scalar wArT = solverPerf.great_;
-    scalar wArTold = wArT;
-
-    // --- Calculate A.psi and T.psi
+    // --- Calculate A.psi
    matrix_.Amul(wA, psi, interfaceBouCoeffs_, interfaces_, cmpt);
-    matrix_.Tmul(wT, psi, interfaceIntCoeffs_, interfaces_, cmpt);

-    // --- Calculate initial residual and transpose residual fields
+    // --- Calculate initial residual field
    scalarField rA(source - wA);
-    scalarField rT(source - wT);
    scalar* __restrict__ rAPtr = rA.begin();
-    scalar* __restrict__ rTPtr = rT.begin();

    // --- Calculate normalisation factor
-    scalar normFactor = this->normFactor(psi, source, wA, pA);
+    const scalar normFactor = this->normFactor(psi, source, wA, pA);

    if (lduMatrix::debug >= 2)
    {
@ -126,6 +114,22 @@ Foam::solverPerformance Foam::PBiCG::solve
     || !solverPerf.checkConvergence(tolerance_, relTol_)
    )
    {
+        scalarField pT(nCells, 0);
+        scalar* __restrict__ pTPtr = pT.begin();
+
+        scalarField wT(nCells);
+        scalar* __restrict__ wTPtr = wT.begin();
+
+        // --- Calculate T.psi
+        matrix_.Tmul(wT, psi, interfaceIntCoeffs_, interfaces_, cmpt);
+
+        // --- Calculate initial transpose residual field
+        scalarField rT(source - wT);
+        scalar* __restrict__ rTPtr = rT.begin();
+
+        // --- Initial value not used
+        scalar wArT = 0;
+
        // --- Select and construct the preconditioner
        autoPtr<lduMatrix::preconditioner> preconPtr =
        lduMatrix::preconditioner::New
@ -138,7 +142,7 @@ Foam::solverPerformance Foam::PBiCG::solve
        do
        {
            // --- Store previous wArT
-            wArTold = wArT;
+            const scalar wArTold = wArT;

            // --- Precondition residuals
            preconPtr->precondition(wA, rA, cmpt);
@ -157,7 +161,7 @@ Foam::solverPerformance Foam::PBiCG::solve
            }
            else
            {
-                scalar beta = wArT/wArTold;
+                const scalar beta = wArT/wArTold;

                for (label cell=0; cell<nCells; cell++)
                {
@ -171,7 +175,7 @@ Foam::solverPerformance Foam::PBiCG::solve
            matrix_.Amul(wA, pA, interfaceBouCoeffs_, interfaces_, cmpt);
            matrix_.Tmul(wT, pT, interfaceIntCoeffs_, interfaces_, cmpt);

-            scalar wApT = gSumProd(wA, pT, matrix().mesh().comm());
+            const scalar wApT = gSumProd(wA, pT, matrix().mesh().comm());

            // --- Test for singularity
            if (solverPerf.checkSingularity(mag(wApT)/normFactor))
@ -182,7 +186,7 @@ Foam::solverPerformance Foam::PBiCG::solve

            // --- Update solution and residual:

-            scalar alpha = wArT/wApT;
+            const scalar alpha = wArT/wApT;

            for (label cell=0; cell<nCells; cell++)
            {
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCGStab/PBiCGStab.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCGStab/PBiCGStab.C
@ -0,0 +1,252 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "PBiCGStab.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(PBiCGStab, 0);
+
+    lduMatrix::solver::addasymMatrixConstructorToTable<PBiCGStab>
+        addPBiCGStabAsymMatrixConstructorToTable_;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::PBiCGStab::PBiCGStab
+(
+    const word& fieldName,
+    const lduMatrix& matrix,
+    const FieldField<Field, scalar>& interfaceBouCoeffs,
+    const FieldField<Field, scalar>& interfaceIntCoeffs,
+    const lduInterfaceFieldPtrsList& interfaces,
+    const dictionary& solverControls
+)
+:
+    lduMatrix::solver
+    (
+        fieldName,
+        matrix,
+        interfaceBouCoeffs,
+        interfaceIntCoeffs,
+        interfaces,
+        solverControls
+    )
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+Foam::solverPerformance Foam::PBiCGStab::solve
+(
+    scalarField& psi,
+    const scalarField& source,
+    const direction cmpt
+) const
+{
+    // --- Setup class containing solver performance data
+    solverPerformance solverPerf
+    (
+        lduMatrix::preconditioner::getName(controlDict_) + typeName,
+        fieldName_
+    );
+
+    const label nCells = psi.size();
+
+    scalar* __restrict__ psiPtr = psi.begin();
+
+    scalarField pA(nCells);
+    scalar* __restrict__ pAPtr = pA.begin();
+
+    scalarField yA(nCells);
+    scalar* __restrict__ yAPtr = yA.begin();
+
+    // --- Calculate A.psi
+    matrix_.Amul(yA, psi, interfaceBouCoeffs_, interfaces_, cmpt);
+
+    // --- Calculate initial residual field
+    scalarField rA(source - yA);
+    scalar* __restrict__ rAPtr = rA.begin();
+
+    // --- Calculate normalisation factor
+    const scalar normFactor = this->normFactor(psi, source, yA, pA);
+
+    if (lduMatrix::debug >= 2)
+    {
+        Info<< "   Normalisation factor = " << normFactor << endl;
+    }
+
+    // --- Calculate normalised residual norm
+    solverPerf.initialResidual() =
+        gSumMag(rA, matrix().mesh().comm())
+       /normFactor;
+    solverPerf.finalResidual() = solverPerf.initialResidual();
+
+    // --- Check convergence, solve if not converged
+    if
+    (
+        minIter_ > 0
+     || !solverPerf.checkConvergence(tolerance_, relTol_)
+    )
+    {
+        scalarField AyA(nCells);
+        scalar* __restrict__ AyAPtr = AyA.begin();
+
+        scalarField sA(nCells);
+        scalar* __restrict__ sAPtr = sA.begin();
+
+        scalarField zA(nCells);
+        scalar* __restrict__ zAPtr = zA.begin();
+
+        scalarField tA(nCells);
+        scalar* __restrict__ tAPtr = tA.begin();
+
+        // --- Store initial residual
+        const scalarField rA0(rA);
+
+        // --- Initial values not used
+        scalar rA0rA = 0;
+        scalar alpha = 0;
+        scalar omega = 0;
+
+        // --- Select and construct the preconditioner
+        autoPtr<lduMatrix::preconditioner> preconPtr =
+        lduMatrix::preconditioner::New
+        (
+            *this,
+            controlDict_
+        );
+
+        // --- Solver iteration
+        do
+        {
+            // --- Store previous rA0rA
+            const scalar rA0rAold = rA0rA;
+
+            rA0rA = gSumProd(rA0, rA, matrix().mesh().comm());
+
+            // --- Test for singularity
+            if (solverPerf.checkSingularity(mag(rA0rA)))
+            {
+                break;
+            }
+
+            // --- Update pA
+            if (solverPerf.nIterations() == 0)
+            {
+                for (label cell=0; cell<nCells; cell++)
+                {
+                    pAPtr[cell] = rAPtr[cell];
+                }
+            }
+            else
+            {
+                // --- Test for singularity
+                if (solverPerf.checkSingularity(mag(omega)))
+                {
+                    break;
+                }
+
+                const scalar beta = (rA0rA/rA0rAold)*(alpha/omega);
+
+                for (label cell=0; cell<nCells; cell++)
+                {
+                    pAPtr[cell] =
+                        rAPtr[cell] + beta*(pAPtr[cell] - omega*AyAPtr[cell]);
+                }
+            }
+
+            // --- Precondition pA
+            preconPtr->precondition(yA, pA, cmpt);
+
+            // --- Calculate AyA
+            matrix_.Amul(AyA, yA, interfaceBouCoeffs_, interfaces_, cmpt);
+
+            const scalar rA0AyA = gSumProd(rA0, AyA, matrix().mesh().comm());
+
+            alpha = rA0rA/rA0AyA;
+
+            // --- Calculate sA
+            for (label cell=0; cell<nCells; cell++)
+            {
+                sAPtr[cell] = rAPtr[cell] - alpha*AyAPtr[cell];
+            }
+
+            // --- Test sA for convergence
+            solverPerf.finalResidual() =
+                gSumMag(sA, matrix().mesh().comm())/normFactor;
+
+            if (solverPerf.checkConvergence(tolerance_, relTol_))
+            {
+                for (label cell=0; cell<nCells; cell++)
+                {
+                    psiPtr[cell] += alpha*yAPtr[cell];
+                }
+
+                solverPerf.nIterations()++;
+
+                return solverPerf;
+            }
+
+            // --- Precondition sA
+            preconPtr->precondition(zA, sA, cmpt);
+
+            // --- Calculate tA
+            matrix_.Amul(tA, zA, interfaceBouCoeffs_, interfaces_, cmpt);
+
+            const scalar tAtA = gSumSqr(tA, matrix().mesh().comm());
+
+            // --- Calculate omega from tA and sA
+            //     (cheaper than using zA with preconditioned tA)
+            omega = gSumProd(tA, sA)/tAtA;
+
+            // --- Update solution and residual
+            for (label cell=0; cell<nCells; cell++)
+            {
+                psiPtr[cell] += alpha*yAPtr[cell] + omega*zAPtr[cell];
+                rAPtr[cell] = sAPtr[cell] - omega*tAPtr[cell];
+            }
+
+            solverPerf.finalResidual() =
+                gSumMag(rA, matrix().mesh().comm())
+               /normFactor;
+        } while
+        (
+            (
+                solverPerf.nIterations()++ < maxIter_
+            && !solverPerf.checkConvergence(tolerance_, relTol_)
+            )
+         || solverPerf.nIterations() < minIter_
+        );
+    }
+
+    return solverPerf;
+}
+
+
+// ************************************************************************* //
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCGStab/PBiCGStab.H
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCGStab/PBiCGStab.H
@ -0,0 +1,123 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::PBiCGStab
+
+Description
+    Preconditioned bi-conjugate gradient stabilized solver for asymmetric
+    lduMatrices using a run-time selectable preconditiioner.
+
+    References:
+    \verbatim
+        Van der Vorst, H. A. (1992).
+        Bi-CGSTAB: A fast and smoothly converging variant of Bi-CG
+        for the solution of nonsymmetric linear systems.
+        SIAM Journal on scientific and Statistical Computing, 13(2), 631-644.
+
+        Barrett, R., Berry, M. W., Chan, T. F., Demmel, J., Donato, J.,
+        Dongarra, J., Eijkhout, V., Pozo, R., Romine, C. & Van der Vorst, H.
+        (1994).
+        Templates for the solution of linear systems:
+        building blocks for iterative methods
+        (Vol. 43). Siam.
+    \endverbatim
+
+SourceFiles
+    PBiCGStab.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef PBiCGStab_H
+#define PBiCGStab_H
+
+#include "lduMatrix.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                           Class PBiCGStab Declaration
+\*---------------------------------------------------------------------------*/
+
+class PBiCGStab
+:
+    public lduMatrix::solver
+{
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct
+        PBiCGStab(const PBiCGStab&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const PBiCGStab&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("PBiCGStab");
+
+
+    // Constructors
+
+        //- Construct from matrix components and solver data stream
+        PBiCGStab
+        (
+            const word& fieldName,
+            const lduMatrix& matrix,
+            const FieldField<Field, scalar>& interfaceBouCoeffs,
+            const FieldField<Field, scalar>& interfaceIntCoeffs,
+            const lduInterfaceFieldPtrsList& interfaces,
+            const dictionary& solverControls
+        );
+
+
+    //- Destructor
+    virtual ~PBiCGStab()
+    {}
+
+
+    // Member Functions
+
+        //- Solve the matrix with this solver
+        virtual solverPerformance solve
+        (
+            scalarField& psi,
+            const scalarField& source,
+            const direction cmpt=0
+        ) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
--- a/tutorials/basic/scalarTransportFoam/pitzDaily/system/fvSolution
+++ b/tutorials/basic/scalarTransportFoam/pitzDaily/system/fvSolution
@ -19,7 +19,7 @@ solvers
 {
    T
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-06;
        relTol          0;
--- a/tutorials/combustion/PDRFoam/flamePropagationWithObstacles/system/fvSolution
+++ b/tutorials/combustion/PDRFoam/flamePropagationWithObstacles/system/fvSolution
@ -49,7 +49,7 @@ solvers

    "(b|Xi|ft|ha|hau|k|epsilon)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-06;
        relTol          0.1;
@ -57,7 +57,7 @@ solvers

    "(b|Xi|ft|ha|hau|k|epsilon)Final"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-06;
        relTol          0;
--- a/tutorials/combustion/XiFoam/ras/moriyoshiHomogeneous/system/fvSolution
+++ b/tutorials/combustion/XiFoam/ras/moriyoshiHomogeneous/system/fvSolution
@ -34,7 +34,7 @@ solvers

    "(U|b|Su|Xi|ha|hau|k|epsilon)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-05;
        relTol          0.1;
@ -42,7 +42,7 @@ solvers

    "(U|b|Su|Xi|ha|hau|k|epsilon)Final"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-05;
        relTol          0;
--- a/tutorials/combustion/chemFoam/gri/system/fvSolution
+++ b/tutorials/combustion/chemFoam/gri/system/fvSolution
@ -19,7 +19,7 @@ solvers
 {
    Yi
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-12;
        relTol          0;
--- a/tutorials/combustion/chemFoam/h2/system/fvSolution
+++ b/tutorials/combustion/chemFoam/h2/system/fvSolution
@ -19,7 +19,7 @@ solvers
 {
    Yi
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-12;
        relTol          0;
--- a/tutorials/combustion/chemFoam/ic8h18/system/fvSolution
+++ b/tutorials/combustion/chemFoam/ic8h18/system/fvSolution
@ -19,7 +19,7 @@ solvers
 {
    Yi
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-12;
        relTol          0;
--- a/tutorials/combustion/chemFoam/ic8h18_TDAC/system/fvSolution
+++ b/tutorials/combustion/chemFoam/ic8h18_TDAC/system/fvSolution
@ -19,7 +19,7 @@ solvers
 {
    Yi
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-12;
        relTol          0;
--- a/tutorials/combustion/chemFoam/nc7h16/system/fvSolution
+++ b/tutorials/combustion/chemFoam/nc7h16/system/fvSolution
@ -19,7 +19,7 @@ solvers
 {
    Yi
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-12;
        relTol          0;
--- a/tutorials/combustion/engineFoam/kivaTest/system/fvSolution
+++ b/tutorials/combustion/engineFoam/kivaTest/system/fvSolution
@ -49,7 +49,7 @@ solvers

    "(U|Xi|hau|eau|ft|b|ha|ea|k|epsilon)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-05;
        relTol          0.1;
--- a/tutorials/combustion/fireFoam/les/flameSpreadWaterSuppressionPanel/system/filmRegion/fvSolution
+++ b/tutorials/combustion/fireFoam/les/flameSpreadWaterSuppressionPanel/system/filmRegion/fvSolution
@ -19,21 +19,21 @@ solvers
 {
    hf
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       0;
        relTol          1e-3;
    }
    "(Uf|deltaf\*rhof)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-10;
        relTol          0;
    }
    deltaf
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-10;
        relTol          0;
--- a/tutorials/combustion/fireFoam/les/flameSpreadWaterSuppressionPanel/system/fvSolution
+++ b/tutorials/combustion/fireFoam/les/flameSpreadWaterSuppressionPanel/system/fvSolution
@ -57,7 +57,7 @@ solvers

    "(U|Yi|h|k)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-8;
        relTol          0.1;
--- a/tutorials/combustion/fireFoam/les/flameSpreadWaterSuppressionPanel/system/pyrolysisRegion/fvSolution
+++ b/tutorials/combustion/fireFoam/les/flameSpreadWaterSuppressionPanel/system/pyrolysisRegion/fvSolution
@ -26,7 +26,7 @@ solvers

    "Yi"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-06;
        relTol          0;
--- a/tutorials/combustion/fireFoam/les/oppositeBurningPanels/system/fvSolution
+++ b/tutorials/combustion/fireFoam/les/oppositeBurningPanels/system/fvSolution
@ -47,7 +47,7 @@ solvers

    "(U|Yi|k|h|omega)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0.1;
--- a/tutorials/combustion/fireFoam/les/oppositeBurningPanels/system/panelRegion/fvSolution
+++ b/tutorials/combustion/fireFoam/les/oppositeBurningPanels/system/panelRegion/fvSolution
@ -26,7 +26,7 @@ solvers

    "Yi"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-06;
        relTol          0;
--- a/tutorials/combustion/fireFoam/les/smallPoolFire2D/system/fvSolution
+++ b/tutorials/combustion/fireFoam/les/smallPoolFire2D/system/fvSolution
@ -56,7 +56,7 @@ solvers

    "(U|Yi|k|h|omega)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0.1;
--- a/tutorials/combustion/fireFoam/les/smallPoolFire3D/system/fvSolution
+++ b/tutorials/combustion/fireFoam/les/smallPoolFire3D/system/fvSolution
@ -48,7 +48,7 @@ solvers

    "(U|Yi|k|h)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0.1;
--- a/tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D/system/fvSolution
+++ b/tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D/system/fvSolution
@ -39,7 +39,7 @@ solvers

    "(U|h|k|epsilon)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0.1;
--- a/tutorials/combustion/reactingFoam/laminar/counterFlowFlame2DLTS/system/fvSolution
+++ b/tutorials/combustion/reactingFoam/laminar/counterFlowFlame2DLTS/system/fvSolution
@ -38,7 +38,7 @@ solvers

    "(U|h|k|epsilon)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0.1;
--- a/tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/system/fvSolution
+++ b/tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/system/fvSolution
@ -39,7 +39,7 @@ solvers

    "(U|h|k|epsilon)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0.1;
--- a/tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/system/fvSolution
+++ b/tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/system/fvSolution
@ -39,7 +39,7 @@ solvers

    "(U|h|k|epsilon)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0.1;
--- a/tutorials/compressible/rhoCentralDyMFoam/movingCone/system/fvSolution
+++ b/tutorials/compressible/rhoCentralDyMFoam/movingCone/system/fvSolution
@ -19,7 +19,7 @@ solvers
 {
    p
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-06;
        relTol          0.01;
--- a/tutorials/compressible/sonicDyMFoam/movingCone/system/fvSolution
+++ b/tutorials/compressible/sonicDyMFoam/movingCone/system/fvSolution
@ -19,7 +19,7 @@ solvers
 {
    p
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-06;
        relTol          0.01;
--- a/tutorials/compressible/sonicFoam/ras/nacaAirfoil/system/fvSolution
+++ b/tutorials/compressible/sonicFoam/ras/nacaAirfoil/system/fvSolution
@ -24,7 +24,7 @@ solvers

    "p.*"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-12;
        relTol          0;
--- a/tutorials/financial/financialFoam/europeanCall/system/fvSolution
+++ b/tutorials/financial/financialFoam/europeanCall/system/fvSolution
@ -19,7 +19,7 @@ solvers
 {
    V
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-05;
        relTol          0;
--- a/tutorials/heatTransfer/buoyantBoussinesqPimpleFoam/hotRoom/system/fvSolution
+++ b/tutorials/heatTransfer/buoyantBoussinesqPimpleFoam/hotRoom/system/fvSolution
@ -33,7 +33,7 @@ solvers

    "(U|T|k|epsilon|R)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0.1;
--- a/tutorials/heatTransfer/buoyantBoussinesqSimpleFoam/hotRoom/system/fvSolution
+++ b/tutorials/heatTransfer/buoyantBoussinesqSimpleFoam/hotRoom/system/fvSolution
@ -27,7 +27,7 @@ solvers

    "(U|T|k|epsilon|R)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-05;
        relTol          0.1;
--- a/tutorials/heatTransfer/buoyantBoussinesqSimpleFoam/iglooWithFridges/system/fvSolution
+++ b/tutorials/heatTransfer/buoyantBoussinesqSimpleFoam/iglooWithFridges/system/fvSolution
@ -27,7 +27,7 @@ solvers

    "(U|T|k|epsilon)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-07;
        relTol          0.1;
--- a/tutorials/heatTransfer/buoyantPimpleFoam/hotRoom/system/fvSolution
+++ b/tutorials/heatTransfer/buoyantPimpleFoam/hotRoom/system/fvSolution
@ -41,7 +41,7 @@ solvers

    "(U|h|e|k|epsilon|R)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0.1;
--- a/tutorials/heatTransfer/buoyantSimpleFoam/buoyantCavity/system/fvSolution
+++ b/tutorials/heatTransfer/buoyantSimpleFoam/buoyantCavity/system/fvSolution
@ -29,7 +29,7 @@ solvers

    "(U|h|k|epsilon|omega)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-8;
        relTol          0.1;
--- a/tutorials/heatTransfer/buoyantSimpleFoam/circuitBoardCooling/system/fvSolution
+++ b/tutorials/heatTransfer/buoyantSimpleFoam/circuitBoardCooling/system/fvSolution
@ -29,7 +29,7 @@ solvers

    "(U|h|k|epsilon|omega)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-7;
        relTol          0.01;
--- a/tutorials/heatTransfer/buoyantSimpleFoam/externalCoupledCavity/system/fvSolution
+++ b/tutorials/heatTransfer/buoyantSimpleFoam/externalCoupledCavity/system/fvSolution
@ -29,7 +29,7 @@ solvers

    "(U|h|k|epsilon|omega)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-8;
        relTol          0.1;
--- a/tutorials/heatTransfer/buoyantSimpleFoam/hotRadiationRoom/system/fvSolution
+++ b/tutorials/heatTransfer/buoyantSimpleFoam/hotRadiationRoom/system/fvSolution
@ -27,7 +27,7 @@ solvers

    "(U|h|k|epsilon)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-05;
        relTol          0.1;
--- a/tutorials/heatTransfer/buoyantSimpleFoam/hotRadiationRoomFvDOM/system/fvSolution
+++ b/tutorials/heatTransfer/buoyantSimpleFoam/hotRadiationRoomFvDOM/system/fvSolution
@ -37,7 +37,7 @@ solvers

    "(U|h|k|epsilon)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-05;
        relTol          0.1;
--- a/tutorials/heatTransfer/chtMultiRegionFoam/multiRegionHeater/system/bottomWater/fvSolution
+++ b/tutorials/heatTransfer/chtMultiRegionFoam/multiRegionHeater/system/bottomWater/fvSolution
@ -50,7 +50,7 @@ solvers

    "(U|h|k|epsilon|R)"
    {
-        solver           PBiCG;
+        solver           PBiCGStab;
        preconditioner   DILU;
        tolerance        1e-7;
        relTol           0.1;
--- a/tutorials/heatTransfer/chtMultiRegionFoam/multiRegionHeater/system/topAir/fvSolution
+++ b/tutorials/heatTransfer/chtMultiRegionFoam/multiRegionHeater/system/topAir/fvSolution
@ -50,7 +50,7 @@ solvers

    "(U|h|k|epsilon|R)"
    {
-        solver           PBiCG;
+        solver           PBiCGStab;
        preconditioner   DILU;
        tolerance        1e-7;
        relTol           0.1;
--- a/tutorials/heatTransfer/chtMultiRegionFoam/snappyMultiRegionHeater/system/bottomAir/fvSolution
+++ b/tutorials/heatTransfer/chtMultiRegionFoam/snappyMultiRegionHeater/system/bottomAir/fvSolution
@ -43,7 +43,7 @@ solvers

    "(U|h|k|epsilon|R)"
    {
-        solver           PBiCG;
+        solver           PBiCGStab;
        preconditioner   DILU;
        tolerance        1e-7;
        relTol           0.1;
--- a/tutorials/heatTransfer/chtMultiRegionFoam/snappyMultiRegionHeater/system/topAir/fvSolution
+++ b/tutorials/heatTransfer/chtMultiRegionFoam/snappyMultiRegionHeater/system/topAir/fvSolution
@ -43,7 +43,7 @@ solvers

    "(U|h|k|epsilon|R)"
    {
-        solver           PBiCG;
+        solver           PBiCGStab;
        preconditioner   DILU;
        tolerance        1e-7;
        relTol           0.1;
--- a/tutorials/heatTransfer/chtMultiRegionSimpleFoam/heatExchanger/system/air/fvSolution
+++ b/tutorials/heatTransfer/chtMultiRegionSimpleFoam/heatExchanger/system/air/fvSolution
@ -31,7 +31,7 @@ solvers

    "(U|h|e|k|epsilon)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0.1;
--- a/tutorials/heatTransfer/chtMultiRegionSimpleFoam/heatExchanger/system/porous/fvSolution
+++ b/tutorials/heatTransfer/chtMultiRegionSimpleFoam/heatExchanger/system/porous/fvSolution
@ -30,7 +30,7 @@ solvers

    "(U|h|e|k|epsilon)"
    {
-        solver           PBiCG;
+        solver           PBiCGStab;
        preconditioner   DILU;
        tolerance        1e-7;
        relTol           0.1;
--- a/tutorials/heatTransfer/chtMultiRegionSimpleFoam/multiRegionHeaterRadiation/system/bottomAir/fvSolution
+++ b/tutorials/heatTransfer/chtMultiRegionSimpleFoam/multiRegionHeaterRadiation/system/bottomAir/fvSolution
@ -36,7 +36,7 @@ solvers

    "(U|h|k|epsilon|G|Ii)"
    {
-        solver           PBiCG;
+        solver           PBiCGStab;
        preconditioner   DILU;
        tolerance        1e-7;
        relTol           0.1;
--- a/tutorials/heatTransfer/chtMultiRegionSimpleFoam/multiRegionHeaterRadiation/system/topAir/fvSolution
+++ b/tutorials/heatTransfer/chtMultiRegionSimpleFoam/multiRegionHeaterRadiation/system/topAir/fvSolution
@ -36,7 +36,7 @@ solvers

    "(U|h|k|epsilon|G|Ii)"
    {
-        solver           PBiCG;
+        solver           PBiCGStab;
        preconditioner   DILU;
        tolerance        1e-7;
        relTol           0.1;
--- a/tutorials/incompressible/simpleFoam/motorBike/system/fvSolution
+++ b/tutorials/incompressible/simpleFoam/motorBike/system/fvSolution
@ -18,10 +18,10 @@ solvers
 {
    p
    {
-        solver           GAMG;
-        tolerance        1e-7;
-        relTol           0.01;
-        smoother         GaussSeidel;
+        solver          GAMG;
+        smoother        GaussSeidel;
+        tolerance       1e-7;
+        relTol          0.01;
    }

    Phi
@ -31,29 +31,29 @@ solvers

    U
    {
-        solver           smoothSolver;
-        smoother         GaussSeidel;
-        tolerance        1e-8;
-        relTol           0.1;
-        nSweeps          1;
+        solver          smoothSolver;
+        smoother        GaussSeidel;
+        tolerance       1e-8;
+        relTol          0.1;
+        nSweeps         1;
    }

    k
    {
-        solver           smoothSolver;
-        smoother         GaussSeidel;
-        tolerance        1e-8;
-        relTol           0.1;
-        nSweeps          1;
+        solver          smoothSolver;
+        smoother        GaussSeidel;
+        tolerance       1e-8;
+        relTol          0.1;
+        nSweeps         1;
    }

    omega
    {
-        solver           smoothSolver;
-        smoother         GaussSeidel;
-        tolerance        1e-8;
-        relTol           0.1;
-        nSweeps          1;
+        solver          smoothSolver;
+        smoother        GaussSeidel;
+        tolerance       1e-8;
+        relTol          0.1;
+        nSweeps         1;
    }
 }

--- a/tutorials/lagrangian/reactingParcelFilmFoam/cylinder/system/fvSolution
+++ b/tutorials/lagrangian/reactingParcelFilmFoam/cylinder/system/fvSolution
@ -70,14 +70,14 @@ solvers

    "(h|Yi|O2|N2|H2O)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0;
    }
    hFinal
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0;
--- a/tutorials/lagrangian/reactingParcelFilmFoam/hotBoxes/system/fvSolution
+++ b/tutorials/lagrangian/reactingParcelFilmFoam/hotBoxes/system/fvSolution
@ -34,7 +34,7 @@ solvers

    "(U|h|k|epsilon)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-5;
        relTol          0.1;
--- a/tutorials/lagrangian/reactingParcelFilmFoam/splashPanel/system/fvSolution
+++ b/tutorials/lagrangian/reactingParcelFilmFoam/splashPanel/system/fvSolution
@ -66,14 +66,14 @@ solvers

    "(h|Yi|O2|N2|H2O)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0;
    }
    hFinal
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0;
--- a/tutorials/lagrangian/reactingParcelFilmFoam/splashPanel/system/wallFilmRegion/fvSolution
+++ b/tutorials/lagrangian/reactingParcelFilmFoam/splashPanel/system/wallFilmRegion/fvSolution
@ -19,14 +19,14 @@ solvers
 {
    "(Uf|hf|deltaf\*rhof)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-10;
        relTol          0;
    }
    deltaf
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-10;
        relTol          0;
--- a/tutorials/lagrangian/reactingParcelFoam/filter/system/fvSolution
+++ b/tutorials/lagrangian/reactingParcelFoam/filter/system/fvSolution
@ -64,7 +64,7 @@ solvers

    "(Yi|O2|N2|H2O)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0;
--- a/tutorials/lagrangian/reactingParcelFoam/parcelInBox/system/fvSolution
+++ b/tutorials/lagrangian/reactingParcelFoam/parcelInBox/system/fvSolution
@ -65,7 +65,7 @@ solvers

    "(Yi|O2|N2|H2O)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0;
--- a/tutorials/lagrangian/reactingParcelFoam/verticalChannel/system/fvSolution
+++ b/tutorials/lagrangian/reactingParcelFoam/verticalChannel/system/fvSolution
@ -69,7 +69,7 @@ solvers

    "(Yi|O2|N2|H2O)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0;
--- a/tutorials/lagrangian/sprayFoam/aachenBomb/system/fvSolution
+++ b/tutorials/lagrangian/sprayFoam/aachenBomb/system/fvSolution
@ -64,7 +64,7 @@ solvers

    "(Yi|O2|N2|H2O)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0;
--- a/tutorials/mesh/moveDynamicMesh/SnakeRiverCanyon/system/fvSolution
+++ b/tutorials/mesh/moveDynamicMesh/SnakeRiverCanyon/system/fvSolution
@ -27,7 +27,7 @@ solvers

    U
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-05;
        relTol          0;
--- a/tutorials/multiphase/driftFluxFoam/ras/dahl/system/fvSolution
+++ b/tutorials/multiphase/driftFluxFoam/ras/dahl/system/fvSolution
@ -59,9 +59,8 @@ solvers

    "(U|k|epsilon)"
    {
-        solver          PBiCG;
+        solver          PBiCGStab;
        preconditioner  DILU;
-        smoother        symGaussSeidel;
        tolerance       1e-7;
        relTol          0.1;
        minIter         1;
--- a/tutorials/multiphase/driftFluxFoam/ras/mixerVessel2D/system/fvSolution
+++ b/tutorials/multiphase/driftFluxFoam/ras/mixerVessel2D/system/fvSolution
@ -67,8 +67,8 @@ solvers

    "(U|k|epsilon)"
    {
-        solver          smoothSolver;
-        smoother        symGaussSeidel;
+        solver          PBiCGStab;
+        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0.1;
    }
--- a/tutorials/multiphase/driftFluxFoam/ras/tank3D/system/fvSolution
+++ b/tutorials/multiphase/driftFluxFoam/ras/tank3D/system/fvSolution
@ -59,8 +59,8 @@ solvers

    "(U|k|epsilon)"
    {
-        solver          smoothSolver;
-        smoother        symGaussSeidel;
+        solver          PBiCGStab;
+        preconditioner  DILU;
        tolerance       1e-7;
        relTol          0.1;
        minIter         1;