conformalVoronoiMesh.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) 2012-2022 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 "conformalVoronoiMesh.H"
27 #include "initialPointsMethod.H"
28 #include "relaxationModel.H"
29 #include "faceAreaWeightModel.H"
30 #include "meshSearch.H"
31 #include "vectorTools.H"
32 #include "IOmanip.H"
33 #include "indexedCellChecks.H"
34 #include "controlMeshRefinement.H"
35 #include "smoothAlignmentSolver.H"
36 #include "OBJstream.H"
37 #include "indexedVertexOps.H"
38 #include "DelaunayMeshTools.H"
39 
40 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
41 
42 namespace Foam
43 {
44  defineTypeNameAndDebug(conformalVoronoiMesh, 0);
45 
46  template<>
47  const char* NamedEnum
48  <
50  5
51  >::names[] =
52  {
53  "internal",
54  "surface",
55  "featureEdge",
56  "featurePoint",
57  "constrained"
58  };
59 }
60 
63 
64 
65 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
66 
67 void Foam::conformalVoronoiMesh::cellSizeMeshOverlapsBackground() const
68 {
69  const cellShapeControlMesh& cellSizeMesh =
70  cellShapeControl_.shapeControlMesh();
71 
72  DynamicList<Foam::point> pts(number_of_vertices());
73 
74  for
75  (
76  Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
77  vit != finite_vertices_end();
78  ++vit
79  )
80  {
81  if (vit->internalOrBoundaryPoint() && !vit->referred())
82  {
83  pts.append(topoint(vit->point()));
84  }
85  }
86 
87  boundBox bb(pts);
88 
89  boundBox cellSizeMeshBb = cellSizeMesh.bounds();
90 
91  bool fullyContained = true;
92 
93  if (!cellSizeMeshBb.contains(bb))
94  {
95  Pout<< "Triangulation not fully contained in cell size mesh."
96  << endl;
97 
98  Pout<< "Cell Size Mesh Bounds = " << cellSizeMesh.bounds() << endl;
99  Pout<< "foamyHexMesh Bounds = " << bb << endl;
100 
101  fullyContained = false;
102  }
103 
104  reduce(fullyContained, andOp<unsigned int>());
105 
106  Info<< "Triangulation is "
107  << (fullyContained ? "fully" : "not fully")
108  << " contained in the cell size mesh"
109  << endl;
110 }
111 
112 
113 void Foam::conformalVoronoiMesh::insertInternalPoints
114 (
115  List<Point>& points,
116  bool distribute
117 )
118 {
119  label nPoints = points.size();
120 
121  if (Pstream::parRun())
122  {
123  reduce(nPoints, sumOp<label>());
124  }
125 
126  Info<< " " << nPoints << " points to insert..." << endl;
127 
128  if (Pstream::parRun() && distribute)
129  {
130  List<Foam::point> transferPoints(points.size());
131 
132  forAll(points, pI)
133  {
134  transferPoints[pI] = topoint(points[pI]);
135  }
136 
137  // Send the points that are not on this processor to the appropriate
138  // place
140  (
141  decomposition_().distributePoints(transferPoints)
142  );
143 
144  transferPoints.clear();
145 
146  map().distribute(points);
147  }
148 
149  label nVert = number_of_vertices();
150 
151  insert(points.begin(), points.end());
152 
153  label nInserted(number_of_vertices() - nVert);
154 
155  if (Pstream::parRun())
156  {
157  reduce(nInserted, sumOp<label>());
158  }
159 
160  Info<< " " << nInserted << " points inserted"
161  << ", failed to insert " << nPoints - nInserted
162  << " ("
163  << 100.0*(nPoints - nInserted)/(nInserted + small)
164  << " %)"<< endl;
165 
166  for
167  (
168  Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
169  vit != finite_vertices_end();
170  ++vit
171  )
172  {
174  {
175  vit->index() = getNewVertexIndex();
176  vit->type() = Vb::vtInternal;
177  }
178  }
179 }
180 
181 
182 Foam::Map<Foam::label> Foam::conformalVoronoiMesh::insertPointPairs
183 (
184  List<Vb>& vertices,
185  bool distribute,
186  bool reIndex
187 )
188 {
189  if (Pstream::parRun() && distribute)
190  {
191  autoPtr<distributionMap> mapDist =
192  decomposition_().distributePoints(vertices);
193 
194  // Re-index the point pairs if one or both have been distributed.
195  // If both, remove
196 
197  // If added a point, then need to know its point pair
198  // If one moved, then update procIndex locally
199 
200  forAll(vertices, vI)
201  {
202  vertices[vI].procIndex() = Pstream::myProcNo();
203  }
204  }
205 
206  label preReinsertionSize(number_of_vertices());
207 
208  Map<label> oldToNewIndices =
209  this->DelaunayMesh<Delaunay>::insertPoints(vertices, reIndex);
210 
211  const label nReinserted = returnReduce
212  (
213  label(number_of_vertices()) - preReinsertionSize,
214  sumOp<label>()
215  );
216 
217  Info<< " Reinserted " << nReinserted << " vertices out of "
218  << returnReduce(vertices.size(), sumOp<label>())
219  << endl;
220 
221  return oldToNewIndices;
222 }
223 
224 
225 void Foam::conformalVoronoiMesh::insertSurfacePointPairs
226 (
227  const pointIndexHitAndFeatureList& surfaceHits,
228  const fileName fName,
229  DynamicList<Vb>& pts
230 )
231 {
232  forAll(surfaceHits, i)
233  {
234  vectorField norm(1);
235 
236  const pointIndexHit surfaceHit = surfaceHits[i].first();
237  const label featureIndex = surfaceHits[i].second();
238 
239  allGeometry_[featureIndex].getNormal
240  (
241  List<pointIndexHit>(1, surfaceHit),
242  norm
243  );
244 
245  const vector& normal = norm[0];
246 
247  const Foam::point& surfacePt(surfaceHit.hitPoint());
248 
250  geometryToConformTo_.meshableSide(featureIndex, surfaceHit);
251 
252  if (meshableSide == extendedFeatureEdgeMesh::BOTH)
253  {
254  createBafflePointPair
255  (
256  pointPairDistance(surfacePt),
257  surfacePt,
258  normal,
259  true,
260  pts
261  );
262  }
263  else if (meshableSide == extendedFeatureEdgeMesh::INSIDE)
264  {
265  createPointPair
266  (
267  pointPairDistance(surfacePt),
268  surfacePt,
269  normal,
270  true,
271  pts
272  );
273  }
274  else if (meshableSide == extendedFeatureEdgeMesh::OUTSIDE)
275  {
276  createPointPair
277  (
278  pointPairDistance(surfacePt),
279  surfacePt,
280  -normal,
281  true,
282  pts
283  );
284  }
285  else
286  {
288  << meshableSide << ", bad"
289  << endl;
290  }
291  }
292 
293  if (foamyHexMeshControls().objOutput() && fName != fileName::null)
294  {
295  DelaunayMeshTools::writeOBJ(time().path()/fName, pts);
296  }
297 }
298 
299 
300 void Foam::conformalVoronoiMesh::insertEdgePointGroups
301 (
302  const pointIndexHitAndFeatureList& edgeHits,
303  const fileName fName,
304  DynamicList<Vb>& pts
305 )
306 {
307  forAll(edgeHits, i)
308  {
309  if (edgeHits[i].first().hit())
310  {
311  const extendedFeatureEdgeMesh& feMesh
312  (
313  geometryToConformTo_.features()[edgeHits[i].second()]
314  );
315 
316 // const bool isBaffle =
317 // geometryToConformTo_.isBaffleFeature(edgeHits[i].second());
318 
319  createEdgePointGroup
320  (
321  feMesh,
322  edgeHits[i].first(),
323  pts
324  );
325  }
326  }
327 
328  if (foamyHexMeshControls().objOutput() && fName != fileName::null)
329  {
330  DelaunayMeshTools::writeOBJ(time().path()/fName, pts);
331  }
332 }
333 
334 
335 bool Foam::conformalVoronoiMesh::nearFeaturePt(const Foam::point& pt) const
336 {
337  scalar exclusionRangeSqr = featurePointExclusionDistanceSqr(pt);
338 
339  pointIndexHit info;
340  label featureHit;
341 
342  geometryToConformTo_.findFeaturePointNearest
343  (
344  pt,
345  exclusionRangeSqr,
346  info,
347  featureHit
348  );
349 
350  return info.hit();
351 }
352 
353 
354 bool Foam::conformalVoronoiMesh::surfacePtNearFeatureEdge
355 (
356  const Foam::point& pt
357 ) const
358 {
359  scalar exclusionRangeSqr = surfacePtExclusionDistanceSqr(pt);
360 
361  pointIndexHit info;
362  label featureHit;
363 
364  geometryToConformTo_.findEdgeNearest
365  (
366  pt,
367  exclusionRangeSqr,
368  info,
369  featureHit
370  );
371 
372  return info.hit();
373 }
374 
375 
376 void Foam::conformalVoronoiMesh::insertInitialPoints()
377 {
378  Info<< nl << "Inserting initial points" << endl;
379 
380  timeCheck("Before initial points call");
381 
382  List<Point> initPts = initialPointsMethod_->initialPoints();
383 
384  timeCheck("After initial points call");
385 
386  // Assume that the initial points method made the correct decision for
387  // which processor each point should be on, so give distribute = false
388  insertInternalPoints(initPts);
389 
390  if (initialPointsMethod_->fixInitialPoints())
391  {
392  for
393  (
394  Finite_vertices_iterator vit = finite_vertices_begin();
395  vit != finite_vertices_end();
396  ++vit
397  )
398  {
399  vit->fixed() = true;
400  }
401  }
402 
403  if (foamyHexMeshControls().objOutput())
404  {
406  (
407  time().path()/"initialPoints.obj",
408  *this,
410  );
411  }
412 }
413 
414 
415 void Foam::conformalVoronoiMesh::distribute()
416 {
417  if (!Pstream::parRun())
418  {
419  return;
420  }
421 
422  DynamicList<Foam::point> points(number_of_vertices());
423  DynamicList<Foam::indexedVertexEnum::vertexType> types
424  (
425  number_of_vertices()
426  );
427  DynamicList<scalar> sizes(number_of_vertices());
428  DynamicList<tensor> alignments(number_of_vertices());
429 
430  for
431  (
432  Finite_vertices_iterator vit = finite_vertices_begin();
433  vit != finite_vertices_end();
434  ++vit
435  )
436  {
437  if (vit->real())
438  {
439  points.append(topoint(vit->point()));
440  types.append(vit->type());
441  sizes.append(vit->targetCellSize());
442  alignments.append(vit->alignment());
443  }
444  }
445 
446  autoPtr<distributionMap> mapDist =
447  DistributedDelaunayMesh<Delaunay>::distribute(decomposition_(), points);
448 
449  mapDist().distribute(types);
450  mapDist().distribute(sizes);
451  mapDist().distribute(alignments);
452 
453  // Reset the entire tessellation
455 
456  Info<< nl << " Inserting distributed tessellation" << endl;
457 
458  // Internal points have to be inserted first
459 
460  DynamicList<Vb> verticesToInsert(points.size());
461 
462  forAll(points, pI)
463  {
464  verticesToInsert.append
465  (
466  Vb
467  (
468  toPoint(points[pI]),
469  -1,
470  types[pI],
472  )
473  );
474 
475  verticesToInsert.last().targetCellSize() = sizes[pI];
476  verticesToInsert.last().alignment() = alignments[pI];
477  }
478 
479  this->rangeInsertWithInfo
480  (
481  verticesToInsert.begin(),
482  verticesToInsert.end(),
483  true
484  );
485 
486  Info<< " Total number of vertices after redistribution "
487  << returnReduce
488  (
489  label(number_of_vertices()), sumOp<label>()
490  )
491  << endl;
492 }
493 
494 
495 void Foam::conformalVoronoiMesh::buildCellSizeAndAlignmentMesh()
496 {
497  controlMeshRefinement meshRefinement
498  (
499  cellShapeControl_
500  );
501 
502  smoothAlignmentSolver meshAlignmentSmoother
503  (
504  cellShapeControl_.shapeControlMesh()
505  );
506 
507  meshRefinement.initialMeshPopulation(decomposition_);
508 
509  cellShapeControlMesh& cellSizeMesh = cellShapeControl_.shapeControlMesh();
510 
511  if (Pstream::parRun())
512  {
513  if (!distributeBackground(cellSizeMesh))
514  {
515  // Synchronise the cell size mesh if it has not been distributed
516  cellSizeMesh.distribute(decomposition_);
517  }
518  }
519 
520  const dictionary& motionControlDict
521  = foamyHexMeshControls().foamyHexMeshDict().subDict("motionControl");
522 
523  const label nMaxIter =
524  motionControlDict.lookup<label>("maxRefinementIterations");
525 
526  Info<< "Maximum number of refinement iterations : " << nMaxIter << endl;
527 
528  for (label i = 0; i < nMaxIter; ++i)
529  {
530  label nAdded = meshRefinement.refineMesh(decomposition_);
531  // cellShapeControl_.refineMesh(decomposition_);
532  reduce(nAdded, sumOp<label>());
533 
534  if (Pstream::parRun())
535  {
536  cellSizeMesh.distribute(decomposition_);
537  }
538 
539  Info<< " Iteration " << i
540  << " Added = " << nAdded << " points"
541  << endl;
542 
543  if (nAdded == 0)
544  {
545  break;
546  }
547  }
548 
549  if (Pstream::parRun())
550  {
551  // Need to distribute the cell size mesh to cover the background mesh
552  if (!distributeBackground(cellSizeMesh))
553  {
554  cellSizeMesh.distribute(decomposition_);
555  }
556  }
557 
558  const label maxSmoothingIterations =
559  motionControlDict.lookup<label>("maxSmoothingIterations");
560 
561  meshAlignmentSmoother.smoothAlignments(maxSmoothingIterations);
562 
563  Info<< "Background cell size and alignment mesh:" << endl;
564  cellSizeMesh.printInfo(Info);
565 
566  Info<< "Triangulation is "
567  << (cellSizeMesh.is_valid() ? "valid" : "not valid!" )
568  << endl;
569 
570  if (foamyHexMeshControls().writeCellShapeControlMesh())
571  {
572  // cellSizeMesh.writeTriangulation();
573  cellSizeMesh.write();
574  }
575 
576  if (foamyHexMeshControls().printVertexInfo())
577  {
578  cellSizeMesh.printVertexInfo(Info);
579  }
580 
581 // Info<< "Estimated number of cells in final mesh = "
582 // << returnReduce
583 // (
584 // cellSizeMesh.estimateCellCount(decomposition_),
585 // sumOp<label>()
586 // )
587 // << endl;
588 }
589 
590 
591 void Foam::conformalVoronoiMesh::setVertexSizeAndAlignment()
592 {
593  Info<< nl << "Calculating target cell alignment and size" << endl;
594 
595  for
596  (
597  Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
598  vit != finite_vertices_end();
599  vit++
600  )
601  {
602  if (vit->internalOrBoundaryPoint())
603  {
604  pointFromPoint pt = topoint(vit->point());
605 
606  cellShapeControls().cellSizeAndAlignment
607  (
608  pt,
609  vit->targetCellSize(),
610  vit->alignment()
611  );
612  }
613  }
614 }
615 
616 
617 Foam::face Foam::conformalVoronoiMesh::buildDualFace
618 (
619  const Delaunay::Finite_edges_iterator& eit
620 ) const
621 {
622  Cell_circulator ccStart = incident_cells(*eit);
623  Cell_circulator cc1 = ccStart;
624  Cell_circulator cc2 = cc1;
625 
626  // Advance the second circulator so that it always stays on the next
627  // cell around the edge;
628  cc2++;
629 
630  DynamicList<label> verticesOnFace;
631 
632  label nUniqueVertices = 0;
633 
634  do
635  {
636  if
637  (
638  cc1->hasFarPoint() || cc2->hasFarPoint()
639  || is_infinite(cc1) || is_infinite(cc2)
640  )
641  {
642  Cell_handle c = eit->first;
643  Vertex_handle vA = c->vertex(eit->second);
644  Vertex_handle vB = c->vertex(eit->third);
645 
646 // DelaunayMeshTools::drawDelaunayCell(Pout, cc1);
647 // DelaunayMeshTools::drawDelaunayCell(Pout, cc2);
648 
650  << "Dual face uses circumcenter defined by a "
651  << "Delaunay tetrahedron with no internal "
652  << "or boundary points. Defining Delaunay edge ends: "
653  << vA->info() << " "
654  << vB->info() << nl
655  <<endl;//<< exit(FatalError);
656  }
657  else
658  {
659  label cc1I = cc1->cellIndex();
660  label cc2I = cc2->cellIndex();
661 
662  if (cc1I != cc2I)
663  {
664  if (findIndex(verticesOnFace, cc1I) == -1)
665  {
666  nUniqueVertices++;
667  }
668 
669  verticesOnFace.append(cc1I);
670  }
671  }
672 
673  cc1++;
674  cc2++;
675 
676  } while (cc1 != ccStart);
677 
678  verticesOnFace.shrink();
679 
680  if (verticesOnFace.size() >= 3 && nUniqueVertices < 3)
681  {
682  // There are not enough unique vertices on this face to
683  // justify its size, it may have a form like:
684 
685  // Vertices:
686  // A B
687  // A B
688 
689  // Face:
690  // ABAB
691 
692  // Setting the size to be below 3, so that it will not be
693  // created
694 
695  verticesOnFace.setSize(nUniqueVertices);
696  }
697 
698  return face(verticesOnFace);
699 }
700 
701 
702 Foam::label Foam::conformalVoronoiMesh::maxFilterCount
703 (
704  const Delaunay::Finite_edges_iterator& eit
705 ) const
706 {
707  Cell_circulator ccStart = incident_cells(*eit);
708  Cell_circulator cc = ccStart;
709 
710  label maxFC = 0;
711 
712  do
713  {
714  if (cc->hasFarPoint())
715  {
716  Cell_handle c = eit->first;
717  Vertex_handle vA = c->vertex(eit->second);
718  Vertex_handle vB = c->vertex(eit->third);
719 
721  << "Dual face uses circumcenter defined by a "
722  << "Delaunay tetrahedron with no internal "
723  << "or boundary points. Defining Delaunay edge ends: "
724  << topoint(vA->point()) << " "
725  << topoint(vB->point()) << nl
726  << exit(FatalError);
727  }
728 
729  if (cc->filterCount() > maxFC)
730  {
731  maxFC = cc->filterCount();
732  }
733 
734  cc++;
735 
736  } while (cc != ccStart);
737 
738  return maxFC;
739 }
740 
741 
742 bool Foam::conformalVoronoiMesh::ownerAndNeighbour
743 (
744  Vertex_handle vA,
745  Vertex_handle vB,
746  label& owner,
747  label& neighbour
748 ) const
749 {
750  bool reverse = false;
751 
752  owner = -1;
753  neighbour = -1;
754 
755  label dualCellIndexA = vA->index();
756 
757  if (!vA->internalOrBoundaryPoint() || vA->referred())
758  {
759  if (!vA->constrained())
760  {
761  dualCellIndexA = -1;
762  }
763  }
764 
765  label dualCellIndexB = vB->index();
766 
767  if (!vB->internalOrBoundaryPoint() || vB->referred())
768  {
769  if (!vB->constrained())
770  {
771  dualCellIndexB = -1;
772  }
773  }
774 
775  if (dualCellIndexA == -1 && dualCellIndexB == -1)
776  {
778  << "Attempting to create a face joining "
779  << "two unindexed dual cells "
780  << exit(FatalError);
781  }
782  else if (dualCellIndexA == -1 || dualCellIndexB == -1)
783  {
784  // boundary face, find which is the owner
785 
786  if (dualCellIndexA == -1)
787  {
788  owner = dualCellIndexB;
789 
790  reverse = true;
791  }
792  else
793  {
794  owner = dualCellIndexA;
795  }
796  }
797  else
798  {
799  // internal face, find the lower cell to be the owner
800 
801  if (dualCellIndexB > dualCellIndexA)
802  {
803  owner = dualCellIndexA;
804  neighbour = dualCellIndexB;
805  }
806  else
807  {
808  owner = dualCellIndexB;
809  neighbour = dualCellIndexA;
810 
811  // reverse face order to correctly orientate normal
812  reverse = true;
813  }
814  }
815 
816  return reverse;
817 }
818 
819 
820 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
821 
823 (
824  const Time& runTime,
825  const dictionary& foamyHexMeshDict
826 )
827 :
828  DistributedDelaunayMesh<Delaunay>(runTime),
829  runTime_(runTime),
830  rndGen_(64293*Pstream::myProcNo()),
831  foamyHexMeshControls_(foamyHexMeshDict),
832  allGeometry_
833  (
834  IOobject
835  (
836  "cvSearchableSurfaces",
837  runTime_.constant(),
838  searchableSurface::geometryDir(runTime),
839  runTime_,
840  IOobject::MUST_READ,
841  IOobject::NO_WRITE
842  ),
843  foamyHexMeshDict.subDict("geometry"),
844  foamyHexMeshDict.lookupOrDefault("singleRegionName", true)
845  ),
846  geometryToConformTo_
847  (
848  runTime_,
849  rndGen_,
850  allGeometry_,
851  foamyHexMeshDict.subDict("surfaceConformation")
852  ),
853  decomposition_
854  (
855  Pstream::parRun()
856  ? new backgroundMeshDecomposition
857  (
858  runTime_,
859  rndGen_,
860  geometryToConformTo_,
861  foamyHexMeshControls().foamyHexMeshDict().subDict
862  (
863  "backgroundMeshDecomposition"
864  )
865  )
866  : nullptr
867  ),
868  cellShapeControl_
869  (
870  runTime_,
871  foamyHexMeshControls_,
872  allGeometry_,
873  geometryToConformTo_
874  ),
875  limitBounds_(),
876  ptPairs_(*this),
877  ftPtConformer_(*this),
878  edgeLocationTreePtr_(),
879  surfacePtLocationTreePtr_(),
880  surfaceConformationVertices_(),
881  initialPointsMethod_
882  (
883  initialPointsMethod::New
884  (
885  foamyHexMeshDict.subDict("initialPoints"),
886  runTime_,
887  rndGen_,
888  geometryToConformTo_,
889  cellShapeControl_,
890  decomposition_
891  )
892  ),
893  relaxationModel_
894  (
895  relaxationModel::New
896  (
897  foamyHexMeshDict.subDict("motionControl"),
898  runTime_
899  )
900  ),
901  faceAreaWeightModel_
902  (
903  faceAreaWeightModel::New
904  (
905  foamyHexMeshDict.subDict("motionControl")
906  )
907  )
908 {}
909 
910 
911 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
912 
914 {}
915 
916 
917 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
918 
920 {
921  if (foamyHexMeshControls().objOutput())
922  {
923  geometryToConformTo_.writeFeatureObj("foamyHexMesh");
924  }
925 
926  buildCellSizeAndAlignmentMesh();
927 
928  insertInitialPoints();
929 
930  insertFeaturePoints(true);
931 
932  setVertexSizeAndAlignment();
933 
934  cellSizeMeshOverlapsBackground();
935 
936  // Improve the guess that the backgroundMeshDecomposition makes with the
937  // initial positions. Use before building the surface conformation to
938  // better balance the surface conformation load.
939  distributeBackground(*this);
940 
941  buildSurfaceConformation();
942 
943  // The introduction of the surface conformation may have distorted the
944  // balance of vertices, distribute if necessary.
945  distributeBackground(*this);
946 
947  if (Pstream::parRun())
948  {
949  sync(decomposition_().procBounds());
950  }
951 
952  // Do not store the surface conformation until after it has been
953  // (potentially) redistributed.
954  storeSurfaceConformation();
955 
956  // Report any Delaunay vertices that do not think that they are in the
957  // domain the processor they are on.
958  // reportProcessorOccupancy();
959 
960  cellSizeMeshOverlapsBackground();
961 
963  {
965  }
966 
967  if (foamyHexMeshControls().objOutput())
968  {
970  (
971  time().path()/"internalPoints_" + time().timeName() + ".obj",
972  *this,
975  );
976  }
977 }
978 
979 
981 {
982  if (Pstream::parRun())
983  {
984  decomposition_.reset
985  (
986  new backgroundMeshDecomposition
987  (
988  runTime_,
989  rndGen_,
990  geometryToConformTo_,
991  foamyHexMeshControls().foamyHexMeshDict().subDict
992  (
993  "backgroundMeshDecomposition"
994  )
995  )
996  );
997  }
998 
999  insertInitialPoints();
1000 
1001  insertFeaturePoints();
1002 
1003  // Improve the guess that the backgroundMeshDecomposition makes with the
1004  // initial positions. Use before building the surface conformation to
1005  // better balance the surface conformation load.
1006  distributeBackground(*this);
1007 
1008  buildSurfaceConformation();
1009 
1010  // The introduction of the surface conformation may have distorted the
1011  // balance of vertices, distribute if necessary.
1012  distributeBackground(*this);
1013 
1014  if (Pstream::parRun())
1015  {
1016  sync(decomposition_().procBounds());
1017  }
1018 
1019  cellSizeMeshOverlapsBackground();
1020 
1022  {
1024  }
1025 }
1026 
1027 
1029 {
1030  timeCheck("Start of move");
1031 
1032  scalar relaxation = relaxationModel_->relaxation();
1033 
1034  Info<< nl << "Relaxation = " << relaxation << endl;
1035 
1036  pointField dualVertices(number_of_finite_cells());
1037 
1038  this->resetCellCount();
1039 
1040  // Find the dual point of each tetrahedron and assign it an index.
1041  for
1042  (
1043  Delaunay::Finite_cells_iterator cit = finite_cells_begin();
1044  cit != finite_cells_end();
1045  ++cit
1046  )
1047  {
1048  cit->cellIndex() = Cb::ctUnassigned;
1049 
1050  if (cit->anyInternalOrBoundaryDualVertex())
1051  {
1052  cit->cellIndex() = getNewCellIndex();
1053 
1054  dualVertices[cit->cellIndex()] = cit->dual();
1055  }
1056 
1057  if (cit->hasFarPoint())
1058  {
1059  cit->cellIndex() = Cb::ctFar;
1060  }
1061  }
1062 
1063  dualVertices.setSize(cellCount());
1064 
1065  setVertexSizeAndAlignment();
1066 
1067  timeCheck("Determined sizes and alignments");
1068 
1069  Info<< nl << "Determining vertex displacements" << endl;
1070 
1071  vectorField cartesianDirections(3);
1072 
1073  cartesianDirections[0] = vector(1, 0, 0);
1074  cartesianDirections[1] = vector(0, 1, 0);
1075  cartesianDirections[2] = vector(0, 0, 1);
1076 
1077  vectorField displacementAccumulator
1078  (
1079  number_of_vertices(),
1080  Zero
1081  );
1082 
1083  PackedBoolList pointToBeRetained
1084  (
1085  number_of_vertices(),
1086  true
1087  );
1088 
1089  DynamicList<Point> pointsToInsert(number_of_vertices());
1090 
1091  for
1092  (
1093  Delaunay::Finite_edges_iterator eit = finite_edges_begin();
1094  eit != finite_edges_end();
1095  ++eit
1096  )
1097  {
1098  Cell_handle c = eit->first;
1099  Vertex_handle vA = c->vertex(eit->second);
1100  Vertex_handle vB = c->vertex(eit->third);
1101 
1102  if
1103  (
1104  (
1105  vA->internalPoint() && !vA->referred()
1106  && vB->internalOrBoundaryPoint()
1107  )
1108  || (
1109  vB->internalPoint() && !vB->referred()
1110  && vA->internalOrBoundaryPoint()
1111  )
1112  )
1113  {
1114  pointFromPoint dVA = topoint(vA->point());
1115  pointFromPoint dVB = topoint(vB->point());
1116 
1117  Field<vector> alignmentDirsA
1118  (
1119  vA->alignment().T() & cartesianDirections
1120  );
1121  Field<vector> alignmentDirsB
1122  (
1123  vB->alignment().T() & cartesianDirections
1124  );
1125 
1126  Field<vector> alignmentDirs(alignmentDirsA);
1127 
1128  forAll(alignmentDirsA, aA)
1129  {
1130  const vector& a = alignmentDirsA[aA];
1131 
1132  scalar maxDotProduct = 0.0;
1133 
1134  forAll(alignmentDirsB, aB)
1135  {
1136  const vector& b = alignmentDirsB[aB];
1137 
1138  const scalar dotProduct = a & b;
1139 
1140  if (mag(dotProduct) > maxDotProduct)
1141  {
1142  maxDotProduct = mag(dotProduct);
1143 
1144  alignmentDirs[aA] = a + sign(dotProduct)*b;
1145 
1146  alignmentDirs[aA] /= mag(alignmentDirs[aA]);
1147  }
1148  }
1149  }
1150 
1151  vector rAB = dVA - dVB;
1152 
1153  scalar rABMag = mag(rAB);
1154 
1155  if (rABMag < small)
1156  {
1157  // Removal of close points
1158 
1159  if
1160  (
1161  vA->internalPoint() && !vA->referred() && !vA->fixed()
1162  && vB->internalPoint() && !vB->referred() && !vB->fixed()
1163  )
1164  {
1165  // Only insert a point at the midpoint of
1166  // the short edge if neither attached
1167  // point has already been identified to be
1168  // removed.
1169 
1170  if
1171  (
1172  pointToBeRetained[vA->index()] == true
1173  && pointToBeRetained[vB->index()] == true
1174  )
1175  {
1176  const Foam::point pt(0.5*(dVA + dVB));
1177 
1178  if (internalPointIsInside(pt))
1179  {
1180  pointsToInsert.append(toPoint(pt));
1181  }
1182  }
1183  }
1184 
1185  if (vA->internalPoint() && !vA->referred() && !vA->fixed())
1186  {
1187  pointToBeRetained[vA->index()] = false;
1188  }
1189 
1190  if (vB->internalPoint() && !vB->referred() && !vB->fixed())
1191  {
1192  pointToBeRetained[vB->index()] = false;
1193  }
1194 
1195  // Do not consider this Delaunay edge any further
1196 
1197  continue;
1198  }
1199 
1200  forAll(alignmentDirs, aD)
1201  {
1202  vector& alignmentDir = alignmentDirs[aD];
1203 
1204  scalar dotProd = rAB & alignmentDir;
1205 
1206  if (dotProd < 0)
1207  {
1208  // swap the direction of the alignment so that has the
1209  // same sense as rAB
1210  alignmentDir *= -1;
1211  dotProd *= -1;
1212  }
1213 
1214  const scalar alignmentDotProd = dotProd/rABMag;
1215 
1216  if
1217  (
1218  alignmentDotProd
1219  > foamyHexMeshControls().cosAlignmentAcceptanceAngle()
1220  )
1221  {
1222  scalar targetCellSize =
1224 
1225  scalar targetFaceArea = sqr(targetCellSize);
1226 
1227  const vector originalAlignmentDir = alignmentDir;
1228 
1229  // Update cell size and face area
1231  (
1232  alignmentDir,
1233  targetFaceArea,
1234  targetCellSize
1235  );
1236 
1237  // Vector to move end points around middle of vector
1238  // to align edge (i.e. dual face normal) with alignment
1239  // directions.
1240  vector delta = alignmentDir - 0.5*rAB;
1241 
1242  face dualFace = buildDualFace(eit);
1243 
1244 // Pout<< dualFace << endl;
1245 // Pout<< " " << vA->info() << endl;
1246 // Pout<< " " << vB->info() << endl;
1247 
1248  const scalar faceArea = dualFace.mag(dualVertices);
1249 
1250  // Update delta vector
1252  (
1253  originalAlignmentDir,
1254  targetCellSize,
1255  rABMag,
1256  delta
1257  );
1258 
1259 // if (targetFaceArea == 0)
1260 // {
1261 // Pout<< vA->info() << vB->info();
1262 //
1263 // Cell_handle ch = locate(vA->point());
1264 // if (is_infinite(ch))
1265 // {
1266 // Pout<< "vA " << vA->targetCellSize() << endl;
1267 // }
1268 //
1269 // ch = locate(vB->point());
1270 // if (is_infinite(ch))
1271 // {
1272 // Pout<< "vB " << vB->targetCellSize() << endl;
1273 // }
1274 // }
1275 
1276  delta *= faceAreaWeightModel_->faceAreaWeight
1277  (
1278  faceArea/targetFaceArea
1279  );
1280 
1281  if
1282  (
1283  (
1284  (vA->internalPoint() && vB->internalPoint())
1285  && (!vA->referred() || !vB->referred())
1286 // ||
1287 // (
1288 // vA->referredInternalPoint()
1289 // && vB->referredInternalPoint()
1290 // )
1291  )
1292  && rABMag
1294  *targetCellSize
1295  && faceArea
1297  *targetFaceArea
1298  && alignmentDotProd
1300  )
1301  {
1302  // Point insertion
1303  if
1304  (
1305  !geometryToConformTo_.findSurfaceAnyIntersection
1306  (
1307  dVA,
1308  dVB
1309  )
1310  )
1311  {
1312  const Foam::point newPt(0.5*(dVA + dVB));
1313 
1314  // Prevent insertions spanning surfaces
1315  if (internalPointIsInside(newPt))
1316  {
1317  if (Pstream::parRun())
1318  {
1319  if
1320  (
1321  decomposition().positionOnThisProcessor
1322  (
1323  newPt
1324  )
1325  )
1326  {
1327  pointsToInsert.append(toPoint(newPt));
1328  }
1329  }
1330  else
1331  {
1332  pointsToInsert.append(toPoint(newPt));
1333  }
1334  }
1335  }
1336  }
1337  else if
1338  (
1339  (
1340  (vA->internalPoint() && !vA->referred())
1341  || (vB->internalPoint() && !vB->referred())
1342  )
1343  && rABMag
1345  *targetCellSize
1346  )
1347  {
1348  // Point removal
1349  if
1350  (
1351  (
1352  vA->internalPoint()
1353  && !vA->referred()
1354  && !vA->fixed()
1355  )
1356  &&
1357  (
1358  vB->internalPoint()
1359  && !vB->referred()
1360  && !vB->fixed()
1361  )
1362  )
1363  {
1364  // Only insert a point at the midpoint of
1365  // the short edge if neither attached
1366  // point has already been identified to be
1367  // removed.
1368  if
1369  (
1370  pointToBeRetained[vA->index()] == true
1371  && pointToBeRetained[vB->index()] == true
1372  )
1373  {
1374  const Foam::point pt(0.5*(dVA + dVB));
1375 
1376  if (internalPointIsInside(pt))
1377  {
1378  pointsToInsert.append(toPoint(pt));
1379  }
1380  }
1381  }
1382 
1383  if
1384  (
1385  vA->internalPoint()
1386  && !vA->referred()
1387  && !vA->fixed()
1388  )
1389  {
1390  pointToBeRetained[vA->index()] = false;
1391  }
1392 
1393  if
1394  (
1395  vB->internalPoint()
1396  && !vB->referred()
1397  && !vB->fixed()
1398  )
1399  {
1400  pointToBeRetained[vB->index()] = false;
1401  }
1402  }
1403  else
1404  {
1405  if
1406  (
1407  vA->internalPoint()
1408  && !vA->referred()
1409  && !vA->fixed()
1410  )
1411  {
1412  if (vB->fixed())
1413  {
1414  displacementAccumulator[vA->index()] += 2*delta;
1415  }
1416  else
1417  {
1418  displacementAccumulator[vA->index()] += delta;
1419  }
1420  }
1421 
1422  if
1423  (
1424  vB->internalPoint()
1425  && !vB->referred()
1426  && !vB->fixed()
1427  )
1428  {
1429  if (vA->fixed())
1430  {
1431  displacementAccumulator[vB->index()] -= 2*delta;
1432  }
1433  else
1434  {
1435  displacementAccumulator[vB->index()] -= delta;
1436  }
1437  }
1438  }
1439  }
1440  }
1441  }
1442  }
1443 
1444  Info<< "Limit displacements" << endl;
1445 
1446  // Limit displacements that pierce, or get too close to the surface
1447  for
1448  (
1449  Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
1450  vit != finite_vertices_end();
1451  ++vit
1452  )
1453  {
1454  if (vit->internalPoint() && !vit->referred() && !vit->fixed())
1455  {
1456  if (pointToBeRetained[vit->index()] == true)
1457  {
1458  limitDisplacement
1459  (
1460  vit,
1461  displacementAccumulator[vit->index()]
1462  );
1463  }
1464  }
1465  }
1466 
1467  vector totalDisp = gSum(displacementAccumulator);
1468  scalar totalDist = gSum(mag(displacementAccumulator));
1469 
1470  displacementAccumulator *= relaxation;
1471 
1472  Info<< "Sum displacements" << endl;
1473 
1474  label nPointsToRetain = 0;
1475  label nPointsToRemove = 0;
1476 
1477  for
1478  (
1479  Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
1480  vit != finite_vertices_end();
1481  ++vit
1482  )
1483  {
1484  if (vit->internalPoint() && !vit->referred() && !vit->fixed())
1485  {
1486  if (pointToBeRetained[vit->index()] == true)
1487  {
1488  // Convert vit->point() to FOAM vector (double) to do addition,
1489  // avoids memory increase because a record of the constructions
1490  // would be kept otherwise.
1491  // See cgal-discuss@lists-sop.inria.fr:
1492  // "Memory issue with openSUSE 11.3, exact kernel, adding
1493  // points/vectors"
1494  // 14/1/2011.
1495  // Only necessary if using an exact constructions kernel
1496  // (extended precision)
1497  Foam::point pt
1498  (
1499  topoint(vit->point())
1500  + displacementAccumulator[vit->index()]
1501  );
1502 
1503  if (internalPointIsInside(pt))
1504  {
1505  pointsToInsert.append(toPoint(pt));
1506  nPointsToRemove++;
1507  }
1508 
1509  nPointsToRetain++;
1510  }
1511  }
1512  }
1513 
1514  pointsToInsert.shrink();
1515 
1516  Info<< returnReduce
1517  (
1518  nPointsToRetain - nPointsToRemove,
1519  sumOp<label>()
1520  )
1521  << " internal points are outside the domain. "
1522  << "They will not be inserted." << endl;
1523 
1524  // Save displacements to file.
1525  if (foamyHexMeshControls().objOutput() && time().writeTime())
1526  {
1527  Info<< "Writing point displacement vectors to file." << endl;
1528  OFstream str
1529  (
1530  time().path()/"displacements_" + runTime_.timeName() + ".obj"
1531  );
1532 
1533  for
1534  (
1535  Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
1536  vit != finite_vertices_end();
1537  ++vit
1538  )
1539  {
1540  if (vit->internalPoint() && !vit->referred())
1541  {
1542  if (pointToBeRetained[vit->index()] == true)
1543  {
1544  meshTools::writeOBJ(str, topoint(vit->point()));
1545 
1546  str << "vn "
1547  << displacementAccumulator[vit->index()][0] << " "
1548  << displacementAccumulator[vit->index()][1] << " "
1549  << displacementAccumulator[vit->index()][2] << " "
1550  << endl;
1551  }
1552  }
1553  }
1554  }
1555 
1556  // Remove the entire tessellation
1558 
1559  timeCheck("Displacement calculated");
1560 
1561  Info<< nl<< "Inserting displaced tessellation" << endl;
1562 
1563  insertFeaturePoints(true);
1564 
1565  insertInternalPoints(pointsToInsert, true);
1566 
1567  // Fix points that have not been significantly displaced
1568 // for
1569 // (
1570 // Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
1571 // vit != finite_vertices_end();
1572 // ++vit
1573 // )
1574 // {
1575 // if (vit->internalPoint())
1576 // {
1577 // if
1578 // (
1579 // mag(displacementAccumulator[vit->index()])
1580 // < 0.1*targetCellSize(topoint(vit->point()))
1581 // )
1582 // {
1583 // vit->setVertexFixed();
1584 // }
1585 // }
1586 // }
1587 
1588  timeCheck("Internal points inserted");
1589 
1590  {
1591  // Check that no index is shared between any of the local points
1592  labelHashSet usedIndices;
1593  for
1594  (
1595  Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
1596  vit != finite_vertices_end();
1597  ++vit
1598  )
1599  {
1600  if (!vit->referred() && !usedIndices.insert(vit->index()))
1601  {
1603  << "Index already used! Could not insert: " << nl
1604  << vit->info()
1605  << abort(FatalError);
1606  }
1607  }
1608  }
1609 
1610  conformToSurface();
1611 
1612  if (foamyHexMeshControls().objOutput())
1613  {
1615  (
1616  time().path()/"internalPoints_" + time().timeName() + ".obj",
1617  *this,
1619  );
1620 
1621  if (reconformToSurface())
1622  {
1624  (
1625  time().path()/"boundaryPoints_" + time().timeName() + ".obj",
1626  *this
1627  );
1628 
1630  (
1631  time().path()/"internalBoundaryPoints_" + time().timeName()
1632  + ".obj",
1633  *this,
1636  );
1637 
1639  (
1640  time().path()/"externalBoundaryPoints_" + time().timeName()
1641  + ".obj",
1642  *this,
1645  );
1646 
1647  OBJstream multipleIntersections
1648  (
1649  "multipleIntersections_"
1650  + time().timeName()
1651  + ".obj"
1652  );
1653 
1654  for
1655  (
1656  Delaunay::Finite_edges_iterator eit = finite_edges_begin();
1657  eit != finite_edges_end();
1658  ++eit
1659  )
1660  {
1661  Cell_handle c = eit->first;
1662  Vertex_handle vA = c->vertex(eit->second);
1663  Vertex_handle vB = c->vertex(eit->third);
1664 
1665  Foam::point ptA(topoint(vA->point()));
1666  Foam::point ptB(topoint(vB->point()));
1667 
1668  List<pointIndexHit> surfHits;
1669  labelList hitSurfaces;
1670 
1672  (
1673  ptA,
1674  ptB,
1675  surfHits,
1676  hitSurfaces
1677  );
1678 
1679  if
1680  (
1681  surfHits.size() >= 2
1682  || (
1683  surfHits.size() == 0
1684  && (
1685  (vA->externalBoundaryPoint()
1686  && vB->internalBoundaryPoint())
1687  || (vB->externalBoundaryPoint()
1688  && vA->internalBoundaryPoint())
1689  )
1690  )
1691  )
1692  {
1693  multipleIntersections.write(linePointRef(ptA, ptB));
1694  }
1695  }
1696  }
1697  }
1698 
1699  timeCheck("After conformToSurface");
1700 
1702  {
1704  }
1705 
1706  if (time().writeTime())
1707  {
1708  writeMesh(time().timeName());
1709  }
1710 
1711  Info<< nl
1712  << "Total displacement = " << totalDisp << nl
1713  << "Total distance = " << totalDist << nl
1714  << endl;
1715 }
1716 
1717 
1718 void Foam::conformalVoronoiMesh::checkCoPlanarCells() const
1719 {
1720  typedef CGAL::Exact_predicates_exact_constructions_kernel Kexact;
1721  typedef CGAL::Point_3<Kexact> PointExact;
1722 
1723  if (!is_valid())
1724  {
1725  Pout<< "Triangulation is invalid!" << endl;
1726  }
1727 
1728  OFstream str("badCells.obj");
1729 
1730  label badCells = 0;
1731 
1732  for
1733  (
1734  Delaunay::Finite_cells_iterator cit = finite_cells_begin();
1735  cit != finite_cells_end();
1736  ++cit
1737  )
1738  {
1739  const scalar quality = foamyHexMeshChecks::coplanarTet(cit, 1e-16);
1740 
1741  if (quality == 0)
1742  {
1743  Pout<< "COPLANAR: " << cit->info() << nl
1744  << " quality = " << quality << nl
1745  << " dual = " << topoint(cit->dual()) << endl;
1746 
1747  DelaunayMeshTools::drawDelaunayCell(str, cit, badCells++);
1748 
1749  FixedList<PointExact, 4> cellVerticesExact(PointExact(0,0,0));
1750  forAll(cellVerticesExact, vI)
1751  {
1752  cellVerticesExact[vI] = PointExact
1753  (
1754  cit->vertex(vI)->point().x(),
1755  cit->vertex(vI)->point().y(),
1756  cit->vertex(vI)->point().z()
1757  );
1758  }
1759 
1760  PointExact synchronisedDual = CGAL::circumcenter<Kexact>
1761  (
1762  cellVerticesExact[0],
1763  cellVerticesExact[1],
1764  cellVerticesExact[2],
1765  cellVerticesExact[3]
1766  );
1767 
1768  Foam::point exactPt
1769  (
1770  CGAL::to_double(synchronisedDual.x()),
1771  CGAL::to_double(synchronisedDual.y()),
1772  CGAL::to_double(synchronisedDual.z())
1773  );
1774 
1775  Info<< "inexact = " << cit->dual() << nl
1776  << "exact = " << exactPt << endl;
1777  }
1778  }
1779 
1780  Pout<< "There are " << badCells << " bad cells out of "
1781  << number_of_finite_cells() << endl;
1782 
1783 
1784  label nNonGabriel = 0;
1785  for
1786  (
1787  Delaunay::Finite_facets_iterator fit = finite_facets_begin();
1788  fit != finite_facets_end();
1789  ++fit
1790  )
1791  {
1792  if (!is_Gabriel(*fit))
1793  {
1794  nNonGabriel++;//Pout<< "Non-gabriel face" << endl;
1795  }
1796  }
1797 
1798  Pout<< "There are " << nNonGabriel << " non-Gabriel faces out of "
1799  << number_of_finite_facets() << endl;
1800 }
1801 
1802 
1803 // ************************************************************************* //
scalar delta
dimensionedScalar sign(const dimensionedScalar &ds)
labelList first(const UList< labelPair > &p)
Definition: patchToPatch.C:38
Delaunay::Cell_handle Cell_handle
CGAL::Delaunay_triangulation_3< K, Tds, CompactLocator > Delaunay
label cellCount() const
Return the cell count (the next unique cell index)
Definition: DelaunayMeshI.H:93
~conformalVoronoiMesh()
Destructor.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
FvWallInfoData< WallInfo, label > label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Delaunay::Vertex_handle Vertex_handle
pointFromPoint topoint(const Point &P)
void writeOBJ(const fileName &fName, const List< Foam::point > &points)
Write list of points to file.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const cellShapeControl & cellShapeControls() const
Return the cellShapeControl object.
static const fileName null
An empty fileName.
Definition: fileName.H:97
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:429
const dimensionedScalar b
Wien displacement law constant: default SI units: [m K].
Definition: createFields.H:27
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
bool distribute(const boundBox &bb)
Holds information (coordinate and normal) regarding nearest wall point.
scalar faceAreaRatioCoeff() const
Return the faceAreaRatioCoeff.
Definition: cvControlsI.H:185
static void timeCheck(const Time &runTime, const string &description=string::null, const bool check=true)
Write the elapsedCpuTime and memory usage, with an optional.
PointIndexHit< point > pointIndexHit
Definition: pointIndexHit.H:42
void writeOBJ(Ostream &os, const point &pt)
Write obj representation of point.
Definition: meshTools.C:203
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
const conformationSurfaces & geometryToConformTo() const
Return the conformationSurfaces object.
void drawDelaunayCell(Ostream &os, const CellHandle &c, label offset=0)
Draws a tet cell to an output stream. The offset is supplied as the tet.
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:51
scalar removalDistCoeff() const
Return removalDistCoeff.
Definition: cvControlsI.H:197
void writeMesh(const fileName &instance)
Prepare data and call writeMesh for polyMesh and.
void printVertexInfo(Ostream &os) const
Write vertex statistics in the form of a table to stream.
void insert(const scalar, DynamicList< floatScalar > &)
Append scalar to given DynamicList.
conformalVoronoiMesh(const Time &runTime, const dictionary &foamyHexMeshDict)
Construct from Time and foamyHexMeshDict.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
An indexed form of CGAL::Triangulation_vertex_base_3<K> used to keep track of the Delaunay vertices i...
Definition: indexedVertex.H:51
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys.
Definition: HashSet.H:211
void move()
Move the vertices according to the controller, re-conforming to the.
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
bool uninitialised(const VertexType &v)
Type gSum(const FieldField< Field, Type > &f)
line< point, const point & > linePointRef
Line using referred points.
Definition: linePointRef.H:45
const pointField & points
void updateCellSizeAndFaceArea(vector &alignmentDir, scalar &targetFaceArea, scalar &targetCellSize) const
static word timeName(const scalar, const int precision=curPrecision_)
Return time name of given scalar time.
Definition: Time.C:666
bool findSurfaceAnyIntersection(const point &start, const point &end) const
const cellAspectRatioControl & aspectRatio() const
sideVolumeType
Normals point to the outside.
void reset()
Clear the entire triangulation.
PointFrompoint toPoint(const Foam::point &p)
word timeName
Definition: getTimeIndex.H:3
List< label > labelList
A List of labels.
Definition: labelList.H:56
void writeFeatureObj(const fileName &prefix) const
Write all components of all the extendedFeatureEdgeMeshes as.
static const zero Zero
Definition: zero.H:97
errorManip< error > abort(error &err)
Definition: errorManip.H:131
Istream and Ostream manipulators taking arguments.
const backgroundMeshDecomposition & decomposition() const
Return the backgroundMeshDecomposition.
void reverse(UList< T > &, const label n)
Definition: UListI.H:334
static const char nl
Definition: Ostream.H:260
defineTypeNameAndDebug(combustionModel, 0)
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurrence of given element and return index,.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
scalar cosInsertionAcceptanceAngle() const
Return the cosInsertionAcceptanceAngle.
Definition: cvControlsI.H:191
void findSurfaceAllIntersections(const point &start, const point &end, List< pointIndexHit > &surfHit, labelList &hitSurface) const
void writeBoundaryPoints(const fileName &fName, const Triangulation &t)
Write the boundary Delaunay points to an OBJ file.
void resetCellCount()
Set the cell count to zero.
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
Foam::scalar averageCellSize(const VertexType &vA, const VertexType &vB)
Return the target cell size from that stored on a pair of Delaunay vertices,.
#define WarningInFunction
Report a warning using Foam::Warning.
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
void updateDeltaVector(const vector &alignmentDir, const scalar targetCellSize, const scalar rABMag, vector &delta) const
label getNewCellIndex() const
Create a new unique cell index and return.
Definition: DelaunayMeshI.H:63
messageStream Info
static const NamedEnum< dualMeshPointType, 5 > dualMeshPointTypeNames_
dimensioned< scalar > mag(const dimensioned< Type > &)
Field< vector > vectorField
Specialisation of Field<T> for vector.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
const cvControls & foamyHexMeshControls() const
Return the foamyHexMeshControls object.
const Time & time() const
Return the Time object.
scalar coplanarTet(Cell &c, const scalar tol=1e-12)
Map< label > insertPoints(const List< Vb > &vertices, const bool reIndex)
Insert the list of vertices (calls rangeInsertWithInfo)
scalar insertionDistCoeff() const
Return the insertionDistCoeff.
Definition: cvControlsI.H:179
void sync(const boundBox &bb)
Refer vertices so that the processor interfaces are consistent.
InfoProxy< IOstream > info() const
Return info proxy.
Definition: IOstream.H:528
Namespace for OpenFOAM.
A HashTable to objects of type <T> with a label key.
Definition: Map.H:49