FreeFOAM The Cross-Platform CFD Toolkit
patchWave.C
Go to the documentation of this file.
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 \*---------------------------------------------------------------------------*/
25 
26 #include "patchWave.H"
27 #include <OpenFOAM/polyMesh.H>
28 #include <meshTools/wallPoint.H>
29 #include <OpenFOAM/MeshWave.H>
31 
32 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
33 
34 void Foam::patchWave::setChangedFaces
35 (
36  const labelHashSet& patchIDs,
37  labelList& changedFaces,
38  List<wallPoint>& faceDist
39 ) const
40 {
41  const polyMesh& mesh = cellDistFuncs::mesh();
42 
43  label nChangedFaces = 0;
44 
45  forAll(mesh.boundaryMesh(), patchI)
46  {
47  if (patchIDs.found(patchI))
48  {
49  const polyPatch& patch = mesh.boundaryMesh()[patchI];
50 
51  forAll(patch.faceCentres(), patchFaceI)
52  {
53  label meshFaceI = patch.start() + patchFaceI;
54 
55  changedFaces[nChangedFaces] = meshFaceI;
56 
57  faceDist[nChangedFaces] =
58  wallPoint
59  (
60  patch.faceCentres()[patchFaceI],
61  0.0
62  );
63 
64  nChangedFaces++;
65  }
66  }
67  }
68 }
69 
70 
71 Foam::label Foam::patchWave::getValues(const MeshWave<wallPoint>& waveInfo)
72 {
73  const List<wallPoint>& cellInfo = waveInfo.allCellInfo();
74  const List<wallPoint>& faceInfo = waveInfo.allFaceInfo();
75 
76  label nIllegal = 0;
77 
78  // Copy cell values
79  distance_.setSize(cellInfo.size());
80 
81  forAll(cellInfo, cellI)
82  {
83  scalar dist = cellInfo[cellI].distSqr();
84 
85  if (cellInfo[cellI].valid())
86  {
87  distance_[cellI] = Foam::sqrt(dist);
88  }
89  else
90  {
91  distance_[cellI] = dist;
92 
93  nIllegal++;
94  }
95  }
96 
97  // Copy boundary values
98  forAll(patchDistance_, patchI)
99  {
100  const polyPatch& patch = mesh().boundaryMesh()[patchI];
101 
102  // Allocate storage for patchDistance
103  scalarField* patchDistPtr = new scalarField(patch.size());
104 
105  patchDistance_.set(patchI, patchDistPtr);
106 
107  scalarField& patchField = *patchDistPtr;
108 
109  forAll(patchField, patchFaceI)
110  {
111  label meshFaceI = patch.start() + patchFaceI;
112 
113  scalar dist = faceInfo[meshFaceI].distSqr();
114 
115  if (faceInfo[meshFaceI].valid())
116  {
117  // Adding SMALL to avoid problems with /0 in the turbulence
118  // models
119  patchField[patchFaceI] = Foam::sqrt(dist) + SMALL;
120  }
121  else
122  {
123  patchField[patchFaceI] = dist;
124 
125  nIllegal++;
126  }
127  }
128  }
129  return nIllegal;
130 }
131 
132 
133 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
134 
135 // Construct from components
137 (
138  const polyMesh& mesh,
139  const labelHashSet& patchIDs,
140  const bool correctWalls
141 )
142 :
143  cellDistFuncs(mesh),
144  patchIDs_(patchIDs),
145  correctWalls_(correctWalls),
146  nUnset_(0),
147  distance_(mesh.nCells()),
148  patchDistance_(mesh.boundaryMesh().size())
149 {
151 }
152 
153 
154 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
155 
157 {}
158 
159 
160 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
161 
162 // Correct for mesh geom/topo changes. Might be more intelligent in the
163 // future (if only small topology change)
165 {
166  //
167  // Set initial changed faces: set wallPoint for wall faces to wall centre
168  //
169 
170  // Count walls
171  label nWalls = sumPatchSize(patchIDs_);
172 
173  List<wallPoint> faceDist(nWalls);
174  labelList changedFaces(nWalls);
175 
176  // Set to faceDist information to facecentre on walls.
177  setChangedFaces(patchIDs_, changedFaces, faceDist);
178 
179 
180  //
181  // Do calculate wall distance by 'growing' from faces.
182  //
183 
184  MeshWave<wallPoint> waveInfo
185  (
186  mesh(),
187  changedFaces,
188  faceDist,
189  mesh().globalData().nTotalCells() // max iterations
190  );
191 
192 
193  //
194  // Copy distance into return field
195  //
196 
197  nUnset_ = getValues(waveInfo);
198 
199  //
200  // Correct wall cells for true distance
201  //
202 
203  if (correctWalls_)
204  {
205  Map<label> nearestFace(2 * nWalls);
206 
207  correctBoundaryFaceCells
208  (
209  patchIDs_,
210  distance_,
211  nearestFace
212  );
213 
214  correctBoundaryPointCells
215  (
216  patchIDs_,
217  distance_,
218  nearestFace
219  );
220  }
221 }
222 
223 
224 // ************************ vim: set sw=4 sts=4 et: ************************ //