src/finiteVolume/fvMesh/extendedStencil/cellToFace/extendedUpwindCellToFaceStencilTemplates.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 "extendedCellToFaceStencil.H"
  27
  28 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
  29
  30 template<class Type>
  31 Foam::tmp<Foam::GeometricField<Type, Foam::fvsPatchField, Foam::surfaceMesh> >
  32 Foam::extendedUpwindCellToFaceStencil::weightedSum
  33 (
  34     const surfaceScalarField& phi,
  35     const GeometricField<Type, fvPatchField, volMesh>& fld,
  36     const List<List<scalar> >& ownWeights,
  37     const List<List<scalar> >& neiWeights
  38 ) const
  39 {
  40     const fvMesh& mesh = fld.mesh();
  41
  42     // Collect internal and boundary values
  43     List<List<Type> > ownFld;
  44     collectData(ownMap(), ownStencil(), fld, ownFld);
  45     List<List<Type> > neiFld;
  46     collectData(neiMap(), neiStencil(), fld, neiFld);
  47
  48     tmp<GeometricField<Type, fvsPatchField, surfaceMesh> > tsfCorr
  49     (
  50         new GeometricField<Type, fvsPatchField, surfaceMesh>
  51         (
  52             IOobject
  53             (
  54                 fld.name(),
  55                 mesh.time().timeName(),
  56                 mesh,
  57                 IOobject::NO_READ,
  58                 IOobject::NO_WRITE,
  59                 false
  60             ),
  61             mesh,
  62             dimensioned<Type>
  63             (
  64                 fld.name(),
  65                 fld.dimensions(),
  66                 pTraits<Type>::zero
  67             )
  68         )
  69     );
  70     GeometricField<Type, fvsPatchField, surfaceMesh>& sf = tsfCorr();
  71
  72     // Internal faces
  73     for (label faceI = 0; faceI < mesh.nInternalFaces(); faceI++)
  74     {
  75         if (phi[faceI] > 0)
  76         {
  77             // Flux out of owner. Use upwind (= owner side) stencil.
  78             const List<Type>& stField = ownFld[faceI];
  79             const List<scalar>& stWeight = ownWeights[faceI];
  80
  81             forAll(stField, i)
  82             {
  83                 sf[faceI] += stField[i]*stWeight[i];
  84             }
  85         }
  86         else
  87         {
  88             const List<Type>& stField = neiFld[faceI];
  89             const List<scalar>& stWeight = neiWeights[faceI];
  90
  91             forAll(stField, i)
  92             {
  93                 sf[faceI] += stField[i]*stWeight[i];
  94             }
  95         }
  96     }
  97
  98     // Boundaries. Either constrained or calculated so assign value
  99     // directly (instead of nicely using operator==)
 100     typename GeometricField<Type, fvsPatchField, surfaceMesh>::
 101         GeometricBoundaryField& bSfCorr = sf.boundaryField();
 102
 103     forAll(bSfCorr, patchi)
 104     {
 105         fvsPatchField<Type>& pSfCorr = bSfCorr[patchi];
 106
 107         if (pSfCorr.coupled())
 108         {
 109             label faceI = pSfCorr.patch().start();
 110
 111             forAll(pSfCorr, i)
 112             {
 113                 if (phi.boundaryField()[patchi][i] > 0)
 114                 {
 115                     // Flux out of owner. Use upwind (= owner side) stencil.
 116                     const List<Type>& stField = ownFld[faceI];
 117                     const List<scalar>& stWeight = ownWeights[faceI];
 118
 119                     forAll(stField, j)
 120                     {
 121                         pSfCorr[i] += stField[j]*stWeight[j];
 122                     }
 123                 }
 124                 else
 125                 {
 126                     const List<Type>& stField = neiFld[faceI];
 127                     const List<scalar>& stWeight = neiWeights[faceI];
 128
 129                     forAll(stField, j)
 130                     {
 131                         pSfCorr[i] += stField[j]*stWeight[j];
 132                     }
 133                 }
 134                 faceI++;
 135             }
 136         }
 137     }
 138
 139     return tsfCorr;
 140 }
 141
 142
 143 // ************************************************************************* //