applications/test/nearWallDist-wave/testWallDist2.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 1991-2010 OpenCFD Ltd.
   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 Description
  25     Wave propagation of nearwall distance through grid. Every iteration
  26     information goes through one layer of cells.
  27
  28 \*---------------------------------------------------------------------------*/
  29
  30 #include "fvCFD.H"
  31 #include "wallFvPatch.H"
  32 #include "MeshWave.H"
  33 #include "wallPoint.H"
  34
  35
  36 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  37 // Main program:
  38
  39 int main(int argc, char *argv[])
  40 {
  41 #   include "setRootCase.H"
  42 #   include "createTime.H"
  43 #   include "createMesh.H"
  44
  45     Info<< "Mesh read in = "
  46         << runTime.cpuTimeIncrement()
  47         << " s\n" << endl << endl;
  48
  49     Info<< "Creating field wDistNC\n" << endl;
  50     volScalarField wallDistUncorrected
  51     (
  52         IOobject
  53         (
  54             "wDistNC",
  55             runTime.timeName(),
  56             mesh,
  57             IOobject::NO_READ,
  58             IOobject::AUTO_WRITE
  59         ),
  60         mesh,
  61         dimensionedScalar
  62         (
  63             "wallDist",
  64             dimensionSet(0, 1, 0, 0, 0),
  65             0.0
  66         )
  67     );
  68
  69     //
  70     // Set initial changed faces: set wallPoint for wall faces to wall centre
  71     //
  72
  73     // Count walls
  74     label nWalls = 0;
  75     forAll (mesh.boundary(), patchI)
  76     {
  77         const fvPatch& patch = mesh.boundary()[patchI];
  78
  79         if (isA<wallFvPatch>(patch))
  80         {
  81             nWalls += patch.size();
  82         }
  83     }
  84
  85     List<wallPoint> faceDist(nWalls);
  86     labelList changedFaces(nWalls);
  87
  88     label nChangedFaces = 0;
  89     forAll (mesh.boundary(), patchI)
  90     {
  91         const fvPatch& patch = mesh.boundary()[patchI];
  92
  93         if (isA<wallFvPatch>(patch))
  94         {
  95             forAll (patch.Cf(), patchFaceI)
  96             {
  97                 const polyPatch& polyPatch = mesh.boundaryMesh()[patchI];
  98
  99                 label meshFaceI = polyPatch.start() + patchFaceI;
 100
 101                 changedFaces[nChangedFaces] = meshFaceI;
 102
 103                 faceDist[nChangedFaces] =
 104                     wallPoint(patch.Cf()[patchFaceI], 0.0);
 105
 106                 nChangedFaces++;
 107             }
 108         }
 109     }
 110
 111
 112     MeshWave<wallPoint> wallDistCalc
 113     (
 114         mesh,
 115         changedFaces,
 116         faceDist,
 117         0           // max iterations
 118     );
 119
 120     Info<< "\nStarting time loop\n" << endl;
 121
 122     while (runTime.loop())
 123     {
 124         Info<< "Time = " << runTime.timeName() << endl;
 125
 126
 127         label nCells = wallDistCalc.faceToCell();
 128
 129         Info<< "    Total changed cells   : " << nCells << endl;
 130
 131         if (nCells == 0)
 132         {
 133             break;
 134         }
 135
 136
 137         label nFaces = wallDistCalc.cellToFace();
 138
 139         Info<< "    Total changed faces   : " << nFaces << endl;
 140
 141         if (nFaces == 0)
 142         {
 143             break;
 144         }
 145
 146
 147         //
 148         // Copy face and cell values into field
 149         //
 150
 151         const List<wallPoint>& cellInfo = wallDistCalc.allCellInfo();
 152         const List<wallPoint>& faceInfo = wallDistCalc.allFaceInfo();
 153
 154         label nIllegal = 0;
 155
 156         // Copy cell values
 157         forAll(cellInfo, cellI)
 158         {
 159             scalar dist = cellInfo[cellI].distSqr();
 160             if (cellInfo[cellI].valid())
 161             {
 162                 wallDistUncorrected[cellI] = Foam::sqrt(dist);
 163             }
 164             else
 165             {
 166                 wallDistUncorrected[cellI] = -1;
 167                 nIllegal++;
 168             }
 169         }
 170
 171         // Copy boundary values
 172         forAll (wallDistUncorrected.boundaryField(), patchI)
 173         {
 174             fvPatchScalarField& patchField =
 175                 wallDistUncorrected.boundaryField()[patchI];
 176
 177             forAll(patchField, patchFaceI)
 178             {
 179                 label meshFaceI =
 180                     patchField.patch().patch().start() + patchFaceI;
 181
 182                 scalar dist = faceInfo[meshFaceI].distSqr();
 183                 if (faceInfo[meshFaceI].valid())
 184                 {
 185                     patchField[patchFaceI] = Foam::sqrt(dist);
 186                 }
 187                 else
 188                 {
 189                     patchField[patchFaceI] = dist;
 190                     nIllegal++;
 191                 }
 192             }
 193         }
 194
 195         Info<< "nIllegal:" << nIllegal << endl;
 196
 197
 198         //
 199         // Write it
 200         //
 201
 202         wallDistUncorrected.write();
 203
 204         Info<< "ExecutionTime = "
 205             << runTime.elapsedCpuTime()
 206             << " s\n" << endl << endl;
 207     }
 208
 209     Info<< "End\n" << endl;
 210
 211     return 0;
 212 }
 213
 214
 215 // ************************************************************************* //