src/coupledMatrix/coupledLduSolver/smoothSolver/coupledSmoothSolver.C

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