Added BlockILU smoother solvers

src/foam/Make/files

@@ -756,6 +756,7 @@ $(BlockLduSolvers)/blockVectorNSolvers.C
 $(BlockLduSolvers)/BlockLduSolver/blockLduSolvers.C
 $(BlockLduSolvers)/BlockDiagonal/blockDiagonalSolvers.C
 $(BlockLduSolvers)/BlockGaussSeidel/blockGaussSeidelSolvers.C
+$(BlockLduSolvers)/BlockILU/blockILUSolvers.C
 $(BlockLduSolvers)/BlockCG/blockCGSolvers.C
 $(BlockLduSolvers)/BlockBiCGStab/blockBiCGStabSolvers.C
 $(BlockLduSolvers)/BlockGMRES/blockGMRESSolvers.C

src/foam/matrices/blockLduMatrix/BlockLduPrecons/BlockCholeskyPrecon/BlockCholeskyPrecon.C

@@ -709,10 +709,16 @@ void Foam::BlockCholeskyPrecon<Type>::diagInterfaceMultiply
         dDiag[fc[coeffI]] +=
             mult.tripleProduct
             (
-                intCoeffs[coeffI],
+                bouCoeffs[coeffI],
                 dDiag[fc[coeffI]],
-                bouCoeffs[coeffI]
+                intCoeffs[coeffI]
             );
+//            mult.tripleProduct
+//            (
+//                intCoeffs[coeffI],
+//                dDiag[fc[coeffI]],
+//                bouCoeffs[coeffI]
+//            );
     }
 }
 

src/foam/matrices/blockLduMatrix/BlockLduSolvers/BlockILU/BlockILUSolver.C

@@ -0,0 +1,125 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | foam-extend: Open Source CFD
+   \\    /   O peration     | Version:     4.0
+    \\  /    A nd           | Web:         http://www.foam-extend.org
+     \\/     M anipulation  | For copyright notice see file Copyright
+-------------------------------------------------------------------------------
+License
+    This file is part of foam-extend.
+
+    foam-extend 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.
+
+    foam-extend 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 foam-extend.  If not, see <http://www.gnu.org/licenses/>.
+
+Description
+    ILU solver for symmetric and asymmetric matrices.  In
+    order to improve efficiency, the residual is evaluated after every
+    nSweeps sweeps.
+
+\*---------------------------------------------------------------------------*/
+
+#include "BlockILUSolver.H"
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+// Construct from matrix and solver data stream
+template<class Type>
+Foam::BlockILUSolver<Type>::BlockILUSolver
+(
+    const word& fieldName,
+    const BlockLduMatrix<Type>& matrix,
+    const dictionary& dict
+)
+:
+    BlockIterativeSolver<Type>
+    (
+        fieldName,
+        matrix,
+        dict
+    ),
+    ilu_(matrix),
+    nSweeps_(readLabel(this->dict().lookup("nSweeps")))
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+template<class Type>
+typename Foam::BlockSolverPerformance<Type>
+Foam::BlockILUSolver<Type>::solve
+(
+    Field<Type>& x,
+    const Field<Type>& b
+)
+{
+    // Create local references to avoid the spread this-> ugliness
+    const BlockLduMatrix<Type>& matrix = this->matrix_;
+
+    // Prepare solver performance
+    BlockSolverPerformance<Type> solverPerf
+    (
+        typeName,
+        this->fieldName()
+    );
+
+    scalar norm = this->normFactor(x, b);
+
+    Field<Type> residual(x.size());
+
+    // Calculate residual.  Note: sign of residual swapped for efficiency
+    matrix.Amul(residual, x);
+
+    forAll (b, i)
+    {
+        residual[i] = b[i] - residual[i];
+    }
+
+    solverPerf.initialResidual() = gSum(cmptMag(residual))/norm;
+    solverPerf.finalResidual() = solverPerf.initialResidual();
+
+    // Check convergence, solve if not converged
+
+    if (!this->stop(solverPerf))
+    {
+        // Iteration loop
+
+        Field<Type> xCorr(x.size());
+
+        do
+        {
+            for (label i = 0; i < nSweeps_; i++)
+            {
+                ilu_.precondition(xCorr, residual);
+
+                x += xCorr;
+
+                solverPerf.nIterations()++;
+            }
+
+            // Re-calculate residual
+            matrix.Amul(residual, x);
+
+            forAll (b, i)
+            {
+                residual[i] = b[i] - residual[i];
+            }
+
+            solverPerf.finalResidual() = gSum(cmptMag(residual))/norm;
+        } while (!this->stop(solverPerf));
+    }
+
+    return solverPerf;
+}
+
+
+// ************************************************************************* //

src/foam/matrices/blockLduMatrix/BlockLduSolvers/BlockILU/BlockILUSolver.H

