src/OpenFOAM/matrices/lduMatrix/smoothers/DIC/DICSmoother.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
   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
  13     the Free Software Foundation, either version 3 of the License, or
  14     (at your 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, see <http://www.gnu.org/licenses/>.
  23
  24 \*---------------------------------------------------------------------------*/
  25
  26 #include "DICSmoother.H"
  27 #include "DICPreconditioner.H"
  28
  29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  30
  31 namespace Foam
  32 {
  33     defineTypeNameAndDebug(DICSmoother, 0);
  34
  35     lduMatrix::smoother::addsymMatrixConstructorToTable<DICSmoother>
  36         addDICSmootherSymMatrixConstructorToTable_;
  37 }
  38
  39
  40 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  41
  42 Foam::DICSmoother::DICSmoother
  43 (
  44     const word& fieldName,
  45     const lduMatrix& matrix,
  46     const FieldField<Field, scalar>& interfaceBouCoeffs,
  47     const FieldField<Field, scalar>& interfaceIntCoeffs,
  48     const lduInterfaceFieldPtrsList& interfaces
  49 )
  50 :
  51     lduMatrix::smoother
  52     (
  53         fieldName,
  54         matrix,
  55         interfaceBouCoeffs,
  56         interfaceIntCoeffs,
  57         interfaces
  58     ),
  59     rD_(matrix_.diag())
  60 {
  61     DICPreconditioner::calcReciprocalD(rD_, matrix_);
  62 }
  63
  64
  65 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
  66
  67 void Foam::DICSmoother::smooth
  68 (
  69     scalarField& psi,
  70     const scalarField& source,
  71     const direction cmpt,
  72     const label nSweeps
  73 ) const
  74 {
  75     const scalar* const __restrict__ rDPtr = rD_.begin();
  76     const scalar* const __restrict__ upperPtr = matrix_.upper().begin();
  77     const label* const __restrict__ uPtr =
  78         matrix_.lduAddr().upperAddr().begin();
  79     const label* const __restrict__ lPtr =
  80         matrix_.lduAddr().lowerAddr().begin();
  81
  82     // Temporary storage for the residual
  83     scalarField rA(rD_.size());
  84     scalar* __restrict__ rAPtr = rA.begin();
  85
  86     for (label sweep=0; sweep<nSweeps; sweep++)
  87     {
  88         matrix_.residual
  89         (
  90             rA,
  91             psi,
  92             source,
  93             interfaceBouCoeffs_,
  94             interfaces_,
  95             cmpt
  96         );
  97
  98         rA *= rD_;
  99
 100         register label nFaces = matrix_.upper().size();
 101         for (register label face=0; face<nFaces; face++)
 102         {
 103             register label u = uPtr[face];
 104             rAPtr[u] -= rDPtr[u]*upperPtr[face]*rAPtr[lPtr[face]];
 105         }
 106
 107         register label nFacesM1 = nFaces - 1;
 108         for (register label face=nFacesM1; face>=0; face--)
 109         {
 110             register label l = lPtr[face];
 111             rAPtr[l] -= rDPtr[l]*upperPtr[face]*rAPtr[uPtr[face]];
 112         }
 113
 114         psi += rA;
 115     }
 116 }
 117
 118
 119 // ************************************************************************* //