fvMeshTools.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) 2012-2016 OpenFOAM Foundation
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 "fvMeshTools.H"
27 
28 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
29 
30 // Adds patch if not yet there. Returns patchID.
32 (
33  fvMesh& mesh,
34  const polyPatch& patch,
35  const dictionary& patchFieldDict,
36  const word& defaultPatchFieldType,
37  const bool validBoundary
38 )
39 {
40  polyBoundaryMesh& polyPatches =
41  const_cast<polyBoundaryMesh&>(mesh.boundaryMesh());
42 
43  label patchi = polyPatches.findPatchID(patch.name());
44  if (patchi != -1)
45  {
46  // Already there
47  return patchi;
48  }
49 
50 
51  // Append at end unless there are processor patches
52  label insertPatchi = polyPatches.size();
53  label startFacei = mesh.nFaces();
54 
55  if (!isA<processorPolyPatch>(patch))
56  {
57  forAll(polyPatches, patchi)
58  {
59  const polyPatch& pp = polyPatches[patchi];
60 
61  if (isA<processorPolyPatch>(pp))
62  {
63  insertPatchi = patchi;
64  startFacei = pp.start();
65  break;
66  }
67  }
68  }
69 
70 
71  // Below is all quite a hack. Feel free to change once there is a better
72  // mechanism to insert and reorder patches.
73 
74  // Clear local fields and e.g. polyMesh parallelInfo.
75  mesh.clearOut();
76 
77  label sz = polyPatches.size();
78 
79  fvBoundaryMesh& fvPatches = const_cast<fvBoundaryMesh&>(mesh.boundary());
80 
81  // Add polyPatch at the end
82  polyPatches.setSize(sz+1);
83  polyPatches.set
84  (
85  sz,
86  patch.clone
87  (
88  polyPatches,
89  insertPatchi, //index
90  0, //size
91  startFacei //start
92  )
93  );
94  fvPatches.setSize(sz+1);
95  fvPatches.set
96  (
97  sz,
99  (
100  polyPatches[sz], // point to newly added polyPatch
101  mesh.boundary()
102  )
103  );
104 
105  addPatchFields<volScalarField>
106  (
107  mesh,
108  patchFieldDict,
109  defaultPatchFieldType,
110  Zero
111  );
112  addPatchFields<volVectorField>
113  (
114  mesh,
115  patchFieldDict,
116  defaultPatchFieldType,
117  Zero
118  );
119  addPatchFields<volSphericalTensorField>
120  (
121  mesh,
122  patchFieldDict,
123  defaultPatchFieldType,
124  Zero
125  );
126  addPatchFields<volSymmTensorField>
127  (
128  mesh,
129  patchFieldDict,
130  defaultPatchFieldType,
131  Zero
132  );
133  addPatchFields<volTensorField>
134  (
135  mesh,
136  patchFieldDict,
137  defaultPatchFieldType,
138  Zero
139  );
140 
141  // Surface fields
142 
143  addPatchFields<surfaceScalarField>
144  (
145  mesh,
146  patchFieldDict,
147  defaultPatchFieldType,
148  Zero
149  );
150  addPatchFields<surfaceVectorField>
151  (
152  mesh,
153  patchFieldDict,
154  defaultPatchFieldType,
155  Zero
156  );
157  addPatchFields<surfaceSphericalTensorField>
158  (
159  mesh,
160  patchFieldDict,
161  defaultPatchFieldType,
162  Zero
163  );
164  addPatchFields<surfaceSymmTensorField>
165  (
166  mesh,
167  patchFieldDict,
168  defaultPatchFieldType,
169  Zero
170  );
171  addPatchFields<surfaceTensorField>
172  (
173  mesh,
174  patchFieldDict,
175  defaultPatchFieldType,
176  Zero
177  );
178 
179  // Create reordering list
180  // patches before insert position stay as is
181  labelList oldToNew(sz+1);
182  for (label i = 0; i < insertPatchi; i++)
183  {
184  oldToNew[i] = i;
185  }
186  // patches after insert position move one up
187  for (label i = insertPatchi; i < sz; i++)
188  {
189  oldToNew[i] = i+1;
190  }
191  // appended patch gets moved to insert position
192  oldToNew[sz] = insertPatchi;
193 
194  // Shuffle into place
195  polyPatches.reorder(oldToNew, validBoundary);
196  fvPatches.reorder(oldToNew);
197 
198  reorderPatchFields<volScalarField>(mesh, oldToNew);
199  reorderPatchFields<volVectorField>(mesh, oldToNew);
200  reorderPatchFields<volSphericalTensorField>(mesh, oldToNew);
201  reorderPatchFields<volSymmTensorField>(mesh, oldToNew);
202  reorderPatchFields<volTensorField>(mesh, oldToNew);
203  reorderPatchFields<surfaceScalarField>(mesh, oldToNew);
204  reorderPatchFields<surfaceVectorField>(mesh, oldToNew);
205  reorderPatchFields<surfaceSphericalTensorField>(mesh, oldToNew);
206  reorderPatchFields<surfaceSymmTensorField>(mesh, oldToNew);
207  reorderPatchFields<surfaceTensorField>(mesh, oldToNew);
208 
209  return insertPatchi;
210 }
211 
212 
213 void Foam::fvMeshTools::setPatchFields
214 (
215  fvMesh& mesh,
216  const label patchi,
217  const dictionary& patchFieldDict
218 )
219 {
220  setPatchFields<volScalarField>(mesh, patchi, patchFieldDict);
221  setPatchFields<volVectorField>(mesh, patchi, patchFieldDict);
222  setPatchFields<volSphericalTensorField>(mesh, patchi, patchFieldDict);
223  setPatchFields<volSymmTensorField>(mesh, patchi, patchFieldDict);
224  setPatchFields<volTensorField>(mesh, patchi, patchFieldDict);
225  setPatchFields<surfaceScalarField>(mesh, patchi, patchFieldDict);
226  setPatchFields<surfaceVectorField>(mesh, patchi, patchFieldDict);
227  setPatchFields<surfaceSphericalTensorField>
228  (
229  mesh,
230  patchi,
231  patchFieldDict
232  );
233  setPatchFields<surfaceSymmTensorField>(mesh, patchi, patchFieldDict);
234  setPatchFields<surfaceTensorField>(mesh, patchi, patchFieldDict);
235 }
236 
237 
239 {
240  setPatchFields<volScalarField>(mesh, patchi, Zero);
241  setPatchFields<volVectorField>(mesh, patchi, Zero);
242  setPatchFields<volSphericalTensorField>
243  (
244  mesh,
245  patchi,
246  Zero
247  );
248  setPatchFields<volSymmTensorField>
249  (
250  mesh,
251  patchi,
252  Zero
253  );
254  setPatchFields<volTensorField>(mesh, patchi, Zero);
255  setPatchFields<surfaceScalarField>(mesh, patchi, Zero);
256  setPatchFields<surfaceVectorField>(mesh, patchi, Zero);
257  setPatchFields<surfaceSphericalTensorField>
258  (
259  mesh,
260  patchi,
261  Zero
262  );
263  setPatchFields<surfaceSymmTensorField>
264  (
265  mesh,
266  patchi,
267  Zero
268  );
269  setPatchFields<surfaceTensorField>(mesh, patchi, Zero);
270 }
271 
272 
273 // Deletes last patch
274 void Foam::fvMeshTools::trimPatches(fvMesh& mesh, const label nPatches)
275 {
276  // Clear local fields and e.g. polyMesh globalMeshData.
277  mesh.clearOut();
278 
279  polyBoundaryMesh& polyPatches =
280  const_cast<polyBoundaryMesh&>(mesh.boundaryMesh());
281  fvBoundaryMesh& fvPatches = const_cast<fvBoundaryMesh&>(mesh.boundary());
282 
283  if (polyPatches.empty())
284  {
286  << "No patches in mesh"
287  << abort(FatalError);
288  }
289 
290  label nFaces = 0;
291  for (label patchi = nPatches; patchi < polyPatches.size(); patchi++)
292  {
293  nFaces += polyPatches[patchi].size();
294  }
295  reduce(nFaces, sumOp<label>());
296 
297  if (nFaces)
298  {
300  << "There are still " << nFaces
301  << " faces in " << polyPatches.size()-nPatches
302  << " patches to be deleted" << abort(FatalError);
303  }
304 
305  // Remove actual patches
306  polyPatches.setSize(nPatches);
307  fvPatches.setSize(nPatches);
308 
309  trimPatchFields<volScalarField>(mesh, nPatches);
310  trimPatchFields<volVectorField>(mesh, nPatches);
311  trimPatchFields<volSphericalTensorField>(mesh, nPatches);
312  trimPatchFields<volSymmTensorField>(mesh, nPatches);
313  trimPatchFields<volTensorField>(mesh, nPatches);
314 
315  trimPatchFields<surfaceScalarField>(mesh, nPatches);
316  trimPatchFields<surfaceVectorField>(mesh, nPatches);
317  trimPatchFields<surfaceSphericalTensorField>(mesh, nPatches);
318  trimPatchFields<surfaceSymmTensorField>(mesh, nPatches);
319  trimPatchFields<surfaceTensorField>(mesh, nPatches);
320 }
321 
322 
324 (
325  fvMesh& mesh,
326  const labelList& oldToNew,
327  const label nNewPatches,
328  const bool validBoundary
329 )
330 {
331  polyBoundaryMesh& polyPatches =
332  const_cast<polyBoundaryMesh&>(mesh.boundaryMesh());
333  fvBoundaryMesh& fvPatches = const_cast<fvBoundaryMesh&>(mesh.boundary());
334 
335  // Shuffle into place
336  polyPatches.reorder(oldToNew, validBoundary);
337  fvPatches.reorder(oldToNew);
338 
339  reorderPatchFields<volScalarField>(mesh, oldToNew);
340  reorderPatchFields<volVectorField>(mesh, oldToNew);
341  reorderPatchFields<volSphericalTensorField>(mesh, oldToNew);
342  reorderPatchFields<volSymmTensorField>(mesh, oldToNew);
343  reorderPatchFields<volTensorField>(mesh, oldToNew);
344  reorderPatchFields<surfaceScalarField>(mesh, oldToNew);
345  reorderPatchFields<surfaceVectorField>(mesh, oldToNew);
346  reorderPatchFields<surfaceSphericalTensorField>(mesh, oldToNew);
347  reorderPatchFields<surfaceSymmTensorField>(mesh, oldToNew);
348  reorderPatchFields<surfaceTensorField>(mesh, oldToNew);
349 
350  // Remove last.
351  trimPatches(mesh, nNewPatches);
352 }
353 
354 
355 // ************************************************************************* //
label nPatches
Definition: readKivaGrid.H:402
void clearOut()
Clear all geometry and addressing.
Definition: fvMesh.C:230
static void reorderPatches(fvMesh &, const labelList &oldToNew, const label nPatches, const bool validBoundary)
Reorder and remove trailing patches. If validBoundary call is parallel.
Definition: fvMeshTools.C:324
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
virtual autoPtr< polyPatch > clone(const polyBoundaryMesh &bm) const
Construct and return a clone, resetting the boundary mesh.
Definition: polyPatch.H:219
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
void reorder(const labelUList &)
Reorders elements. Ordering does not have to be done in.
Definition: PtrList.C:197
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:421
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
static void zeroPatchFields(fvMesh &mesh, const label patchi)
Change patchField to zero on registered fields.
Definition: fvMeshTools.C:238
dynamicFvMesh & mesh
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:300
A class for handling words, derived from string.
Definition: word.H:59
static label addPatch(fvMesh &mesh, const polyPatch &patch, const dictionary &patchFieldDict, const word &defaultPatchFieldType, const bool validBoundary)
Add patch. Inserts patch before all processor patches.
Definition: fvMeshTools.C:32
bool set(const label) const
Is element set.
Definition: PtrListI.H:65
static const zero Zero
Definition: zero.H:91
const word & name() const
Return name.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
void setSize(const label)
Reset size of PtrList. If extending the PtrList, new entries are.
Definition: PtrList.C:131
Foam::polyBoundaryMesh.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
label nFaces() const
bool empty() const
Return true if the UPtrList is empty (ie, size() is zero)
Definition: UPtrListI.H:36
label patchi
Foam::fvBoundaryMesh.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
label findPatchID(const word &patchName) const
Find patch index given a name.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:545
void reorder(const labelUList &, const bool validBoundary)
Reorders patches. Ordering does not have to be done in.
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29