createPolyBoundary.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2023 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 Description
25  boundary faces
26  - use pointCells when searching for connectivity
27  - initialise the cell connectivity with '-1'
28  - find both cell faces corresponding to the baffles and mark them
29  to prevent a connection
30  - standard connectivity checks
31 
32  - added baffle and monitoring support
33 
34 \*---------------------------------------------------------------------------*/
35 
36 #include "meshReader.H"
37 #include "Time.H"
38 #include "polyPatch.H"
39 #include "emptyPolyPatch.H"
40 #include "preservePatchTypes.H"
41 
42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43 
44 void Foam::meshReader::addPolyBoundaryFace
45 (
46  const label cellId,
47  const label cellFaceId,
48  const label nCreatedFaces
49 )
50 {
51 #ifdef DEBUG_BOUNDARY
52  Info<< nCreatedFaces
53  << " add bnd for cell " << cellId
54  << " face " << cellFaceId
55  << " (original cell " << origCellId_[cellId] << ")"
56  << endl;
57 #endif
58 
59  // standard case: volume cells
60  const face& thisFace = cellFaces_[cellId][cellFaceId];
61 
62  // Debugging
63  if (cellPolys_[cellId][cellFaceId] > nInternalFaces_)
64  {
66  << "Problem with face: " << thisFace << endl
67  << "Probably multiple definitions "
68  << "of a single boundary face." << endl
69  << endl;
70  }
71  else if (cellPolys_[cellId][cellFaceId] >= 0)
72  {
74  << "Problem with face: " << thisFace << endl
75  << "Probably trying to define a boundary face "
76  << "on a previously matched internal face." << endl
77  << "Internal face: "
78  << meshFaces_[cellPolys_[cellId][cellFaceId]]
79  << endl;
80  }
81 
82  meshFaces_[nCreatedFaces] = thisFace;
83  cellPolys_[cellId][cellFaceId] = nCreatedFaces;
84 }
85 
86 
87 void Foam::meshReader::addPolyBoundaryFace
88 (
89  const cellFaceIdentifier& identifier,
90  const label nCreatedFaces
91 )
92 {
93  addPolyBoundaryFace(identifier.cell, identifier.face, nCreatedFaces);
94 }
95 
96 
97 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
98 
99 void Foam::meshReader::createPolyBoundary()
100 {
101  label nBoundaryFaces = 0;
102  label nMissingFaces = 0;
103  label nInterfaces = 0;
104 
105  const faceListList& cFaces = cellFaces();
106 
107  // determine number of non-patched faces:
108  forAll(cellPolys_, celli)
109  {
110  cell& curCell = cellPolys_[celli];
111 
112  forAll(curCell, fI)
113  {
114  if (curCell[fI] < 0)
115  {
116  nMissingFaces++;
117  }
118  }
119  }
120 
121  forAll(boundaryIds_, patchi)
122  {
123  nBoundaryFaces += boundaryIds_[patchi].size();
124  }
125 
126  Info<< nl
127  << "There are " << nMissingFaces
128  << " faces to be patched and " << nBoundaryFaces
129  << " specified - collect missed boundaries to final patch" << endl;
130 
131  patchStarts_.setSize(boundaryIds_.size());
132  patchSizes_.setSize(boundaryIds_.size());
133 
134  label nCreatedFaces = nInternalFaces_;
135  label baffleOffset = cFaces.size();
136  interfaces_.setSize(baffleIds_.size());
137  nBoundaryFaces = 0;
138 
139  forAll(boundaryIds_, patchi)
140  {
141  const List<cellFaceIdentifier>& idList = boundaryIds_[patchi];
142 
143  patchStarts_[patchi] = nCreatedFaces;
144 
145  // write each baffle side separately
146  if (patchPhysicalTypes_[patchi] == "baffle")
147  {
148  label count = 0;
149 
150  for (label side = 0; side < 2; ++side)
151  {
152  label position = nInterfaces;
153 
154  forAll(idList, bndI)
155  {
156  label baffleI = idList[bndI].cell - baffleOffset;
157 
158  if
159  (
160  baffleI >= 0
161  && baffleI < baffleFaces_.size()
162  && baffleIds_[baffleI].size()
163  )
164  {
165  addPolyBoundaryFace
166  (
167  baffleIds_[baffleI][side],
168  nCreatedFaces
169  );
170 
171  // remove applied boundaries (2nd pass)
172  if (side == 1)
173  {
174  baffleIds_[baffleI].clear();
175  }
176 
177  interfaces_[position][side] = nCreatedFaces;
178 
179  nBoundaryFaces++;
180  nCreatedFaces++;
181  position++;
182  count++;
183  }
184  }
185  }
186 
187  nInterfaces += (count - (count % 2)) / 2;
188  }
189  else if (patchPhysicalTypes_[patchi] == "monitoring")
190  {
191  // translate the "monitoring" pseudo-boundaries to face sets
192  List<label> monitoring(idList.size());
193 
194  label monitorI = 0;
195  forAll(idList, bndI)
196  {
197  label cellId = idList[bndI].cell;
198  label faceId = idList[bndI].face;
199 
200  // standard case: volume cells
201  if (cellId < baffleOffset)
202  {
203  label faceNr = cellPolys_[cellId][faceId];
204  if (faceNr >= 0)
205  {
206  monitoring[monitorI++] = faceNr;
207  }
208  }
209  }
210 
211  monitoringSets_.insert(patchNames_[patchi], monitoring);
212  }
213  else
214  {
215  forAll(idList, bndI)
216  {
217  // standard case: volume cells
218  if (idList[bndI].cell < baffleOffset)
219  {
220  addPolyBoundaryFace
221  (
222  idList[bndI],
223  nCreatedFaces
224  );
225 
226  nBoundaryFaces++;
227  nCreatedFaces++;
228  }
229  }
230  }
231 
232  patchSizes_[patchi] = nCreatedFaces - patchStarts_[patchi];
233  }
234 
235  // add in missing faces
236  Info<< "Missing faces added to patch after face "
237  << nCreatedFaces << ":" <<endl;
238  nMissingFaces = 0;
239 
240  // look for baffles first - keep them together at the start of the patch
241  for (label side = 0; side < 2; ++side)
242  {
243  label position = nInterfaces;
244 
245  forAll(baffleIds_, baffleI)
246  {
247  if (baffleIds_[baffleI].size())
248  {
249  // add each side for each baffle
250  addPolyBoundaryFace
251  (
252  baffleIds_[baffleI][side],
253  nCreatedFaces
254  );
255 
256  interfaces_[position][side] = nCreatedFaces;
257 
258  // remove applied boundaries (2nd pass)
259  if (side == 1)
260  {
261  baffleIds_[baffleI].clear();
262  }
263 
264  nMissingFaces++;
265  nCreatedFaces++;
266  position++;
267  }
268  }
269  }
270 
271  nInterfaces += (nMissingFaces - (nMissingFaces % 2)) / 2;
272 
273  // scan for any other missing faces
274  forAll(cellPolys_, celli)
275  {
276  const labelList& curFaces = cellPolys_[celli];
277 
278  forAll(curFaces, cellFacei)
279  {
280  if (curFaces[cellFacei] < 0)
281  {
282  // just report the first few
283  if (nMissingFaces < 4)
284  {
285  const face& thisFace = cFaces[celli][cellFacei];
286 
287  Info<< " cell " << celli << " face " << cellFacei
288  << " (original cell " << origCellId_[celli] << ")"
289  << " face: " << thisFace
290  << endl;
291  }
292  else if (nMissingFaces == 5)
293  {
294  Info<< " ..." << nl << endl;
295  }
296 
297  addPolyBoundaryFace(celli, cellFacei, nCreatedFaces);
298  nMissingFaces++;
299  nCreatedFaces++;
300  }
301  }
302  }
303 
304  Info<< "Added " << nMissingFaces << " unmatched faces" << endl;
305 
306  // Add missing faces to last patch ('Default_Empty' etc.)
307  if (nMissingFaces > 0)
308  {
309  patchSizes_.last() = nMissingFaces;
310  }
311 
312 
313  // reset the size of the face list
314  meshFaces_.setSize(nCreatedFaces);
315 
316  // check the mesh for face mismatch
317  // (faces addressed once or more than twice)
318  labelList markupFaces(meshFaces_.size(), 0);
319 
320  forAll(cellPolys_, celli)
321  {
322  const labelList& curFaces = cellPolys_[celli];
323 
324  forAll(curFaces, facei)
325  {
326  markupFaces[curFaces[facei]]++;
327  }
328  }
329 
330  for (label i = nInternalFaces_; i < markupFaces.size(); i++)
331  {
332  markupFaces[i]++;
333  }
334 
335  label nProblemFaces = 0;
336 
337  forAll(markupFaces, facei)
338  {
339  if (markupFaces[facei] != 2)
340  {
341  const face& problemFace = meshFaces_[facei];
342 
344  << "Problem with face " << facei << ": addressed "
345  << markupFaces[facei] << " times (should be 2!). Face: "
346  << problemFace << endl;
347 
348  nProblemFaces++;
349  }
350  }
351 
352  if (nProblemFaces > 0)
353  {
354  Info<< "Number of incorrectly matched faces: "
355  << nProblemFaces << endl;
356  }
357 
358  // adjust for missing members
359  if (nInterfaces < interfaces_.size())
360  {
361  interfaces_.setSize(nInterfaces);
362  }
363 
364  Info<< "Number of boundary faces: " << nBoundaryFaces << nl
365  << "Total number of faces: " << nCreatedFaces << nl
366  << "Number of interfaces: " << nInterfaces << endl;
367 }
368 
369 
370 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
371 
373 Foam::meshReader::polyBoundaryPatches(const polyMesh& mesh)
374 {
375  label nUsed = 0, nEmpty = 0;
376  label nPatches = patchStarts_.size();
377 
378  // avoid empty patches - move to the end of the lists and truncate
379  labelList oldToNew = identityMap(nPatches);
380  forAll(patchSizes_, patchi)
381  {
382  if (patchSizes_[patchi] > 0)
383  {
384  oldToNew[patchi] = nUsed++;
385  }
386  else
387  {
388  nEmpty++;
389  oldToNew[patchi] = nPatches - nEmpty;
390  }
391  }
392 
393  nPatches = nUsed;
394 
395  if (nEmpty)
396  {
397  Info<< "Removing " << nEmpty << " empty patches" << endl;
398 
399  inplaceReorder(oldToNew, patchTypes_);
400  inplaceReorder(oldToNew, patchNames_);
401  inplaceReorder(oldToNew, patchStarts_);
402  inplaceReorder(oldToNew, patchSizes_);
403  }
404 
405  patchTypes_.setSize(nPatches);
406  patchNames_.setSize(nPatches);
407  patchStarts_.setSize(nPatches);
408  patchSizes_.setSize(nPatches);
409 
410 
411  List<polyPatch*> p(nPatches);
412 
413  // All patch dictionaries
414  PtrList<dictionary> patchDicts(patchNames_.size());
415  // Default boundary patch types
416  word defaultFacesType(emptyPolyPatch::typeName);
417 
418  // we could consider dropping this entirely
420  (
421  mesh,
422  mesh.instance(),
423  mesh.meshDir(),
424  patchNames_,
425  patchDicts,
426  "defaultFaces",
428  );
430  {
431  if (!patchDicts.set(patchi))
432  {
433  patchDicts.set(patchi, new dictionary());
434  }
435  dictionary& patchDict = patchDicts[patchi];
436 
437  // add but not overwrite type
438  patchDict.add("type", patchTypes_[patchi], false);
439  if (patchPhysicalTypes_.size() && patchPhysicalTypes_[patchi].size())
440  {
441  patchDict.add("startFace", patchPhysicalTypes_[patchi], false);
442  }
443 
444  // overwrite sizes and start
445  patchDict.add("nFaces", patchSizes_[patchi], true);
446  patchDict.add("startFace", patchStarts_[patchi], true);
447  }
448 
449 
450  forAll(patchStarts_, patchi)
451  {
453  (
454  patchNames_[patchi],
456  patchi,
457  mesh.boundaryMesh()
458  ).ptr();
459  }
460 
461  return p;
462 }
463 
464 
465 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:91
labelList origCellId_
Lookup original Cell number for a given cell.
Definition: meshReader.H:252
faceListList cellFaces_
List of faces for every cell.
Definition: meshReader.H:268
static autoPtr< polyPatch > New(const word &patchType, const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm)
Return a pointer to a new patch created on freestore from.
Definition: polyPatchNew.C:32
label patchi
const label cellId
label faceId(-1)
#define InfoInFunction
Report an information message using Foam::Info.
List< label > labelList
A List of labels.
Definition: labelList.H:56
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
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
messageStream Info
List< faceList > faceListList
Definition: faceListFwd.H:43
void preservePatchTypes(const objectRegistry &obr, const word &meshInstance, const fileName &meshDir, const wordList &patchNames, PtrList< dictionary > &patchDicts, const word &defaultFacesName, word &defaultFacesType)
Preserve patch types.
label count(const ListType &l, typename ListType::const_reference x)
Count the number of occurrences of a value in a list.
labelList identityMap(const label len)
Create identity map (map[i] == i) of given length.
Definition: ListOps.C:104
static const char nl
Definition: Ostream.H:266
void inplaceReorder(const labelUList &oldToNew, ListType &)
Inplace reorder the elements of a list.
preservePatchTypes
PtrList< dictionary > patchDicts
Definition: readKivaGrid.H:531
word defaultFacesType
Definition: readKivaGrid.H:455
label nPatches
Definition: readKivaGrid.H:396
volScalarField & p