partTetMeshAddressing.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 #include "demandDrivenData.H"
29 #include "polyMeshGenModifier.H"
30 #include "partTetMesh.H"
31 #include "tetMatcher.H"
32 #include "polyMeshGenAddressing.H"
33 #include "helperFunctionsPar.H"
34 
35 //#define DEBUGSmooth
36 
37 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
38 
39 namespace Foam
40 {
41 
42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43 
44 void partTetMesh::createPointsAndTets
45 (
46  const List<direction>& useCell,
47  const boolList& lockedPoints
48 )
49 {
50  const pointFieldPMG& points = origMesh_.points();
51  const faceListPMG& faces = origMesh_.faces();
52  const cellListPMG& cells = origMesh_.cells();
53  const labelList& owner = origMesh_.owner();
54  const labelList& neighbour = origMesh_.neighbour();
55  const PtrList<boundaryPatch>& boundaries = origMesh_.boundaries();
56  const PtrList<processorBoundaryPatch>& procBoundaries =
57  origMesh_.procBoundaries();
58  const label nInternalFaces = origMesh_.nInternalFaces();
59 
60  //- check how many neighbours of a face are marked for smoothing
61  labelList usedFace(faces.size(), 0);
62 
63  //- mark faces
64  forAll(faces, faceI)
65  {
66  if( useCell[owner[faceI]] )
67  ++usedFace[faceI];
68 
69  if( neighbour[faceI] < 0 )
70  continue;
71 
72  if( useCell[neighbour[faceI]] )
73  ++usedFace[faceI];
74  }
75 
76  //- send data at processor boundaries
77  forAll(procBoundaries, patchI)
78  {
79  const label start = procBoundaries[patchI].patchStart();
80  const label size = procBoundaries[patchI].patchSize();
81 
82  labelLongList dataToSend;
83  for(label faceI=0;faceI<size;++faceI)
84  {
85  if( usedFace[start+faceI] )
86  dataToSend.append(faceI);
87  }
88 
89  OPstream toOtherProc
90  (
91  Pstream::blocking,
92  procBoundaries[patchI].neiProcNo(),
93  dataToSend.byteSize()
94  );
95 
96  toOtherProc << dataToSend;
97  }
98 
99  //- receive data at proc boundaries
100  forAll(procBoundaries, patchI)
101  {
102  labelLongList receivedData;
103 
104  IPstream fromOtherProc
105  (
106  Pstream::blocking,
107  procBoundaries[patchI].neiProcNo()
108  );
109 
110  fromOtherProc >> receivedData;
111 
112  const label start = procBoundaries[patchI].patchStart();
113  forAll(receivedData, faceI)
114  ++usedFace[start+receivedData[faceI]];
115  }
116 
117  const vectorField& faceCentres = origMesh_.addressingData().faceCentres();
118  const vectorField& cellCentres = origMesh_.addressingData().cellCentres();
119 
120  labelLongList nodeLabelForPoint(points.size(), -1);
121  labelLongList nodeLabelForFace(faces.size(), -1);
122  labelLongList nodeLabelForCell(cells.size(), -1);
123 
124  points_.clear();
125  smoothVertex_.clear();
126 
127  //- create BOUNDARY points
128  forAll(boundaries, patchI)
129  {
130  const boundaryPatch& patch = boundaries[patchI];
131  const label start = patch.patchStart();
132  const label end = start + patch.patchSize();
133 
134  for(label faceI=start;faceI<end;++faceI)
135  {
136  if( !usedFace[faceI] )
137  continue;
138 
139  const face& f = faces[faceI];
140 
141  if( f.size() > 3 )
142  {
143  //- create face centre
144  nodeLabelForFace[faceI] = points_.size();
145  points_.append(faceCentres[faceI]);
146  smoothVertex_.append(FACECENTRE);
147  }
148 
149  //- add face corners
150  forAll(f, pI)
151  {
152  const label pointI = f[pI];
153  if( nodeLabelForPoint[pointI] == -1 )
154  {
155  nodeLabelForPoint[pointI] = points_.size();
156  points_.append(points[pointI]);
157 
158  smoothVertex_.append(BOUNDARY);
159  }
160  }
161  }
162  }
163 
164  //- create points at processor boundaries
165  forAll(procBoundaries, patchI)
166  {
167  const processorBoundaryPatch& patch = procBoundaries[patchI];
168  const label start = patch.patchStart();
169  const label end = start + patch.patchSize();
170 
171  for(label faceI=start;faceI<end;++faceI)
172  {
173  if( !usedFace[faceI] )
174  continue;
175 
176  const face& f = faces[faceI];
177 
178  if( f.size() > 3 )
179  {
180  //- create face centre
181  nodeLabelForFace[faceI] = points_.size();
182  points_.append(faceCentres[faceI]);
183  smoothVertex_.append(FACECENTRE);
184  }
185 
186  //- add face corners
187  const direction vType = usedFace[faceI]==2?SMOOTH:NONE;
188  forAll(f, pI)
189  {
190  const label pointI = f[pI];
191  if( nodeLabelForPoint[pointI] == -1 )
192  {
193  nodeLabelForPoint[pointI] = points_.size();
194  points_.append(points[pointI]);
195 
196  smoothVertex_.append(vType);
197  }
198  else if( vType == NONE )
199  {
200  smoothVertex_[nodeLabelForPoint[pointI]] = NONE;
201  }
202  }
203  }
204  }
205 
206  //- create points for internal faces
207  for(label faceI=0;faceI<nInternalFaces;++faceI)
208  {
209  if( usedFace[faceI] )
210  {
211  const face& f = faces[faceI];
212 
213  if( f.size() > 3 )
214  {
215  //- create face centre
216  nodeLabelForFace[faceI] = points_.size();
217  points_.append(faceCentres[faceI]);
218  smoothVertex_.append(FACECENTRE);
219  }
220 
221  //- add face corners
222  const direction vType = usedFace[faceI]==2?SMOOTH:NONE;
223  forAll(f, pI)
224  {
225  const label pointI = f[pI];
226  if( nodeLabelForPoint[pointI] == -1 )
227  {
228  nodeLabelForPoint[pointI] = points_.size();
229  points_.append(points[pointI]);
230 
231  smoothVertex_.append(vType);
232  }
233  else if( vType == NONE )
234  {
235  smoothVertex_[nodeLabelForPoint[pointI]] = NONE;
236  }
237  }
238  }
239  }
240 
241  //- create tets
242  tetMatcher tet;
243  forAll(useCell, cI)
244  if( useCell[cI] )
245  {
246  const cell& c = cells[cI];
247 
248  if( tet.matchShape(false, faces, owner, cI, cells[cI]) )
249  {
250  const labelList& tVrt = tet.vertLabels();
251 
252  //- add tet
253  tets_.append
254  (
255  partTet
256  (
257  nodeLabelForPoint[tVrt[0]],
258  nodeLabelForPoint[tVrt[1]],
259  nodeLabelForPoint[tVrt[2]],
260  nodeLabelForPoint[tVrt[3]]
261  )
262  );
263 
264  continue;
265  }
266 
267  nodeLabelForCell[cI] = points_.size();
268  const label centreLabel = points_.size();
269  points_.append(cellCentres[cI]);
270  smoothVertex_.append(CELLCENTRE);
271 
272  forAll(c, fI)
273  {
274  const face& f = faces[c[fI]];
275 
276  if( owner[c[fI]] == cI )
277  {
278  if( f.size() == 3 )
279  {
280  partTet tet
281  (
282  nodeLabelForPoint[f[0]],
283  nodeLabelForPoint[f[2]],
284  nodeLabelForPoint[f[1]],
285  centreLabel
286  );
287 
288  # ifdef DEBUGSmooth
289  Info << "1.1 Tet " << tets_.size() << " is "
290  << tet << endl;
291  # endif
292 
293  tets_.append(tet);
294  }
295  else
296  {
297  forAll(f, pI)
298  {
299  partTet tet
300  (
301  nodeLabelForPoint[f[pI]],
302  nodeLabelForPoint[f.prevLabel(pI)],
303  nodeLabelForFace[c[fI]],
304  centreLabel
305  );
306 
307  # ifdef DEBUGSmooth
308  Info << "1.2 Tet " << tets_.size() << " is "
309  << tet << endl;
310  # endif
311 
312  tets_.append(tet);
313  }
314  }
315  }
316  else
317  {
318  if( f.size() == 3 )
319  {
320  partTet tet
321  (
322  nodeLabelForPoint[f[0]],
323  nodeLabelForPoint[f[1]],
324  nodeLabelForPoint[f[2]],
325  centreLabel
326  );
327 
328  # ifdef DEBUGSmooth
329  Info << "2.1 Tet " << tets_.size() << " is "
330  << tet << endl;
331  # endif
332 
333  tets_.append(tet);
334  }
335  else
336  {
337  forAll(f, pI)
338  {
339  partTet tet
340  (
341  nodeLabelForPoint[f[pI]],
342  nodeLabelForPoint[f.nextLabel(pI)],
343  nodeLabelForFace[c[fI]],
344  centreLabel
345  );
346 
347  # ifdef DEBUGSmooth
348  Info << "2.2 Tet " << tets_.size() << " is "
349  << tet << endl;
350  # endif
351 
352  tets_.append(tet);
353  }
354  }
355  }
356  }
357  }
358 
359  //- create node labels in origMesh_
360  nodeLabelInOrigMesh_.setSize(points_.size());
361  nodeLabelInOrigMesh_ = -1;
362  forAll(nodeLabelForPoint, pI)
363  if( nodeLabelForPoint[pI] != -1 )
364  {
365  //- lock mesh vertices
366  if( lockedPoints[pI] )
367  smoothVertex_[nodeLabelForPoint[pI]] |= LOCKED;
368 
369  nodeLabelInOrigMesh_[nodeLabelForPoint[pI]] = pI;
370  }
371 
372 
373  //- create pointTets_
374  pointTets_.reverseAddressing(points_.size(), tets_);
375 
376  //- create addressing for parallel runs
377  if( Pstream::parRun() )
378  {
379  createParallelAddressing
380  (
381  nodeLabelForPoint,
382  nodeLabelForFace,
383  nodeLabelForCell
384  );
385 
386  createBufferLayers();
387  }
388 
389  # ifdef DEBUGSmooth
390  forAll(nodeLabelInOrigMesh_, pI)
391  if(
392  (nodeLabelInOrigMesh_[pI] != -1) &&
393  (mag(points_[pI] - points[nodeLabelInOrigMesh_[pI]]) > SMALL)
394  )
395  FatalErrorIn
396  (
397  "void partTetMesh::createPointsAndTets"
398  "(const boolList& useCell)"
399  ) << "Node " << pI << " is dislocated" << abort(FatalError);
400  # endif
401 }
402 
403 void partTetMesh::createSMOOTHPointsOrdering() const
404 {
405  internalPointsOrderPtr_ = new VRWGraph();
406  VRWGraph& internalPointsOrder = *internalPointsOrderPtr_;
407 
408  internalPointsOrder.setSize(0);
409  labelLongList order(points_.size(), -1);
410  boolList helper(points_.size());
411 
412  bool found;
413  do
414  {
415  found = false;
416  helper = false;
417  labelLongList selectedPoints;
418 
419  forAll(points_, nodeI)
420  {
421  if( smoothVertex_[nodeI] & SMOOTH )
422  {
423  if( helper[nodeI] )
424  continue;
425  if( order[nodeI] != -1 )
426  continue;
427 
428  //- find neighbouring FACECENTRE and CELLCENTRE points
429  DynList<label, 64> neiCentrePoints, neiSmoothPoints;
430  forAllRow(pointTets_, nodeI, ptI)
431  {
432  const partTet& tet = tets_[pointTets_(nodeI, ptI)];
433 
434  for(label i=0;i<4;++i)
435  if( smoothVertex_[tet[i]] & (FACECENTRE+CELLCENTRE) )
436  {
437  neiCentrePoints.appendIfNotIn(tet[i]);
438  }
439  else if( smoothVertex_[tet[i]] & SMOOTH )
440  {
441  neiSmoothPoints.appendIfNotIn(tet[i]);
442  }
443  }
444 
445  //- find neighbouring SMOOTH points
446  forAll(neiCentrePoints, ncI)
447  {
448  const label centreI = neiCentrePoints[ncI];
449 
450  forAllRow(pointTets_, centreI, ptI)
451  {
452  const partTet& tet = tets_[pointTets_(centreI, ptI)];
453 
454  for(label i=0;i<4;++i)
455  if( smoothVertex_[tet[i]] & SMOOTH )
456  neiSmoothPoints.appendIfNotIn(tet[i]);
457  }
458  }
459 
460  //- select the point and mark neighbouring SMOOTH points
461  selectedPoints.append(nodeI);
462  order[nodeI] = internalPointsOrder.size();
463 
464  forAll(neiSmoothPoints, i)
465  helper[neiSmoothPoints[i]] = true;
466  }
467  }
468 
469  if( selectedPoints.size() != 0 )
470  {
471  internalPointsOrder.appendList(selectedPoints);
472  found = true;
473  }
474 
475  } while( found );
476 }
477 
478 void partTetMesh::createBOUNDARYPointsOrdering() const
479 {
480  boundaryPointsOrderPtr_ = new VRWGraph();
481  VRWGraph& boundaryPointsOrder = *boundaryPointsOrderPtr_;
482 
483  boundaryPointsOrder.setSize(0);
484  labelLongList order(points_.size(), -1);
485  boolList helper(points_.size());
486 
487  bool found;
488  do
489  {
490  found = false;
491  helper = false;
492 
493  labelLongList selectedPoints;
494  forAll(points_, nodeI)
495  {
496  if( smoothVertex_[nodeI] & BOUNDARY )
497  {
498  if( helper[nodeI] )
499  continue;
500  if( order[nodeI] != -1 )
501  continue;
502 
503  //- find neighbouring FACECENTRE and CELLCENTRE points
504  DynList<label, 64> neiCentrePoints, neiSmoothPoints;
505  forAllRow(pointTets_, nodeI, ptI)
506  {
507  const partTet& tet = tets_[pointTets_(nodeI, ptI)];
508 
509  for(label i=0;i<4;++i)
510  if( smoothVertex_[tet[i]] & (FACECENTRE+CELLCENTRE) )
511  {
512  neiCentrePoints.appendIfNotIn(tet[i]);
513  }
514  else if( smoothVertex_[tet[i]] & BOUNDARY )
515  {
516  neiSmoothPoints.appendIfNotIn(tet[i]);
517  }
518  }
519 
520  //- find neighbouring BOUNDARY points
521  forAll(neiCentrePoints, ncI)
522  {
523  const label centreI = neiCentrePoints[ncI];
524 
525  forAllRow(pointTets_, centreI, ptI)
526  {
527  const partTet& tet = tets_[pointTets_(centreI, ptI)];
528 
529  for(label i=0;i<4;++i)
530  if( smoothVertex_[tet[i]] & BOUNDARY )
531  neiSmoothPoints.appendIfNotIn(tet[i]);
532  }
533  }
534 
535  //- select the point and mark neighbouring BOUNDARY points
536  selectedPoints.append(nodeI);
537  order[nodeI] = boundaryPointsOrder.size();
538 
539  forAll(neiSmoothPoints, i)
540  helper[neiSmoothPoints[i]] = true;
541  }
542  }
543 
544  if( selectedPoints.size() != 0 )
545  {
546  boundaryPointsOrder.appendList(selectedPoints);
547  found = true;
548  }
549 
550  } while( found );
551 }
552 
553 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
554 
555 } // End namespace Foam
556 
557 // ************************************************************************* //
Foam::tetMatcher::matchShape
virtual bool matchShape(const bool checkOnly, const faceList &faces, const labelList &faceOwner, const label cellI, const labelList &myFaces)
Low level shape recognition. Return true if matches.
Definition: tetMatcher.C:64
Foam::LongList::append
void append(const T &e)
Append an element at the end of the list.
Definition: LongListI.H:265
Foam::partTetMesh::FACECENTRE
@ FACECENTRE
Definition: partTetMesh.H:162
Foam::partTetMesh::smoothVertex_
LongList< direction > smoothVertex_
shall a node be used for smoothing or not
Definition: partTetMesh.H:75
Foam::boundaryPatchBase::patchSize
label patchSize() const
Definition: boundaryPatchBase.H:161
helperFunctionsPar.H
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:406
Foam::partTetMesh::boundaryPointsOrderPtr_
VRWGraph * boundaryPointsOrderPtr_
boundary point ordering for parallel runs
Definition: partTetMesh.H:84
Foam::OPstream
Output inter-processor communications stream.
Definition: OPstream.H:50
demandDrivenData.H
Template functions to aid in the implementation of demand driven data.
Foam::UPstream::parRun
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:377
polyMeshGenModifier.H
Foam::tetMatcher
A cellMatcher for tet cells.
Definition: tetMatcher.H:51
Foam::partTet
Definition: partTet.H:53
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Foam::partTetMesh::tets_
LongList< partTet > tets_
tetrahedra making the mesh
Definition: partTetMesh.H:69
Foam::LongList::size
label size() const
Size of the active part of the list.
Definition: LongListI.H:203
Foam::mag
dimensioned< scalar > mag(const dimensioned< Type > &)
Foam::cellListPMG
Definition: cellListPMG.H:49
Foam::cellMatcher::vertLabels
const labelList & vertLabels() const
Definition: cellMatcherI.H:74
Foam::boundaryPatch
Like polyPatch but without reference to mesh. patchIdentifier::index is not used. Used in boundaryMes...
Definition: boundaryPatch.H:50
Foam::LongList< label >
Foam::UPstream::blocking
@ blocking
Definition: UPstream.H:66
Foam::partTetMesh::internalPointsOrderPtr_
VRWGraph * internalPointsOrderPtr_
internal point ordering for parallel runs
Definition: partTetMesh.H:81
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::Field
Pre-declare SubField and related Field type.
Definition: Field.H:57
Foam::List::append
void append(const T &)
Append an element at the end of the list.
Foam::Info
messageStream Info
Foam::VRWGraph::size
label size() const
Returns the number of rows.
Definition: VRWGraphI.H:122
Foam::DynList::appendIfNotIn
void appendIfNotIn(const T &e)
Definition: DynListI.H:324
Foam::partTetMesh::createSMOOTHPointsOrdering
void createSMOOTHPointsOrdering() const
create ordering of internal points for parallel execution
Definition: partTetMeshAddressing.C:403
partTetMesh.H
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::partTetMesh::createPointsAndTets
void createPointsAndTets(const List< direction > &usedCells, const boolList &lockedPoints)
create points and tets
Definition: partTetMeshAddressing.C:45
Foam::FatalError
error FatalError
Foam::partTetMesh::SMOOTH
@ SMOOTH
Definition: partTetMesh.H:161
Foam
Namespace for OpenFOAM.
Definition: combustionModel.C:30
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:131
Foam::processorBoundaryPatch
Definition: processorBoundaryPatch.H:47
Foam::DynList< label, 64 >
found
bool found
Definition: TABSMDCalcMethod2.H:32
Foam::partTetMesh::BOUNDARY
@ BOUNDARY
Definition: partTetMesh.H:165
Foam::partTetMesh::createBOUNDARYPointsOrdering
void createBOUNDARYPointsOrdering() const
create order of boundary points for parallel execution
Definition: partTetMeshAddressing.C:478
Foam::partTetMesh::CELLCENTRE
@ CELLCENTRE
Definition: partTetMesh.H:163
Foam::partTetMesh::pointTets_
VRWGraph pointTets_
addressing data
Definition: partTetMesh.H:78
Foam::LongList::byteSize
label byteSize() const
Return the binary size in number of characters of the UList.
Definition: LongListI.H:209
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::boundaryPatchBase::patchStart
label patchStart() const
Definition: boundaryPatchBase.H:151
Foam::List< direction >
points
const pointField & points
Definition: gmvOutputHeader.H:1
tetMatcher.H
Foam::direction
unsigned char direction
Definition: direction.H:43
Foam::constant::universal::c
const dimensionedScalar c
Speed of light in a vacuum.
Foam::partTetMesh::points_
LongList< point > points_
points in the tet mesh
Definition: partTetMesh.H:66
Foam::IPstream
Input inter-processor communications stream.
Definition: IPstream.H:50
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
cells
const cellShapeList & cells
Definition: gmvOutputHeader.H:3
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

Copyright © 2011-2018  OpenFOAM | OPENFOAM® is a registered trademark of OpenCFD Ltd.

Generated by   doxygen 1.8.17