@@ -0,0 +1,125 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | foam-extend: Open Source CFD
+   \\    /   O peration     | Version:     4.0
+    \\  /    A nd           | Web:         http://www.foam-extend.org
+     \\/     M anipulation  | For copyright notice see file Copyright
+-------------------------------------------------------------------------------
+License
+    This file is part of foam-extend.
+
+    foam-extend 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.
+
+    foam-extend 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 foam-extend.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    BlockILUSolver
+
+Description
+    ILU solver
+
+Author
+    Hrvoje Jasak, Wikki Ltd.  All rights reserved
+
+SourceFiles
+    BlockILUSolver.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef BlockILUSolver_H
+#define BlockILUSolver_H
+
+#include "blockLduSolvers.H"
+#include "BlockIterativeSolver.H"
+#include "BlockSolverPerformance.H"
+#include "blockCholeskyPrecons.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                       Class BlockILUSolver Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class Type>
+class BlockILUSolver
+:
+    public BlockIterativeSolver<Type>
+{
+    // Private data
+
+        //- Gauss-Seidel preconditioner
+        BlockCholeskyPrecon<Type> ilu_;
+
+        //- Number of sweeps before evaluating residual
+        label nSweeps_;
+
+
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct
+        BlockILUSolver(const BlockILUSolver<Type>&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const BlockILUSolver<Type>&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("ILU");
+
+
+    // Constructors
+
+        //- Construct from matrix components and solver data stream
+        BlockILUSolver
+        (
+            const word& fieldName,
+            const BlockLduMatrix<Type>& matrix,
+            const dictionary& dict
+        );
+
+
+        //- Destructor
+        virtual ~BlockILUSolver()
+        {}
+
+
+    // Member Functions
+
+        //- Solve the matrix with this solver
+        virtual BlockSolverPerformance<Type> solve
+        (
+            Field<Type>& x,
+            const Field<Type>& b
+        );
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+#   include "BlockILUSolver.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/foam/matrices/blockLduMatrix/BlockLduSolvers/BlockILU/blockILUSolvers.C

@@ -0,0 +1,48 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | foam-extend: Open Source CFD
+   \\    /   O peration     | Version:     4.0
+    \\  /    A nd           | Web:         http://www.foam-extend.org
+     \\/     M anipulation  | For copyright notice see file Copyright
+-------------------------------------------------------------------------------
+License
+    This file is part of foam-extend.
+
+    foam-extend 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.
+
+    foam-extend 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 foam-extend.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "blockLduMatrices.H"
+#include "blockILUSolvers.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+makeBlockSolvers(blockILUSolver);
+
+addSymSolverToBlockMatrix(blockILUSolver);
+
+addAsymSolverToBlockMatrix(blockILUSolver);
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// ************************************************************************* //

src/foam/matrices/blockLduMatrix/BlockLduSolvers/BlockILU/blockILUSolvers.H

@@ -0,0 +1,63 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | foam-extend: Open Source CFD
+   \\    /   O peration     | Version:     4.0
+    \\  /    A nd           | Web:         http://www.foam-extend.org
+     \\/     M anipulation  | For copyright notice see file Copyright
+-------------------------------------------------------------------------------
+License
+    This file is part of foam-extend.
+
+    foam-extend 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.
+
+    foam-extend 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 foam-extend.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    BlockILUSolvers
+
+Description
+    Typedefs for BlockILUSolvers
+
+Author
+    Hrvoje Jasak, Wikki Ltd.  All rights reserved
+
+SourceFiles
+    blockILUSolvers.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef blockILUSolvers_H
+#define blockILUSolvers_H
+
+#include "BlockILUSolver.H"
+#include "fieldTypes.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+typedef BlockILUSolver<scalar> blockILUSolverScalar;
+typedef BlockILUSolver<vector> blockILUSolverVector;
+typedef BlockILUSolver<tensor> blockILUSolverTensor;
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/foam/matrices/blockLduMatrix/BlockLduSolvers/blockVectorNSolvers.C

@@ -43,6 +43,7 @@ License
 #include "BlockBiCGStabSolver.H"
 #include "BlockCGSolver.H"
 #include "BlockGaussSeidelSolver.H"
+#include "BlockILUSolver.H"
 #include "BlockGMRESSolver.H"
 
 // KRJ: 2012-12-15: Multigrid solver
@@ -138,6 +139,11 @@ makeBlockSolverTypeName(block##Type##GaussSeidelSolver);                      \
 addSolverToBlockMatrix(Type, block##Type##GaussSeidelSolver, symMatrix);      \
 addSolverToBlockMatrix(Type, block##Type##GaussSeidelSolver, asymMatrix);     \
                                                                               \
+typedef BlockILUSolver<type > block##Type##ILUSolver;                         \
+makeBlockSolverTypeName(block##Type##ILUSolver);                              \
+addSolverToBlockMatrix(Type, block##Type##ILUSolver, symMatrix);              \
+addSolverToBlockMatrix(Type, block##Type##ILUSolver, asymMatrix);             \
+                                                                              \
 typedef BlockGMRESSolver<type > block##Type##GMRESSolver;                     \
 makeBlockSolverTypeName(block##Type##GMRESSolver);                            \
 addSolverToBlockMatrix(Type, block##Type##GMRESSolver, symMatrix);            \