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