conformalVoronoiMeshI.H
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2012-2015 OpenFOAM Foundation
9 Copyright (C) 2021 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "indexedVertexOps.H"
30#include "indexedCellOps.H"
31
32// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
33
34inline Foam::scalar Foam::conformalVoronoiMesh::defaultCellSize() const
35{
37}
38
39
41(
42 const Foam::point& pt
43) const
44{
45 return cellShapeControls().cellSize(pt);
46}
47
48
49inline Foam::scalar Foam::conformalVoronoiMesh::averageAnyCellSize
50(
51 const Vertex_handle& vA,
52 const Vertex_handle& vB
53) const
54{
55 if
56 (
57 (!vA->internalOrBoundaryPoint() || vA->referred())
58 && (!vB->internalOrBoundaryPoint() || vB->referred())
59 )
60 {
61 // There are no internalOrBoundaryPoints available, determine
62 // size from scratch
63
64 // Geometric mean
65 return sqrt
66 (
67 targetCellSize(topoint(vA->point()))
68 *targetCellSize(topoint(vB->point()))
69 );
70 }
71 else if (!vB->internalOrBoundaryPoint() || vB->referred())
72 {
73 return vA->targetCellSize();
74 }
75 else if (!vA->internalOrBoundaryPoint() || vA->referred())
76 {
77 return vB->targetCellSize();
78 }
79
81}
82
83
84inline Foam::scalar Foam::conformalVoronoiMesh::averageAnyCellSize
85(
86 const Delaunay::Finite_facets_iterator& fit
87) const
88{
89 // Arithmetic mean
90
91 scalar sizeSum = 0;
92
93 label nProducts = 0;
94
95 const Cell_handle c(fit->first);
96 const label oppositeVertex = fit->second;
97
98 for (label i = 0; i < 3; i++)
99 {
100 Vertex_handle v = c->vertex(vertex_triple_index(oppositeVertex, i));
101
102 if (v->internalOrBoundaryPoint() && !v->referred())
103 {
104
105 sizeSum += v->targetCellSize();
106
107 nProducts++;
108 }
109 }
110
111 if (nProducts < 1)
112 {
113 // There are no internalOrBoundaryPoints available, determine
114 // size from scratch
115
116 for (label i = 0; i < 3; i++)
117 {
118 Vertex_handle v = c->vertex(vertex_triple_index(oppositeVertex, i));
119
120 sizeSum += targetCellSize(topoint(v->point()));
121 }
122
123 nProducts = 3;
124 }
125
126 if (sizeSum < 0)
127 {
129 << "sizeSum = " << sizeSum
130 << endl;
131
132 return 0;
133 }
134
135 return pow(sizeSum, (1.0/nProducts));
136}
137
138
140(
141 const Foam::point& pt
142) const
143{
144 return targetCellSize(pt)*foamyHexMeshControls().pointPairDistanceCoeff();
145}
146
147
149(
150 const Foam::point& pt
151) const
152{
153 return
154 pointPairDistance(pt)
155 *foamyHexMeshControls().mixedFeaturePointPPDistanceCoeff();
156}
157
158
160(
161 const Foam::point& pt
162) const
163{
164 return
165 sqr
166 (
167 targetCellSize(pt)
168 *foamyHexMeshControls().featurePointExclusionDistanceCoeff()
169 );
170}
171
172
174(
175 const Foam::point& pt
176) const
177{
178 return
179 sqr
180 (
181 targetCellSize(pt)
182 *foamyHexMeshControls().featureEdgeExclusionDistanceCoeff()
183 );
184}
185
186
188(
189 const Foam::point& pt
190) const
191{
192 return
193 sqr
194 (
195 targetCellSize(pt)
196 *foamyHexMeshControls().surfacePtExclusionDistanceCoeff()
197 );
198}
199
200
202(
203 const Foam::point& pt
204) const
205{
206 return
207 sqr
208 (
209 targetCellSize(pt)
210 *foamyHexMeshControls().surfaceSearchDistanceCoeff()
211 );
212}
213
214
216(
217 const Foam::point& pt
218) const
219{
220 return
221 targetCellSize(pt)
222 *foamyHexMeshControls().maxSurfaceProtrusionCoeff();
223}
224
225
226inline void Foam::conformalVoronoiMesh::createPointPair
227(
228 const scalar ppDist,
229 const Foam::point& surfPt,
230 const vector& n,
231 const bool ptPair,
232 DynamicList<Vb>& pts
233) const
234{
235 vector ppDistn = ppDist*n;
236
237// const Foam::point internalPt = surfPt - ppDistn;
238// const Foam::point externalPt = surfPt + ppDistn;
239
240// bool internalInside = geometryToConformTo_.inside(internalPt);
241// bool externalOutside = geometryToConformTo_.outside(externalPt);
242
243// if (internalInside && externalOutside)
244 {
245 pts.append
246 (
247 Vb
248 (
249 surfPt - ppDistn,
250 vertexCount() + pts.size(),
253 )
254 );
255
256 pts.append
257 (
258 Vb
259 (
260 surfPt + ppDistn,
261 vertexCount() + pts.size(),
264 )
265 );
266
267 if (ptPair)
268 {
269 ptPairs_.addPointPair
270 (
271 pts[pts.size() - 2].index(),
272 pts[pts.size() - 1].index() // external 0 -> slave
273 );
274 }
275 }
276// else
277// {
278// Info<< "Warning: point pair not inside/outside" << nl
279// << " surfPt = " << surfPt << nl
280// << " internal = " << internalPt << " " << internalInside << nl
281// << " external = " << externalPt << " " << externalOutside
282// << endl;
283// }
284}
285
286
287inline Foam::point Foam::conformalVoronoiMesh::perturbPoint
288(
289 const Foam::point& pt
290) const
291{
292 Foam::point perturbedPt(pt);
293
294// vector delta(xR/ni, yR/nj, zR/nk);
295// scalar pert = randomPerturbationCoeff*cmptMin(delta);
296
297 scalar pert = 1e-12*defaultCellSize();
298
299 perturbedPt.x() += pert*(rndGen_.sample01<scalar>() - 0.5);
300 perturbedPt.y() += pert*(rndGen_.sample01<scalar>() - 0.5);
301 perturbedPt.z() += pert*(rndGen_.sample01<scalar>() - 0.5);
302
303 return perturbedPt;
304}
305
306
307inline void Foam::conformalVoronoiMesh::createBafflePointPair
308(
309 const scalar ppDist,
310 const Foam::point& surfPt,
311 const vector& n,
312 const bool ptPair,
313 DynamicList<Vb>& pts
314) const
315{
316 vector ppDistn = ppDist*n;
317
318 pts.append
319 (
320 Vb
321 (
322 surfPt - ppDistn,
323 vertexCount() + pts.size(),
326 )
327 );
328
329 pts.append
330 (
331 Vb
332 (
333 surfPt + ppDistn,
334 vertexCount() + pts.size(),
337 )
338 );
339
340 if (ptPair)
341 {
342 ptPairs_.addPointPair
343 (
344 pts[pts.size() - 2].index(), // external 0 -> slave
345 pts[pts.size() - 1].index()
346 );
347 }
348}
349
350
351inline bool Foam::conformalVoronoiMesh::internalPointIsInside
352(
353 const Foam::point& pt
354) const
355{
356 if
357 (
358 !geometryToConformTo_.globalBounds().contains(pt)
359 || !geometryToConformTo_.inside(pt)
360 )
361 {
362 return false;
363 }
364
365 return true;
366}
367
368
369inline bool Foam::conformalVoronoiMesh::isBoundaryDualFace
370(
371 const Delaunay::Finite_edges_iterator& eit
372) const
373{
374 Cell_handle c = eit->first;
375 Vertex_handle vA = c->vertex(eit->second);
376 Vertex_handle vB = c->vertex(eit->third);
377
378// if (vA->internalBoundaryPoint() && vB->externalBoundaryPoint())
379// {
380// if (vA->index() == vB->index() - 1)
381// {
382// return true;
383// }
384// }
385// else if (vA->externalBoundaryPoint() && vB->internalBoundaryPoint())
386// {
387// if (vA->index() == vB->index() + 1)
388// {
389// return true;
390// }
391// }
392//
393// return false;
394
395 // A dual face on the boundary will result from one Dv inside and
396 // one outside
397 return
398 (
399 (
400 (vA->internalOrBoundaryPoint() && !vA->referred())
401 || (vB->internalOrBoundaryPoint() && !vB->referred())
402 )
403 && (
404 !vA->internalOrBoundaryPoint()
405 || !vB->internalOrBoundaryPoint()
406 )
407 );
408}
409
410
411inline Foam::List<bool> Foam::conformalVoronoiMesh::dualFaceBoundaryPoints
412(
413 const Delaunay::Finite_edges_iterator& eit
414) const
415{
416 Cell_circulator ccStart = incident_cells(*eit);
417 Cell_circulator cc1 = ccStart;
418 Cell_circulator cc2 = cc1;
419
420 // Advance the second circulator so that it always stays on the next
421 // cell around the edge;
422 cc2++;
423
424 DynamicList<bool> tmpFaceBoundaryPoints;
425
426 do
427 {
428 label cc1I = cc1->cellIndex();
429
430 label cc2I = cc2->cellIndex();
431
432 if (cc1I != cc2I)
433 {
434 if (cc1->boundaryDualVertex())
435 {
436 tmpFaceBoundaryPoints.append(true);
437 }
438 else
439 {
440 tmpFaceBoundaryPoints.append(false);
441 }
442 }
443
444 cc1++;
445
446 cc2++;
447
448 } while (cc1 != ccStart);
449
450 return tmpFaceBoundaryPoints;
451}
452
453
454inline Foam::List<Foam::label> Foam::conformalVoronoiMesh::processorsAttached
455(
456 const Delaunay::Finite_facets_iterator& fit
457) const
458{
459 DynamicList<label> procsAttached(8);
460
461 const Cell_handle c1(fit->first);
462 const label oppositeVertex = fit->second;
463 const Cell_handle c2(c1->neighbor(oppositeVertex));
464
465 FixedList<label, 4> c1Procs(CGAL::indexedCellOps::processorsAttached(c1));
466
467 FixedList<label, 4> c2Procs(CGAL::indexedCellOps::processorsAttached(c2));
468
469 forAll(c1Procs, aPI)
470 {
471 procsAttached.appendUniq(c1Procs[aPI]);
472 procsAttached.appendUniq(c2Procs[aPI]);
473 }
474
475 return List<label>(procsAttached);
476}
477
478
479inline bool Foam::conformalVoronoiMesh::isParallelDualEdge
480(
481 const Delaunay::Finite_facets_iterator& fit
482) const
483{
484 const Cell_handle c1(fit->first);
485 const label oppositeVertex = fit->second;
486
487 return
488 (
489 c1->vertex(vertex_triple_index(oppositeVertex, 0))->referred()
490 || c1->vertex(vertex_triple_index(oppositeVertex, 1))->referred()
491 || c1->vertex(vertex_triple_index(oppositeVertex, 2))->referred()
492 );
493}
494
495
496inline bool Foam::conformalVoronoiMesh::isProcBoundaryEdge
497(
498 const Delaunay::Finite_edges_iterator& eit
499) const
500{
501 bool isProcBoundaryEdge = false;
502
503 Cell_handle c = eit->first;
504 Vertex_handle vA = c->vertex(eit->second);
505 Vertex_handle vB = c->vertex(eit->third);
506
507 if
508 (
509 (
510 (vA->referred() && !vB->referred())
511 || (vB->referred() && !vA->referred())
512 )
513 && vA->internalOrBoundaryPoint()
514 && vB->internalOrBoundaryPoint()
515 )
516 {
517 isProcBoundaryEdge = true;
518 }
519
520 return isProcBoundaryEdge;
521}
522
523
524// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
525
527{
528 return runTime_;
529}
530
531
533{
534 return rndGen_;
535}
536
537
538inline const Foam::searchableSurfaces&
540{
541 return allGeometry_;
542}
543
544
545inline const Foam::conformationSurfaces&
547{
548 return geometryToConformTo_;
549}
550
551
554{
555 if (!Pstream::parRun())
556 {
558 << "The backgroundMeshDecomposition cannot be asked for in serial."
559 << exit(FatalError) << endl;
560 }
561
562 return *decomposition_;
563}
564
565
566inline const Foam::cellShapeControl&
568{
569 return cellShapeControl_;
570}
571
572
573inline const Foam::cvControls&
575{
576 return foamyHexMeshControls_;
577}
578
579
580// ************************************************************************* //
label n
An indexed form of CGAL::Triangulation_vertex_base_3<K> used to keep track of the Delaunay vertices i...
Definition: indexedVertex.H:88
Foam::scalar & targetCellSize()
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
Definition: DynamicList.H:72
void append(const T &val)
Copy append an element to the end of this list.
Definition: DynamicListI.H:503
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:77
Random number generator.
Definition: Random.H:60
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:80
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
static bool & parRun() noexcept
Test if this a parallel run.
Definition: UPstream.H:433
Store a background polyMesh to use for the decomposition of space and queries for parallel conformalV...
scalar featurePointExclusionDistanceSqr(const Foam::point &pt) const
Return the square of the local feature point exclusion distance.
const conformationSurfaces & geometryToConformTo() const
Return the conformationSurfaces object.
const Time & time() const
Return the Time object.
scalar mixedFeaturePointDistance(const Foam::point &pt) const
Return the local mixed feature point placement distance.
const cvControls & foamyHexMeshControls() const
Return the foamyHexMeshControls object.
const backgroundMeshDecomposition & decomposition() const
Return the backgroundMeshDecomposition.
scalar pointPairDistance(const Foam::point &pt) const
Return the local point pair separation at the given location.
const searchableSurfaces & allGeometry() const
Return the allGeometry object.
scalar surfacePtExclusionDistanceSqr(const Foam::point &pt) const
Return the square of the local surface point exclusion distance.
const cellShapeControl & cellShapeControls() const
Return the cellShapeControl object.
Random & rndGen() const
Return the random number generator.
scalar surfaceSearchDistanceSqr(const Foam::point &pt) const
Return the square of the local surface search distance.
scalar maxSurfaceProtrusion(const Foam::point &pt) const
Return the local maximum surface protrusion distance.
scalar featureEdgeExclusionDistanceSqr(const Foam::point &pt) const
Return the square of the local feature edge exclusion distance.
Controls for the conformalVoronoiMesh mesh generator.
Definition: cvControls.H:56
scalar defaultCellSize() const
Return the defaultCellSize.
Definition: cvControlsI.H:139
int myProcNo() const noexcept
Return processor number.
Container for searchableSurfaces. The collection is specified as a dictionary. For example,...
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
#define WarningInFunction
Report a warning using Foam::Warning.
Foam::FixedList< Foam::label, 4 > processorsAttached(const CellType &c)
Foam::scalar averageCellSize(const VertexType &vA, const VertexType &vB)
Return the target cell size from that stored on a pair of Delaunay vertices,.
const dimensionedScalar c2
Second radiation constant: default SI units: [m.K].
const dimensionedScalar c1
First radiation constant: default SI units: [W/m2].
const dimensionedScalar c
Speed of light in a vacuum.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
pointFromPoint topoint(const Point &P)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
dimensionedScalar sqrt(const dimensionedScalar &ds)
error FatalError
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
volScalarField & e
Definition: createFields.H:11
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333