FreeFOAM The Cross-Platform CFD Toolkit
patchDataWave.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 "patchDataWave.H"
27 #include <OpenFOAM/MeshWave.H>
28 
29 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
30 
31 // Set initial set of changed faces (= all wall faces)
32 template<class TransferType>
34 (
35  const labelHashSet& patchIDs,
36  labelList& changedFaces,
37  List<TransferType>& faceDist
38 ) const
39 {
40  const polyMesh& mesh = cellDistFuncs::mesh();
41 
42  label nChangedFaces = 0;
43 
44  forAll(mesh.boundaryMesh(), patchI)
45  {
46  if (patchIDs.found(patchI))
47  {
48  const polyPatch& patch = mesh.boundaryMesh()[patchI];
49 
50  const Field<Type>& patchField = initialPatchValuePtrs_[patchI];
51 
52  forAll(patch.faceCentres(), patchFaceI)
53  {
54  label meshFaceI = patch.start() + patchFaceI;
55 
56  changedFaces[nChangedFaces] = meshFaceI;
57 
58  faceDist[nChangedFaces] =
59  TransferType
60  (
61  patch.faceCentres()[patchFaceI],
62  patchField[patchFaceI],
63  0.0
64  );
65 
66  nChangedFaces++;
67  }
68  }
69  }
70 }
71 
72 
73 // Copy from MeshWave data into *this (distance) and field_ (transported data)
74 template<class TransferType>
76 (
77  const MeshWave<TransferType>& waveInfo
78 )
79 {
80  const polyMesh& mesh = cellDistFuncs::mesh();
81 
82  const List<TransferType>& cellInfo = waveInfo.allCellInfo();
83  const List<TransferType>& faceInfo = waveInfo.allFaceInfo();
84 
85  label nIllegal = 0;
86 
87  // Copy cell values
88  distance_.setSize(cellInfo.size());
89 
90  forAll(cellInfo, cellI)
91  {
92  const TransferType & wpn = cellInfo[cellI];
93 
94  scalar dist = wpn.distSqr();
95 
96  if (cellInfo[cellI].valid())
97  {
98  distance_[cellI] = Foam::sqrt(dist);
99 
100  cellData_[cellI] = cellInfo[cellI].data();
101  }
102  else
103  {
104  // Illegal/unset value. What to do with data?
105 
106  distance_[cellI] = dist;
107 
108  //cellData_[cellI] = wallPoint::greatPoint;
109  cellData_[cellI] = cellInfo[cellI].data();
110 
111  nIllegal++;
112  }
113  }
114 
115  // Copy boundary values
116  forAll(patchDistance_, patchI)
117  {
118  const polyPatch& patch = mesh.boundaryMesh()[patchI];
119 
120  // Allocate storage for patchDistance
121  scalarField* patchFieldPtr = new scalarField(patch.size());
122 
123  patchDistance_.set(patchI, patchFieldPtr);
124 
125  scalarField& patchField = *patchFieldPtr;
126 
127  // Allocate storage for patchData
128  Field<Type>* patchDataFieldPtr = new Field<Type>(patch.size());
129 
130  patchData_.set(patchI, patchDataFieldPtr);
131 
132  Field<Type>& patchDataField = *patchDataFieldPtr;
133 
134  // Copy distance and data
135  forAll(patchField, patchFaceI)
136  {
137  label meshFaceI = patch.start() + patchFaceI;
138 
139  scalar dist = faceInfo[meshFaceI].distSqr();
140 
141  if (faceInfo[meshFaceI].valid())
142  {
143  // Adding SMALL to avoid problems with /0 in the turbulence
144  // models
145  patchField[patchFaceI] = Foam::sqrt(dist) + SMALL;
146 
147  patchDataField[patchFaceI] = faceInfo[meshFaceI].data();
148  }
149  else
150  {
151  // Illegal/unset value. What to do with data?
152 
153  patchField[patchFaceI] = dist;
154 
155  //patchDataField[patchFaceI] = wallPoint::greatPoint;
156  patchDataField[patchFaceI] = faceInfo[meshFaceI].data();
157 
158  nIllegal++;
159  }
160  }
161  }
162 
163  return nIllegal;
164 }
165 
166 
167 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
168 
169 // Construct from components
170 template<class TransferType>
172 (
173  const polyMesh& mesh,
174  const labelHashSet& patchIDs,
175  const UPtrList<Field<Type> >& initialPatchValuePtrs,
176  const bool correctWalls
177 )
178 :
179  cellDistFuncs(mesh),
180  patchIDs_(patchIDs),
181  initialPatchValuePtrs_(initialPatchValuePtrs),
182  correctWalls_(correctWalls),
183  nUnset_(0),
184  distance_(mesh.nCells()),
185  patchDistance_(mesh.boundaryMesh().size()),
186  cellData_(mesh.nCells()),
187  patchData_(mesh.boundaryMesh().size())
188 {
190 }
191 
192 
193 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
194 
195 template<class TransferType>
197 {}
198 
199 
200 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
201 
202 // Correct for mesh geom/topo changes
203 template<class TransferType>
205 {
206  //
207  // Set initial changed faces: set TransferType for wall faces
208  // to wall centre.
209  //
210 
211  // Count walls
212  label nWalls = sumPatchSize(patchIDs_);
213 
214  List<TransferType> faceDist(nWalls);
215  labelList changedFaces(nWalls);
216 
217  setChangedFaces(patchIDs_, changedFaces, faceDist);
218 
219  //
220  // Do calculate wall distance by 'growing' from faces.
221  //
222 
223  MeshWave<TransferType> waveInfo
224  (
225  mesh(),
226  changedFaces,
227  faceDist,
228  mesh().globalData().nTotalCells() // max iterations
229  );
230 
231 
232  //
233  // Copy distance into return field
234  //
235 
236  nUnset_ = getValues(waveInfo);
237 
238  //
239  // Correct wall cells for true distance
240  //
241 
242  if (correctWalls_)
243  {
244  Map<label> nearestFace(2 * nWalls);
245 
246  // Get distance and indices of nearest face
247  correctBoundaryFaceCells
248  (
249  patchIDs_,
250  distance_,
251  nearestFace
252  );
253 
254  correctBoundaryPointCells
255  (
256  patchIDs_,
257  distance_,
258  nearestFace
259  );
260 
261  // Transfer data from nearest face to cell
262  const List<TransferType>& faceInfo = waveInfo.allFaceInfo();
263 
264  const labelList wallCells(nearestFace.toc());
265 
266  forAll(wallCells, wallCellI)
267  {
268  label cellI = wallCells[wallCellI];
269 
270  label faceI = nearestFace[cellI];
271 
272  cellData_[cellI] = faceInfo[faceI].data();
273  }
274  }
275 }
276 
277 
278 // ************************ vim: set sw=4 sts=4 et: ************************ //