src/coupledMatrix/coupledLduSolver/smoothSolver/coupledSmoothSolver.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | foam-extend: Open Source CFD
   4    \\    /   O peration     | Version:     3.2
   5     \\  /    A nd           | Web:         http://www.foam-extend.org
   6      \\/     M anipulation  | For copyright notice see file Copyright
   7 -------------------------------------------------------------------------------
   8 License
   9     This file is part of foam-extend.
  10
  11     foam-extend is free software: you can redistribute it and/or modify it
  12     under the terms of the GNU General Public License as published by the
  13     Free Software Foundation, either version 3 of the License, or (at your
  14     option) any later version.
  15
  16     foam-extend is distributed in the hope that it will be useful, but
  17     WITHOUT ANY WARRANTY; without even the implied warranty of
  18     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  19     General Public License for more details.
  20
  21     You should have received a copy of the GNU General Public License
  22     along with foam-extend.  If not, see <http://www.gnu.org/licenses/>.
  23
  24 Description
  25     Smoother-solver for coupled diagonal lduMatrices.
  26
  27 Author
  28     Hrvoje Jasak, Wikki Ltd.  All rights reserved
  29
  30 \*---------------------------------------------------------------------------*/
  31
  32 #include "coupledSmoothSolver.H"
  33 #include "addToRunTimeSelectionTable.H"
  34 #include "coupledLduSmoother.H"
  35
  36 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  37
  38 namespace Foam
  39 {
  40     defineTypeNameAndDebug(coupledSmoothSolver, 0);
  41
  42     coupledLduSolver::addsymMatrixConstructorToTable<coupledSmoothSolver>
  43         addGaussSeidelSolverSymMatrixConstructorToTable_;
  44
  45     coupledLduSolver::addasymMatrixConstructorToTable<coupledSmoothSolver>
  46         addGaussSeidelSolverAsymMatrixConstructorToTable_;
  47 }
  48
  49
  50 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  51
  52 //- Construct from matrix
  53 Foam::coupledSmoothSolver::coupledSmoothSolver
  54 (
  55     const word& fieldName,
  56     const coupledLduMatrix& matrix,
  57     const PtrList<FieldField<Field, scalar> >& bouCoeffs,
  58     const PtrList<FieldField<Field, scalar> >& intCoeffs,
  59     const lduInterfaceFieldPtrsListList& interfaces,
  60     const dictionary& solverData
  61 )
  62 :
  63     coupledIterativeSolver
  64     (
  65         fieldName,
  66         matrix,
  67         bouCoeffs,
  68         intCoeffs,
  69         interfaces,
  70         solverData
  71     ),
  72     nSweeps_(1)
  73 {
  74     readControls();
  75 }
  76
  77
  78 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  79
  80 void Foam::coupledSmoothSolver::readControls()
  81 {
  82     coupledIterativeSolver::readControls();
  83     dict().readIfPresent("nSweeps", nSweeps_);
  84 }
  85
  86
  87 Foam::coupledSolverPerformance Foam::coupledSmoothSolver::solve
  88 (
  89     FieldField<Field, scalar>& x,
  90     const FieldField<Field, scalar>& b,
  91     const direction cmpt
  92 ) const
  93 {
  94     // Prepare solver performance
  95     coupledSolverPerformance solverPerf(typeName, fieldName());
  96
  97     // Do a minimum number of sweeps
  98     // HJ, 19/Jan/2009
  99     if (minIter() > 0)
 100     {
 101         autoPtr<coupledLduSmoother> smootherPtr = coupledLduSmoother::New
 102         (
 103             matrix_,
 104             bouCoeffs_,
 105             intCoeffs_,
 106             interfaces_,
 107             dict()
 108         );
 109
 110         smootherPtr->smooth
 111         (
 112             x,
 113             b,
 114             cmpt,
 115             minIter()
 116         );
 117
 118         solverPerf.nIterations() += minIter();
 119     }
 120
 121     // Now do normal sweeps.  HJ, 19/Jan/2009
 122
 123     FieldField<Field, scalar> Ax(x.size());
 124     FieldField<Field, scalar> temp(x.size());
 125
 126     forAll (x, rowI)
 127     {
 128         Ax.set(rowI, new scalarField(x[rowI].size(), 0));
 129         temp.set(rowI, new scalarField(x[rowI].size(), 0));
 130     }
 131
 132     // Calculate initial residual. Note: for efficiency, swapping sign
 133     matrix_.Amul(Ax, x, bouCoeffs_, interfaces_, cmpt);
 134
 135     scalar normFactor = this->normFactor(x, b, Ax, temp, cmpt);
 136
 137     Ax -= b;
 138
 139     solverPerf.initialResidual() = gSumMag(Ax)/normFactor;
 140     solverPerf.finalResidual() = solverPerf.initialResidual();
 141
 142     if (!solverPerf.checkConvergence(tolerance_, relTolerance_))
 143     {
 144         autoPtr<coupledLduSmoother> smootherPtr =
 145             coupledLduSmoother::New
 146             (
 147                 matrix_,
 148                 bouCoeffs_,
 149                 intCoeffs_,
 150                 interfaces_,
 151                 dict()
 152             );
 153
 154         // Smoothing loop
 155         do
 156         {
 157             smootherPtr->smooth(x, b, cmpt, nSweeps_);
 158
 159             // Re-calculate residual
 160             matrix_.Amul(Ax, x, bouCoeffs_, interfaces_, cmpt);
 161             Ax -= b;
 162
 163             solverPerf.finalResidual() = gSumMag(Ax)/normFactor;
 164
 165             solverPerf.nIterations() += nSweeps_;
 166         } while (!stop(solverPerf));
 167     }
 168
 169     return solverPerf;
 170 }
 171
 172
 173 // ************************************************************************* //