conformalVoronoiMeshI.H
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-2019 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 "indexedVertexOps.H"
27 #include "indexedCellOps.H"
28 
29 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
30 
31 inline Foam::scalar Foam::conformalVoronoiMesh::defaultCellSize() const
32 {
34 }
35 
36 
37 inline Foam::scalar Foam::conformalVoronoiMesh::targetCellSize
38 (
39  const Foam::point& pt
40 ) const
41 {
42  return cellShapeControls().cellSize(pt);
43 }
44 
45 
46 inline Foam::scalar Foam::conformalVoronoiMesh::averageAnyCellSize
47 (
48  const Vertex_handle& vA,
49  const Vertex_handle& vB
50 ) const
51 {
52  if
53  (
54  (!vA->internalOrBoundaryPoint() || vA->referred())
55  && (!vB->internalOrBoundaryPoint() || vB->referred())
56  )
57  {
58  // There are no internalOrBoundaryPoints available, determine
59  // size from scratch
60 
61  // Geometric mean
62  return sqrt
63  (
64  targetCellSize(topoint(vA->point()))
65  *targetCellSize(topoint(vB->point()))
66  );
67  }
68  else if (!vB->internalOrBoundaryPoint() || vB->referred())
69  {
70  return vA->targetCellSize();
71  }
72  else if (!vA->internalOrBoundaryPoint() || vA->referred())
73  {
74  return vB->targetCellSize();
75  }
76 
78 }
79 
80 
81 inline Foam::scalar Foam::conformalVoronoiMesh::averageAnyCellSize
82 (
83  const Delaunay::Finite_facets_iterator& fit
84 ) const
85 {
86  // Arithmetic mean
87 
88  scalar sizeSum = 0;
89 
90  label nProducts = 0;
91 
92  const Cell_handle c(fit->first);
93  const label oppositeVertex = fit->second;
94 
95  for (label i = 0; i < 3; i++)
96  {
97  Vertex_handle v = c->vertex(vertex_triple_index(oppositeVertex, i));
98 
99  if (v->internalOrBoundaryPoint() && !v->referred())
100  {
101 
102  sizeSum += v->targetCellSize();
103 
104  nProducts++;
105  }
106  }
107 
108  if (nProducts < 1)
109  {
110  // There are no internalOrBoundaryPoints available, determine
111  // size from scratch
112 
113  for (label i = 0; i < 3; i++)
114  {
115  Vertex_handle v = c->vertex(vertex_triple_index(oppositeVertex, i));
116 
117  sizeSum += targetCellSize(topoint(v->point()));
118  }
119 
120  nProducts = 3;
121  }
122 
123  if (sizeSum < 0)
124  {
126  << "sizeSum = " << sizeSum
127  << endl;
128 
129  return 0;
130  }
131 
132  return pow(sizeSum, (1.0/nProducts));
133 }
134 
135 
137 (
138  const Foam::point& pt
139 ) const
140 {
141  return targetCellSize(pt)*foamyHexMeshControls().pointPairDistanceCoeff();
142 }
143 
144 
146 (
147  const Foam::point& pt
148 ) const
149 {
150  return
153 }
154 
155 
157 (
158  const Foam::point& pt
159 ) const
160 {
161  return
162  sqr
163  (
164  targetCellSize(pt)
165  *foamyHexMeshControls().featurePointExclusionDistanceCoeff()
166  );
167 }
168 
169 
171 (
172  const Foam::point& pt
173 ) const
174 {
175  return
176  sqr
177  (
178  targetCellSize(pt)
179  *foamyHexMeshControls().featureEdgeExclusionDistanceCoeff()
180  );
181 }
182 
183 
185 (
186  const Foam::point& pt
187 ) const
188 {
189  return
190  sqr
191  (
192  targetCellSize(pt)
193  *foamyHexMeshControls().surfacePtExclusionDistanceCoeff()
194  );
195 }
196 
197 
199 (
200  const Foam::point& pt
201 ) const
202 {
203  return
204  sqr
205  (
206  targetCellSize(pt)
207  *foamyHexMeshControls().surfaceSearchDistanceCoeff()
208  );
209 }
210 
211 
213 (
214  const Foam::point& pt
215 ) const
216 {
217  return
218  targetCellSize(pt)
220 }
221 
222 
223 inline void Foam::conformalVoronoiMesh::createPointPair
224 (
225  const scalar ppDist,
226  const Foam::point& surfPt,
227  const vector& n,
228  const bool ptPair,
229  DynamicList<Vb>& pts
230 ) const
231 {
232  vector ppDistn = ppDist*n;
233 
234 // const Foam::point internalPt = surfPt - ppDistn;
235 // const Foam::point externalPt = surfPt + ppDistn;
236 
237 // bool internalInside = geometryToConformTo_.inside(internalPt);
238 // bool externalOutside = geometryToConformTo_.outside(externalPt);
239 
240 // if (internalInside && externalOutside)
241  {
242  pts.append
243  (
244  Vb
245  (
246  surfPt - ppDistn,
247  vertexCount() + pts.size(),
250  )
251  );
252 
253  pts.append
254  (
255  Vb
256  (
257  surfPt + ppDistn,
258  vertexCount() + pts.size(),
261  )
262  );
263 
264  if (ptPair)
265  {
266  ptPairs_.addPointPair
267  (
268  pts[pts.size() - 2].index(),
269  pts[pts.size() - 1].index() // external 0 -> slave
270  );
271  }
272  }
273 // else
274 // {
275 // Info<< "Warning: point pair not inside/outside" << nl
276 // << " surfPt = " << surfPt << nl
277 // << " internal = " << internalPt << " " << internalInside << nl
278 // << " external = " << externalPt << " " << externalOutside
279 // << endl;
280 // }
281 }
282 
283 
284 inline Foam::point Foam::conformalVoronoiMesh::perturbPoint
285 (
286  const Foam::point& pt
287 ) const
288 {
289  Foam::point perturbedPt(pt);
290 
291 // vector delta(xR/ni, yR/nj, zR/nk);
292 // scalar pert = randomPerturbationCoeff*cmptMin(delta);
293 
294  scalar pert = 1e-12*defaultCellSize();
295 
296  perturbedPt.x() += pert*(rndGen_.scalar01() - 0.5);
297  perturbedPt.y() += pert*(rndGen_.scalar01() - 0.5);
298  perturbedPt.z() += pert*(rndGen_.scalar01() - 0.5);
299 
300  return perturbedPt;
301 }
302 
303 
304 inline void Foam::conformalVoronoiMesh::createBafflePointPair
305 (
306  const scalar ppDist,
307  const Foam::point& surfPt,
308  const vector& n,
309  const bool ptPair,
310  DynamicList<Vb>& pts
311 ) const
312 {
313  vector ppDistn = ppDist*n;
314 
315  pts.append
316  (
317  Vb
318  (
319  surfPt - ppDistn,
320  vertexCount() + pts.size(),
323  )
324  );
325 
326  pts.append
327  (
328  Vb
329  (
330  surfPt + ppDistn,
331  vertexCount() + pts.size(),
334  )
335  );
336 
337  if (ptPair)
338  {
339  ptPairs_.addPointPair
340  (
341  pts[pts.size() - 2].index(), // external 0 -> slave
342  pts[pts.size() - 1].index()
343  );
344  }
345 }
346 
347 
348 inline bool Foam::conformalVoronoiMesh::internalPointIsInside
349 (
350  const Foam::point& pt
351 ) const
352 {
353  if
354  (
355  !geometryToConformTo_.globalBounds().contains(pt)
356  || !geometryToConformTo_.inside(pt)
357  )
358  {
359  return false;
360  }
361 
362  return true;
363 }
364 
365 
366 inline bool Foam::conformalVoronoiMesh::isBoundaryDualFace
367 (
368  const Delaunay::Finite_edges_iterator& eit
369 ) const
370 {
371  Cell_handle c = eit->first;
372  Vertex_handle vA = c->vertex(eit->second);
373  Vertex_handle vB = c->vertex(eit->third);
374 
375 // if (vA->internalBoundaryPoint() && vB->externalBoundaryPoint())
376 // {
377 // if (vA->index() == vB->index() - 1)
378 // {
379 // return true;
380 // }
381 // }
382 // else if (vA->externalBoundaryPoint() && vB->internalBoundaryPoint())
383 // {
384 // if (vA->index() == vB->index() + 1)
385 // {
386 // return true;
387 // }
388 // }
389 //
390 // return false;
391 
392  // A dual face on the boundary will result from one Dv inside and
393  // one outside
394  return
395  (
396  (
397  (vA->internalOrBoundaryPoint() && !vA->referred())
398  || (vB->internalOrBoundaryPoint() && !vB->referred())
399  )
400  && (
401  !vA->internalOrBoundaryPoint()
402  || !vB->internalOrBoundaryPoint()
403  )
404  );
405 }
406 
407 
408 inline Foam::List<bool> Foam::conformalVoronoiMesh::dualFaceBoundaryPoints
409 (
410  const Delaunay::Finite_edges_iterator& eit
411 ) const
412 {
413  Cell_circulator ccStart = incident_cells(*eit);
414  Cell_circulator cc1 = ccStart;
415  Cell_circulator cc2 = cc1;
416 
417  // Advance the second circulator so that it always stays on the next
418  // cell around the edge;
419  cc2++;
420 
421  DynamicList<bool> tmpFaceBoundaryPoints;
422 
423  do
424  {
425  label cc1I = cc1->cellIndex();
426 
427  label cc2I = cc2->cellIndex();
428 
429  if (cc1I != cc2I)
430  {
431  if (cc1->boundaryDualVertex())
432  {
433  tmpFaceBoundaryPoints.append(true);
434  }
435  else
436  {
437  tmpFaceBoundaryPoints.append(false);
438  }
439  }
440 
441  cc1++;
442 
443  cc2++;
444 
445  } while (cc1 != ccStart);
446 
447  return Foam::move(tmpFaceBoundaryPoints);
448 }
449 
450 
451 inline Foam::List<Foam::label> Foam::conformalVoronoiMesh::processorsAttached
452 (
453  const Delaunay::Finite_facets_iterator& fit
454 ) const
455 {
456  DynamicList<label> procsAttached(8);
457 
458  const Cell_handle c1(fit->first);
459  const label oppositeVertex = fit->second;
460  const Cell_handle c2(c1->neighbor(oppositeVertex));
461 
463 
465 
466  forAll(c1Procs, aPI)
467  {
468  if (findIndex(procsAttached, c1Procs[aPI] == -1))
469  {
470  procsAttached.append(c1Procs[aPI]);
471  }
472 
473  if (findIndex(procsAttached, c2Procs[aPI] == -1))
474  {
475  procsAttached.append(c2Procs[aPI]);
476  }
477  }
478 
479  return List<label>(procsAttached);
480 }
481 
482 
483 inline bool Foam::conformalVoronoiMesh::isParallelDualEdge
484 (
485  const Delaunay::Finite_facets_iterator& fit
486 ) const
487 {
488  const Cell_handle c1(fit->first);
489  const label oppositeVertex = fit->second;
490 
491  return
492  (
493  c1->vertex(vertex_triple_index(oppositeVertex, 0))->referred()
494  || c1->vertex(vertex_triple_index(oppositeVertex, 1))->referred()
495  || c1->vertex(vertex_triple_index(oppositeVertex, 2))->referred()
496  );
497 }
498 
499 
500 inline bool Foam::conformalVoronoiMesh::isProcBoundaryEdge
501 (
502  const Delaunay::Finite_edges_iterator& eit
503 ) const
504 {
505  bool isProcBoundaryEdge = false;
506 
507  Cell_handle c = eit->first;
508  Vertex_handle vA = c->vertex(eit->second);
509  Vertex_handle vB = c->vertex(eit->third);
510 
511  if
512  (
513  (
514  (vA->referred() && !vB->referred())
515  || (vB->referred() && !vA->referred())
516  )
517  && vA->internalOrBoundaryPoint()
518  && vB->internalOrBoundaryPoint()
519  )
520  {
521  isProcBoundaryEdge = true;
522  }
523 
524  return isProcBoundaryEdge;
525 }
526 
527 
528 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
529 
531 {
532  return runTime_;
533 }
534 
535 
537 {
538  return rndGen_;
539 }
540 
541 
542 inline const Foam::searchableSurfaces&
544 {
545  return allGeometry_;
546 }
547 
548 
549 inline const Foam::conformationSurfaces&
551 {
552  return geometryToConformTo_;
553 }
554 
555 
558 {
559  if (!Pstream::parRun())
560  {
562  << "The backgroundMeshDecomposition cannot be asked for in serial."
563  << exit(FatalError) << endl;
564  }
565 
566  return decomposition_();
567 }
568 
569 
570 inline const Foam::cellShapeControl&
572 {
573  return cellShapeControl_;
574 }
575 
576 
577 inline const Foam::cvControls&
579 {
580  return foamyHexMeshControls_;
581 }
582 
583 
584 // ************************************************************************* //
Field< bool > inside(const pointField &samplePts) const
Check if points are inside surfaces to conform to.
Delaunay::Cell_handle Cell_handle
#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
const searchableSurfaces & allGeometry() const
Return the allGeometry object.
scalar surfaceSearchDistanceSqr(const Foam::point &pt) const
Return the square of the local surface search distance.
pointFromPoint topoint(const Point &P)
scalar pointPairDistance(const Foam::point &pt) const
Return the local point pair separation at the given location.
const treeBoundBox & globalBounds() const
Return the global bounds.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
scalar pointPairDistanceCoeff() const
Return the pointPairDistanceCoeff.
Definition: cvControlsI.H:34
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.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:429
scalar mixedFeaturePointDistance(const Foam::point &pt) const
Return the local mixed feature point placement distance.
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
const Cmpt & z() const
Definition: VectorI.H:87
const conformationSurfaces & geometryToConformTo() const
Return the conformationSurfaces object.
bool addPointPair(const labelPair &vA, const labelPair &vB)
scalar defaultCellSize() const
Return the defaultCellSize.
Definition: cvControlsI.H:137
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:69
const dimensionedScalar c
Speed of light in a vacuum.
const dimensionedScalar c2
Second radiation constant: default SI units: [m K].
const Cmpt & y() const
Definition: VectorI.H:81
Random & rndGen() const
Return the random number generator.
An indexed form of CGAL::Triangulation_vertex_base_3<K> used to keep track of the Delaunay vertices i...
Definition: indexedVertex.H:51
scalar maxSurfaceProtrusionCoeff() const
Return the maxSurfaceProtrusionCoeff.
Definition: cvControlsI.H:95
const dimensionedScalar c1
First radiation constant: default SI units: [W/m^2].
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
Definition: DynamicList.H:56
Foam::FixedList< Foam::label, 4 > processorsAttached(const CellType &c)
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:296
Container for searchableSurfaces.
scalar surfacePtExclusionDistanceSqr(const Foam::point &pt) const
Return the square of the local surface point exclusion distance.
const Cmpt & x() const
Definition: VectorI.H:75
Random number generator.
Definition: Random.H:57
Controls for the conformalVoronoiMesh mesh generator.
Definition: cvControls.H:53
Templated 3D Vector derived from VectorSpace adding construction from 3 components, element access using x(), y() and z() member functions and the inner-product (dot-product) and cross product operators.
Definition: Vector.H:57
bool contains(const vector &dir, const point &) const
Contains point (inside or on edge) and moving in direction.
Definition: treeBoundBox.C:394
const backgroundMeshDecomposition & decomposition() const
Return the backgroundMeshDecomposition.
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurrence of given element and return index,.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
scalar maxSurfaceProtrusion(const Foam::point &pt) const
Return the local maximum surface protrusion distance.
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
label vertexCount() const
Return the vertex count (the next unique vertex index)
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.
label n
Store a background polyMesh to use for the decomposition of space and queries for parallel conformalV...
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 featureEdgeExclusionDistanceSqr(const Foam::point &pt) const
Return the square of the local feature edge exclusion distance.
scalar featurePointExclusionDistanceSqr(const Foam::point &pt) const
Return the square of the local feature point exclusion distance.
scalar scalar01()
Advance the state and return a scalar sample from a uniform.
Definition: RandomI.H:57
scalar mixedFeaturePointPPDistanceCoeff() const
Return the mixedFeaturePointPPDistanceCoeff.
Definition: cvControlsI.H:40
scalar cellSize(const point &pt) const
Return the cell size at the given location.