src/fvMotionSolver/motionDiffusivity/inversePointDistance/inversePointDistanceDiffusivity.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 \*---------------------------------------------------------------------------*/
  26
  27 #include "inversePointDistanceDiffusivity.H"
  28 #include "addToRunTimeSelectionTable.H"
  29 #include "HashSet.H"
  30 #include "pointEdgePoint.H"
  31 #include "PointEdgeWave.H"
  32
  33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  34
  35 namespace Foam
  36 {
  37     defineTypeNameAndDebug(inversePointDistanceDiffusivity, 0);
  38
  39     addToRunTimeSelectionTable
  40     (
  41         motionDiffusivity,
  42         inversePointDistanceDiffusivity,
  43         Istream
  44     );
  45 }
  46
  47
  48 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  49
  50 Foam::inversePointDistanceDiffusivity::inversePointDistanceDiffusivity
  51 (
  52     const fvMotionSolver& mSolver,
  53     Istream& mdData
  54 )
  55 :
  56     uniformDiffusivity(mSolver, mdData),
  57     patchNames_(mdData)
  58 {
  59     correct();
  60 }
  61
  62
  63 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
  64
  65 Foam::inversePointDistanceDiffusivity::~inversePointDistanceDiffusivity()
  66 {}
  67
  68
  69 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
  70
  71 void Foam::inversePointDistanceDiffusivity::correct()
  72 {
  73     const polyMesh& mesh = mSolver().mesh();
  74     const polyBoundaryMesh& bdry = mesh.boundaryMesh();
  75
  76     labelHashSet patchSet(bdry.size());
  77
  78     label nPatchEdges = 0;
  79
  80     forAll (patchNames_, i)
  81     {
  82         label pID = bdry.findPatchID(patchNames_[i]);
  83
  84         if (pID > -1)
  85         {
  86             patchSet.insert(pID);
  87             nPatchEdges += bdry[pID].nEdges();
  88         }
  89     }
  90
  91     // Distance to wall on points and edges.
  92     List<pointEdgePoint> pointWallDist(mesh.nPoints());
  93     List<pointEdgePoint> edgeWallDist(mesh.nEdges());
  94
  95     {
  96         // Seeds
  97         List<pointEdgePoint> seedInfo(nPatchEdges);
  98         labelList seedPoints(nPatchEdges);
  99
 100         nPatchEdges = 0;
 101
 102         forAllConstIter(labelHashSet, patchSet, iter)
 103         {
 104             const polyPatch& patch = bdry[iter.key()];
 105
 106             const labelList& meshPoints = patch.meshPoints();
 107
 108             forAll(meshPoints, i)
 109             {
 110                 label pointI = meshPoints[i];
 111
 112                 if (!pointWallDist[pointI].valid())
 113                 {
 114                     // Not yet seeded
 115                     seedInfo[nPatchEdges] = pointEdgePoint
 116                     (
 117                         mesh.points()[pointI],
 118                         0.0
 119                     );
 120                     seedPoints[nPatchEdges] = pointI;
 121                     pointWallDist[pointI] = seedInfo[nPatchEdges];
 122
 123                     nPatchEdges++;
 124                 }
 125             }
 126         }
 127         seedInfo.setSize(nPatchEdges);
 128         seedPoints.setSize(nPatchEdges);
 129
 130         // Do calculations
 131         PointEdgeWave<pointEdgePoint> waveInfo
 132         (
 133             mesh,
 134             seedPoints,
 135             seedInfo,
 136
 137             pointWallDist,
 138             edgeWallDist,
 139             mesh.globalData().nTotalPoints() // max iterations
 140         );
 141     }
 142
 143
 144     for (label faceI=0; faceI<mesh.nInternalFaces(); faceI++)
 145     {
 146         const face& f = mesh.faces()[faceI];
 147
 148         scalar dist = 0;
 149
 150         forAll(f, fp)
 151         {
 152             dist += sqrt(pointWallDist[f[fp]].distSqr());
 153         }
 154         dist /= f.size();
 155
 156         faceDiffusivity_[faceI] = 1.0/dist;
 157     }
 158
 159     forAll(faceDiffusivity_.boundaryField(), patchI)
 160     {
 161         fvsPatchScalarField& bfld = faceDiffusivity_.boundaryField()[patchI];
 162
 163         if (patchSet.found(patchI))
 164         {
 165             const unallocLabelList& faceCells = bfld.patch().faceCells();
 166
 167             forAll(bfld, i)
 168             {
 169                 const cell& ownFaces = mesh.cells()[faceCells[i]];
 170
 171                 labelHashSet cPoints(4*ownFaces.size());
 172
 173                 scalar dist = 0;
 174
 175                 forAll(ownFaces, ownFaceI)
 176                 {
 177                     const face& f = mesh.faces()[ownFaces[ownFaceI]];
 178
 179                     forAll(f, fp)
 180                     {
 181                         if (cPoints.insert(f[fp]))
 182                         {
 183                             dist += sqrt(pointWallDist[f[fp]].distSqr());
 184                         }
 185                     }
 186                 }
 187                 dist /= cPoints.size();
 188
 189                 bfld[i] = 1.0/dist;
 190             }
 191         }
 192         else
 193         {
 194             label start = bfld.patch().patch().start();
 195
 196             forAll(bfld, i)
 197             {
 198                 const face& f = mesh.faces()[start+i];
 199
 200                 scalar dist = 0;
 201
 202                 forAll(f, fp)
 203                 {
 204                     dist += sqrt(pointWallDist[f[fp]].distSqr());
 205                 }
 206                 dist /= f.size();
 207
 208                 bfld[i] = 1.0/dist;
 209             }
 210         }
 211     }
 212 }
 213
 214
 215 // ************************************************************************* //