decomposeFaces.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | cfMesh: A library for mesh generation
4  \\ / O peration |
5  \\ / A nd | Author: Franjo Juretic (franjo.juretic@c-fields.com)
6  \\/ M anipulation | Copyright (C) Creative Fields, Ltd.
7 -------------------------------------------------------------------------------
8 License
9  This file is part of cfMesh.
10 
11  cfMesh is free software; you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by the
13  Free Software Foundation; either version 3 of the License, or (at your
14  option) any later version.
15 
16  cfMesh 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 cfMesh. If not, see <http://www.gnu.org/licenses/>.
23 
24 Description
25 
26 \*---------------------------------------------------------------------------*/
27 
28 
29 #include "decomposeFaces.H"
30 #include "boolList.H"
31 
32 # ifdef USE_OMP
33 #include <omp.h>
34 # endif
35 
36 //#define DEBUGDec
37 
38 # ifdef DEBUGDec
39 #include "polyMeshGenAddressing.H"
40 # endif
41 
42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43 
44 namespace Foam
45 {
46 
47 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
48 //- Constructor
50 :
51  mesh_(mesh),
52  newFacesForFace_(),
53  done_(false)
54 {}
55 
56 //- Destructor
58 {}
59 
60 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
61 
62 void decomposeFaces::decomposeMeshFaces(const boolList& decomposeFace)
63 {
64  done_ = false;
65 
69 
70  const label nIntFaces = mesh_.nInternalFaces();
71  label nFaces(0), nPoints = mesh_.points().size();
72  point p;
73 
75  forAll(decomposeFace, fI)
76  if( decomposeFace[fI] )
77  ++nFaces;
78 
79  points.setSize(nPoints + nFaces);
80 
81  polyMeshGenModifier meshModifier(mesh_);
82  faceListPMG& faces = meshModifier.facesAccess();
83 
84  if( decomposeFace.size() != faces.size() )
86  (
87  "void decomposeFaces::decomposeMeshFaces(const boolList&)"
88  ) << "Incorrect size of the decomposeFace list!" << abort(FatalError);
89 
90  nFaces = 0;
91  VRWGraph newFaces;
92 
93  //- decompose internal faces
94  for(label faceI=0;faceI<nIntFaces;++faceI)
95  {
96  const face& f = faces[faceI];
97 
98  if( decomposeFace[faceI] )
99  {
100  # ifdef DEBUGDec
101  Info << "Decomposing internal face " << faceI << " with nodes "
102  << f << endl;
103  # endif
104 
105  FixedList<label, 3> newF;
106 
107  forAll(f, pI)
108  {
109  newF[0] = f[pI];
110  newF[1] = f.nextLabel(pI);
111  newF[2] = nPoints;
112 
113  # ifdef DEBUGDec
114  Info << "Storing face " << newF << " with label "
115  << nFaces << endl;
116  # endif
117 
118  newFaces.appendList(newF);
119  newFacesForFace_.append(faceI, nFaces++);
120  }
121 
122  p = f.centre(points);
123  points[nPoints] = p;
124  ++nPoints;
125  }
126  else
127  {
128  # ifdef DEBUGDec
129  Info << "Storing internal face " << faceI << " with nodes "
130  << f << " as new face " << faceI << endl;
131  # endif
132 
133  newFaces.appendList(f);
134  newFacesForFace_.append(faceI, nFaces++);
135  }
136  }
137 
138  //- decompose boundary faces
139  PtrList<boundaryPatch>& boundaries = meshModifier.boundariesAccess();
140  forAll(boundaries, patchI)
141  {
142  const label start = boundaries[patchI].patchStart();
143  const label end = start + boundaries[patchI].patchSize();
144 
145  boundaries[patchI].patchStart() = nFaces;
146 
147  for(label bfI=start;bfI<end;++bfI)
148  {
149  const face& f = faces[bfI];
150  if( decomposeFace[bfI] )
151  {
152  # ifdef DEBUGDec
153  Info << "Decomposing boundary face " << bfI
154  << " with nodes " << f << endl;
155  # endif
156 
157  FixedList<label, 3> newF;
158 
159  forAll(f, pI)
160  {
161  newF[0] = f[pI];
162  newF[1] = f.nextLabel(pI);
163  newF[2] = nPoints;
164 
165  # ifdef DEBUGDec
166  Info << "Storing face " << newF << " with label "
167  << nFaces << endl;
168  # endif
169 
170  newFaces.appendList(newF);
171  newFacesForFace_.append(bfI, nFaces++);
172  }
173 
174  p = f.centre(points);
175  points[nPoints++] = p;
176  }
177  else
178  {
179  # ifdef DEBUGDec
180  Info << "Storing boundary face " << bfI << " in patch"
181  << patchI << " as new face " << bfI << endl;
182  # endif
183 
184  newFaces.appendList(f);
185  newFacesForFace_.append(bfI, nFaces++);
186  }
187  }
188 
189  boundaries[patchI].patchSize() =
190  nFaces - boundaries[patchI].patchStart();
191  }
192 
193  //- decompose processor faces
194  if( Pstream::parRun() )
195  {
196  PtrList<processorBoundaryPatch>& procBoundaries =
197  meshModifier.procBoundariesAccess();
198 
199  forAll(procBoundaries, patchI)
200  {
201  const label start = procBoundaries[patchI].patchStart();
202  const label end = start + procBoundaries[patchI].patchSize();
203 
204  const bool own = procBoundaries[patchI].owner();
205 
206  procBoundaries[patchI].patchStart() = nFaces;
207 
208  for(label bfI=start;bfI<end;++bfI)
209  {
210  face f;
211  if( own )
212  {
213  f = faces[bfI];
214  }
215  else
216  {
217  f = faces[bfI].reverseFace();
218  }
219 
220  if( decomposeFace[bfI] )
221  {
222  # ifdef DEBUGDec
223  Info << "Decomposing processor boundary face " << bfI
224  << " with nodes " << f << endl;
225  # endif
226 
227  face newF(3);
228 
229  forAll(f, pI)
230  {
231  newF[0] = f[pI];
232  newF[1] = f.nextLabel(pI);
233  newF[2] = nPoints;
234 
235  # ifdef DEBUGDec
236  Info << "Storing face " << newF << " with label "
237  << nFaces << endl;
238  # endif
239 
240  if( own )
241  {
242  newFaces.appendList(newF);
243  }
244  else
245  {
246  newFaces.appendList(newF.reverseFace());
247  }
248  newFacesForFace_.append(bfI, nFaces++);
249  }
250 
251  p = f.centre(points);
252  points[nPoints++] = p;
253  }
254  else
255  {
256  # ifdef DEBUGDec
257  Info << "Storing boundary face " << bfI << " in patch"
258  << patchI << " as new face " << bfI << endl;
259  # endif
260 
261  if( own )
262  {
263  newFaces.appendList(f);
264  }
265  else
266  {
267  newFaces.appendList(f.reverseFace());
268  }
269  newFacesForFace_.append(bfI, nFaces++);
270  }
271  }
272 
273  procBoundaries[patchI].patchSize() =
274  nFaces - procBoundaries[patchI].patchStart();
275  }
276  }
277 
278  //- store the faces back into their list
279  faces.setSize(nFaces);
280  forAll(faces, faceI)
281  {
282  face& f = faces[faceI];
283  f.setSize(newFaces.sizeOfRow(faceI));
284 
285  forAll(f, pI)
286  f[pI] = newFaces(faceI, pI);
287  }
288  newFaces.setSize(0);
289 
290  //- update subsets
292 
293  //- change the mesh
294  cellListPMG& cells = meshModifier.cellsAccess();
295 
296  # ifdef USE_OMP
297  # pragma omp parallel for schedule(dynamic, 40)
298  # endif
299  forAll(cells, cellI)
300  {
301  cell& c = cells[cellI];
302 
303  DynList<label> newC;
304 
305  forAll(c, fJ)
306  {
307  const label faceI = c[fJ];
308  forAllRow(newFacesForFace_, faceI, nfI)
309  newC.append(newFacesForFace_(faceI, nfI));
310  }
311 
312  # ifdef DEBUGDec
313  Info << "Cell " << cellI << " with faces " << c
314  << " is changed into " << newC << endl;
315  # endif
316 
317  c.setSize(newC.size());
318  forAll(newC, fJ)
319  c[fJ] = newC[fJ];
320  }
321 
322  meshModifier.clearAll();
323 
324  done_ = true;
325 
326  # ifdef DEBUGDec
327  Info << "New cells are " << cells << endl;
328  mesh_.write();
329  # endif
330 }
331 
333 (
334  const boolList& concaveVertex
335 )
336 {
337  if( Pstream::parRun() )
338  {
340  (
341  "void decomposeFaces::decomposeConcaveInternalFaces"
342  "(const boolList& concaveVertex)"
343  ) << "This procedure is not parallelised!" << exit(FatalError);
344  }
345 
346  done_ = false;
347 
348  newFacesForFace_.setSize(mesh_.faces().size());
349  forAll(newFacesForFace_, fI)
350  newFacesForFace_.setRowSize(fI, 0);
351 
352  const label nIntFaces = mesh_.nInternalFaces();
353 
354  polyMeshGenModifier meshModifier(mesh_);
355  pointFieldPMG& points = meshModifier.pointsAccess();
356  faceListPMG& faces = meshModifier.facesAccess();
357 
358  if( concaveVertex.size() != mesh_.points().size() )
360  (
361  "void decomposeFaces::decomposeMeshFaces(const boolList&)"
362  ) << "Incorrect size of the concaveVertex list!" << abort(FatalError);
363 
364  VRWGraph newFaces;
365  DynList<label> newF;
366  newF.setSize(3);
367 
368  # ifdef DEBUGDec
369  const label id = mesh_.addFaceSubset("decomposedFaces");
370  # endif
371 
372  //- decompose internal faces
373  for(label faceI=0;faceI<nIntFaces;++faceI)
374  {
375  const face& f = faces[faceI];
376 
377  DynList<label> concavePos;
378  forAll(f, pI)
379  if( concaveVertex[f[pI]] )
380  {
381  concavePos.append(pI);
382  }
383 
384  if( concavePos.size() == 1 )
385  {
386  # ifdef DEBUGDec
387  Info << "1. Decomposing internal face " << faceI << " with nodes "
388  << f << endl;
389  mesh_.addFaceToSubset(id, faceI);
390  # endif
391 
392  newF[0] = f[concavePos[0]];
393 
394  for(label pI=1;pI<(f.size()-1);++pI)
395  {
396  const label pJ = (concavePos[0] + pI) % f.size();
397  newF[1] = f[pJ];
398  newF[2] = f.nextLabel(pJ);
399 
400  # ifdef DEBUGDec
401  Info << "Storing face " << newF << " with label "
402  << newFaces.size() << endl;
403  # endif
404 
405  newFacesForFace_.append(faceI, newFaces.size());
406  newFaces.appendList(newF);
407  }
408  }
409  else if( concavePos.size() > 1 )
410  {
411  # ifdef DEBUGDec
412  Info << "2. Decomposing internal face " << faceI << " with nodes "
413  << f << endl;
414  mesh_.addFaceToSubset(id, faceI);
415  # endif
416 
417  newF[0] = points.size();
418  forAll(f, pI)
419  {
420  newF[1] = f[pI];
421  newF[2] = f.nextLabel(pI);
422 
423  # ifdef DEBUGDec
424  Info << "2. Storing face " << newF << " with label "
425  << newFaces.size() << endl;
426  # endif
427 
428  newFacesForFace_.append(faceI, newFaces.size());
429  newFaces.appendList(newF);
430  }
431 
432  const point fCent = f.centre(points);
433  points.append(fCent);
434  }
435  else
436  {
437  # ifdef DEBUGDec
438  Info << "Storing internal face " << faceI << " with nodes "
439  << f << " as new face " << newFaces.size() << endl;
440  # endif
441 
442  newFacesForFace_.append(faceI, newFaces.size());
443  newFaces.appendList(f);
444  }
445  }
446 
447  PtrList<boundaryPatch>& boundaries = meshModifier.boundariesAccess();
448  forAll(boundaries, patchI)
449  {
450  const label start = boundaries[patchI].patchStart();
451  const label end = start + boundaries[patchI].patchSize();
452 
453  //- set new patch start
454  boundaries[patchI].patchStart() = newFaces.size();
455 
456  //- store faces into newFaces
457  for(label bfI=start;bfI<end;++bfI)
458  {
459  newFacesForFace_.append(bfI, newFaces.size());
460  newFaces.appendList(faces[bfI]);
461  }
462  }
463 
464  //- copy new faces into the faceListPMG
465  faces.setSize(newFaces.size());
466  forAll(newFaces, faceI)
467  {
468  faces[faceI].setSize(newFaces.sizeOfRow(faceI));
469 
470  forAllRow(newFaces, faceI, pI)
471  faces[faceI][pI] = newFaces(faceI, pI);
472  }
473 
474  newFaces.setSize(0);
475 
476  //- update cells
477  cellListPMG& cells = meshModifier.cellsAccess();
478 
479  # ifdef USE_OMP
480  # pragma omp parallel for schedule(dynamic, 40)
481  # endif
482  forAll(cells, cellI)
483  {
484  cell& c = cells[cellI];
485 
486  DynList<label, 24> newC;
487 
488  forAll(c, fJ)
489  {
490  const label faceI = c[fJ];
491  forAllRow(newFacesForFace_, faceI, nfI)
492  newC.append(newFacesForFace_(faceI, nfI));
493  }
494 
495  # ifdef DEBUGDec
496  Info << "Cell " << cellI << " with faces " << c
497  << " is changed into " << newC << endl;
498  # endif
499 
500  c.setSize(newC.size());
501  forAll(newC, fJ)
502  c[fJ] = newC[fJ];
503  }
504 
505  meshModifier.clearAll();
506 
507  //- update subsets
508  mesh_.updateFaceSubsets(newFacesForFace_);
509 
510  # ifdef DEBUGDec
511  Info << "New cells are " << cells << endl;
512  mesh_.write();
513  ::exit(1);
514  # endif
515 
516  meshModifier.removeUnusedVertices();
517 
518  done_ = true;
519 }
520 
522 {
523  if( !done_ )
524  WarningIn
525  (
526  "const VRWGraph& decomposeFaces::newFacesForFace() const"
527  ) << "Decomposition is not yet performed!" << endl;
528 
529  return newFacesForFace_;
530 }
531 
532 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *//
533 
534 } // End namespace Foam
535 
536 // ************************************************************************* //
Foam::decomposeFaces::~decomposeFaces
~decomposeFaces()
Destructor.
Definition: decomposeFaces.C:57
Foam::face::reverseFace
face reverseFace() const
Return face with reverse direction.
Definition: face.C:611
boolList.H
p
p
Definition: pEqn.H:62
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:406
Foam::decomposeFaces::newFacesForFace
const VRWGraph & newFacesForFace() const
Definition: decomposeFaces.C:521
Foam::UPstream::parRun
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:377
decomposeFaces.H
Foam::polyMeshGen
Definition: polyMeshGen.H:46
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Foam::polyMeshGenPoints::points
const pointFieldPMG & points() const
access to points
Definition: polyMeshGenPointsI.H:44
Foam::cellListPMG
Definition: cellListPMG.H:49
Foam::polyMeshGenFaces::updateFaceSubsets
void updateFaceSubsets(const ListType &)
Definition: polyMeshGenFacesI.H:194
Foam::polyMeshGenFaces::faces
const faceListPMG & faces() const
access to faces
Definition: polyMeshGenFacesI.H:43
Foam::polyMeshGenFaces::nInternalFaces
label nInternalFaces() const
return number of internal faces
Definition: polyMeshGenFacesI.H:48
nPoints
label nPoints
Definition: gmvOutputHeader.H:2
Foam::polyMeshGenModifier::boundariesAccess
PtrList< boundaryPatch > & boundariesAccess()
access to boundary data
Definition: polyMeshGenModifier.H:131
Foam::polyMeshGenModifier::procBoundariesAccess
PtrList< processorBoundaryPatch > & procBoundariesAccess()
access to processor boundary data
Definition: polyMeshGenModifier.H:125
Foam::polyMeshGenModifier::cellsAccess
cellListPMG & cellsAccess()
access to cells
Definition: polyMeshGenModifier.H:119
Foam::label
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
Foam::Ostream::write
virtual Ostream & write(const token &)=0
Write next token to stream.
Foam::Info
messageStream Info
Foam::decomposeFaces::decomposeConcaveInternalFaces
void decomposeConcaveInternalFaces(const boolList &concaveVertex)
decompose internal faces containing concave nodes
Definition: decomposeFaces.C:333
Foam::VRWGraph::size
label size() const
Returns the number of rows.
Definition: VRWGraphI.H:122
Foam::faceListPMG::setSize
void setSize(const label nElmts)
set the number of used elements
Definition: faceListPMGI.H:78
Foam::VRWGraph::setSize
void setSize(const label)
Reset the number of rows.
Definition: VRWGraphI.H:132
Foam::PtrList
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:61
forAllRow
#define forAllRow(graph, rowI, index)
Definition: VRWGraph.H:277
Foam::pointFieldPMG::size
label size() const
return the number of used elements
Definition: pointFieldPMGI.H:71
Foam::polyMeshGenModifier::removeUnusedVertices
void removeUnusedVertices()
remove unused vertices
Definition: polyMeshGenModifierRemoveUnusedVertices.C:37
Foam::VRWGraph::sizeOfRow
label sizeOfRow(const label rowI) const
Returns the number of elements in the given row.
Definition: VRWGraphI.H:127
Foam::FatalError
error FatalError
Foam::polyMeshGenModifier::pointsAccess
pointFieldPMG & pointsAccess()
access to mesh points
Definition: polyMeshGenModifier.H:107
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:18
Foam::decomposeFaces::done_
bool done_
is decomposition performed
Definition: decomposeFaces.H:57
Foam
Namespace for OpenFOAM.
Definition: combustionModel.C:30
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:131
Foam::DynList< label >
Foam::decomposeFaces::decomposeMeshFaces
void decomposeMeshFaces(const boolList &decomposeFace)
decompose selected faces into triangles using midnode subdivision
Definition: decomposeFaces.C:62
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
Foam::polyMeshGenModifier
Definition: polyMeshGenModifier.H:52
Foam::decomposeFaces::mesh_
polyMeshGen & mesh_
reference to the mesh
Definition: decomposeFaces.H:51
f
labelList f(nPoints)
Foam::faceListPMG::size
label size() const
return the number of used elements
Definition: faceListPMGI.H:73
Foam::VRWGraph::appendList
void appendList(const ListType &l)
Append a list as a row at the end of the graph.
Definition: VRWGraphI.H:286
Foam::Vector< scalar >
Foam::List
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
Foam::polyMeshGen::write
void write() const
Definition: polyMeshGen.C:126
Foam::DynList::setSize
void setSize(const label)
Reset size of List.
Definition: DynListI.H:263
Foam::FixedList< label, 3 >
points
const pointField & points
Definition: gmvOutputHeader.H:1
Foam::VRWGraph::append
void append(const label rowI, const label)
Append an element to the given row.
Definition: VRWGraphI.H:303
Foam::polyMeshGenModifier::clearAll
void clearAll()
clear out all allocated data
Definition: polyMeshGenModifier.H:200
Foam::constant::universal::c
const dimensionedScalar c
Speed of light in a vacuum.
Foam::decomposeFaces::newFacesForFace_
VRWGraph newFacesForFace_
number of points
Definition: decomposeFaces.H:54
FatalErrorIn
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:313
Foam::face
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
Foam::DynList::size
label size() const
Definition: DynListI.H:235
cells
const cellShapeList & cells
Definition: gmvOutputHeader.H:3
WarningIn
#define WarningIn(functionName)
Report a warning using Foam::Warning.
Definition: messageStream.H:254
Foam::polyMeshGenModifier::facesAccess
faceListPMG & facesAccess()
access to mesh faces
Definition: polyMeshGenModifier.H:113
Foam::List::size
void size(const label)
Override size to be inconsistent with allocated storage.
Foam::faceListPMG
Definition: faceListPMG.H:50
Foam::pointFieldPMG
Definition: pointFieldPMG.H:50
Foam::VRWGraph
Definition: VRWGraph.H:101
polyMeshGenAddressing.H
Foam::cell
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:56
Foam::VRWGraph::setRowSize
void setRowSize(const label rowI, const label newSize)
Reset the size of the given row.
Definition: VRWGraphI.H:204
Foam::DynList::append
void append(const T &e)
Append an element at the end of the list.
Definition: DynListI.H:304
Foam::decomposeFaces::decomposeFaces
decomposeFaces(const decomposeFaces &)
copy constructor