indexedCellI.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-2016 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
29template<class Gt, class Cb>
31(
32 const Foam::globalIndex& globalDelaunayVertexIndices
33) const
34{
35 Foam::tetCell tVGI;
36
37 for (int i = 0; i < 4; i++)
38 {
39 Vertex_handle v = this->vertex(i);
40
41 // Finding the global index of each Delaunay vertex
42 tVGI[i] = globalDelaunayVertexIndices.toGlobal
43 (
45 v->index()
46 );
47 }
48
49 return tVGI;
50}
51
52
53// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
54
55template<class Gt, class Cb>
57:
58 Cb(),
59 index_(ctUnassigned),
60 filterCount_(0)
61{}
62
63
64template<class Gt, class Cb>
66(
68)
69:
70 Cb(v0, v1, v2, v3),
71 index_(ctUnassigned),
72 filterCount_(0)
73{}
74
75
76template<class Gt, class Cb>
78(
83 Cell_handle n0,
84 Cell_handle n1,
85 Cell_handle n2,
87)
88:
89 Cb(v0, v1, v2, v3, n0, n1, n2, n3),
90 index_(ctUnassigned),
91 filterCount_(0)
92{}
93
94
95// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
96
97template<class Gt, class Cb>
99{
100 return index_;
101}
102
103
104template<class Gt, class Cb>
106{
107 return index_;
108}
109
110
111#ifdef CGAL_INEXACT
112
113 template<class Gt, class Cb>
115 {
116 return reinterpret_cast<const Foam::point&>(this->circumcenter());
117 }
118
119#else
120
121 template<class Gt, class Cb>
123 {
124 const typename Gt::Point_3& P = this->circumcenter();
125
126 return Foam::point
127 (
128 CGAL::to_double(P.x()),
129 CGAL::to_double(P.y()),
130 CGAL::to_double(P.z())
131 );
132 }
133
134#endif
135
136
137template<class Gt, class Cb>
139{
140 return index_ == ctUnassigned;
141}
142
143
144template<class Gt, class Cb>
146{
147 return filterCount_;
148}
149
150
151template<class Gt, class Cb>
153{
154 return filterCount_;
155}
156
157
158template<class Gt, class Cb>
160{
161 return
162 (
163 (
164 this->vertex(0)->real()
165 || this->vertex(1)->real()
166 || this->vertex(2)->real()
167 || this->vertex(3)->real()
168 )
169 &&
170 !(
171 this->vertex(0)->farPoint()
172 || this->vertex(1)->farPoint()
173 || this->vertex(2)->farPoint()
174 || this->vertex(3)->farPoint()
175 )
176 );
177}
178
179
180template<class Gt, class Cb>
182{
183 return
184 (
185 this->vertex(0)->farPoint()
186 || this->vertex(1)->farPoint()
187 || this->vertex(2)->farPoint()
188 || this->vertex(3)->farPoint()
189 );
190}
191
192
193template<class Gt, class Cb>
195{
196 return
197 (
198 this->vertex(0)->referred()
199 || this->vertex(1)->referred()
200 || this->vertex(2)->referred()
201 || this->vertex(3)->referred()
202 );
203}
204
205
206template<class Gt, class Cb>
208{
209 return
210 (
211 this->vertex(0)->featurePoint()
212 || this->vertex(1)->featurePoint()
213 || this->vertex(2)->featurePoint()
214 || this->vertex(3)->featurePoint()
215 );
216}
217
218
219template<class Gt, class Cb>
221{
222 return
223 (
224 this->vertex(0)->seedPoint()
225 || this->vertex(1)->seedPoint()
226 || this->vertex(2)->seedPoint()
227 || this->vertex(3)->seedPoint()
228 );
229}
230
231
232template<class Gt, class Cb>
234{
235 return
236 (
237 this->vertex(0)->internalPoint()
238 || this->vertex(1)->internalPoint()
239 || this->vertex(2)->internalPoint()
240 || this->vertex(3)->internalPoint()
241 );
242}
243
244
245template<class Gt, class Cb>
247{
248 return
249 (
250 this->vertex(0)->boundaryPoint()
251 || this->vertex(1)->boundaryPoint()
252 || this->vertex(2)->boundaryPoint()
253 || this->vertex(3)->boundaryPoint()
254 );
255}
256
257
258template<class Gt, class Cb>
260{
261 return
262 (
263 this->vertex(0)->constrained()
264 || this->vertex(1)->constrained()
265 || this->vertex(2)->constrained()
266 || this->vertex(3)->constrained()
267 );
268}
269
270
271template<class Gt, class Cb>
273{
274 return
275 (
276 !this->hasFarPoint()
277 &&
278 (
279 this->vertex(0)->referred()
280 || this->vertex(1)->referred()
281 || this->vertex(2)->referred()
282 || this->vertex(3)->referred()
283 )
284 &&
285 (
286 this->vertex(0)->real()
287 || this->vertex(1)->real()
288 || this->vertex(2)->real()
289 || this->vertex(3)->real()
290 )
291 );
292}
293
294
295template<class Gt, class Cb>
297{
298 Foam::label lowestProc = -1;
299
300 for (int i = 0; i < 4; ++i)
301 {
302 if (this->vertex(i)->referred())
303 {
304 lowestProc = min(lowestProc, this->vertex(i)->procIndex());
305 }
306 }
307
308 return lowestProc;
309}
310
311
312template<class Gt, class Cb>
314(
315 const Foam::globalIndex& globalDelaunayVertexIndices
316) const
317{
318 // tetVertexGlobalIndices
319 Foam::tetCell tVGI
320 = unsortedVertexGlobalIndices(globalDelaunayVertexIndices);
321
322 // bubble sort
323 for (int i = 0; i < tVGI.size(); i++)
324 {
325 for (int j = tVGI.size() - 1 ; j > i; j--)
326 {
327 if (tVGI[j - 1] > tVGI[j])
328 {
329 std::swap(tVGI[j - 1], tVGI[j]);
330 }
331 }
332 }
333
334 return tVGI;
335}
336
337
338template<class Gt, class Cb>
341(
342 const Foam::globalIndex& globalDelaunayVertexIndices
343) const
344{
345 // tetVertexGlobalIndices
346 Foam::tetCell tVGI
347 = unsortedVertexGlobalIndices(globalDelaunayVertexIndices);
348
350
351 // bubble sort
352 for (int i = 0; i < tVGI.size(); i++)
353 {
354 for (int j = tVGI.size() - 1 ; j > i; j--)
355 {
356 if (tVGI[j - 1] > tVGI[j])
357 {
358 std::swap(tVGI[j - 1], tVGI[j]);
359 std::swap(vertexMap[j - 1], vertexMap[j]);
360 }
361 }
362 }
363
364 for (int i = 0; i < 4; i++)
365 {
366 tVGI[i] = vertexMap[i];
367 }
368
369 return std::move(tVGI);
370}
371
372
373template<class Gt, class Cb>
375{
376 return
377 (
378 this->vertex(0)->internalOrBoundaryPoint()
379 || this->vertex(1)->internalOrBoundaryPoint()
380 || this->vertex(2)->internalOrBoundaryPoint()
381 || this->vertex(3)->internalOrBoundaryPoint()
382 );
383}
384
385
386template<class Gt, class Cb>
388{
389 return
390 (
391 this->vertex(0)->internalOrBoundaryPoint()
392 || this->vertex(0)->externalBoundaryPoint()
393 || this->vertex(1)->internalOrBoundaryPoint()
394 || this->vertex(1)->externalBoundaryPoint()
395 || this->vertex(2)->internalOrBoundaryPoint()
396 || this->vertex(2)->externalBoundaryPoint()
397 || this->vertex(3)->internalOrBoundaryPoint()
398 || this->vertex(3)->externalBoundaryPoint()
399 );
400}
401
402
403template<class Gt, class Cb>
405{
406// return
407// (
408// this->vertex(0)->boundaryPoint()
409// && this->vertex(1)->boundaryPoint()
410// && this->vertex(2)->boundaryPoint()
411// && this->vertex(3)->boundaryPoint()
412// );
413
414 return
415 (
416 (
417 this->vertex(0)->internalBoundaryPoint()
418 || this->vertex(1)->internalBoundaryPoint()
419 || this->vertex(2)->internalBoundaryPoint()
420 || this->vertex(3)->internalBoundaryPoint()
421 )
422 && (
423 this->vertex(0)->externalBoundaryPoint()
424 || this->vertex(1)->externalBoundaryPoint()
425 || this->vertex(2)->externalBoundaryPoint()
426 || this->vertex(3)->externalBoundaryPoint()
427 )
428 );
429
430// Foam::label nBoundaryPoints = 0;
431//
432// for (Foam::label i = 0; i < 4; ++i)
433// {
434// Vertex_handle v = this->vertex(i);
435//
436// if (v->boundaryPoint())
437// {
438// nBoundaryPoints++;
439// }
440// }
441//
442// return (nBoundaryPoints > 1);
443}
444
445
446template<class Gt, class Cb>
448{
449 return
450 (
451 (
452 this->vertex(0)->internalBaffleSurfacePoint()
453 || this->vertex(1)->internalBaffleSurfacePoint()
454 || this->vertex(2)->internalBaffleSurfacePoint()
455 || this->vertex(3)->internalBaffleSurfacePoint()
456 )
457 && (
458 this->vertex(0)->externalBaffleSurfacePoint()
459 || this->vertex(1)->externalBaffleSurfacePoint()
460 || this->vertex(2)->externalBaffleSurfacePoint()
461 || this->vertex(3)->externalBaffleSurfacePoint()
462 )
463 );
464}
465
466
467template<class Gt, class Cb>
469{
470 return
471 (
472 (
473 this->vertex(0)->internalBaffleEdgePoint()
474 || this->vertex(1)->internalBaffleEdgePoint()
475 || this->vertex(2)->internalBaffleEdgePoint()
476 || this->vertex(3)->internalBaffleEdgePoint()
477 )
478 && (
479 this->vertex(0)->externalBaffleEdgePoint()
480 || this->vertex(1)->externalBaffleEdgePoint()
481 || this->vertex(2)->externalBaffleEdgePoint()
482 || this->vertex(3)->externalBaffleEdgePoint()
483 )
484 );
485}
486
487
488template<class Gt, class Cb>
490{
491 return
492 (
493 this->vertex(0)->featureEdgePoint()
494 && this->vertex(1)->featureEdgePoint()
495 && this->vertex(2)->featureEdgePoint()
496 && this->vertex(3)->featureEdgePoint()
497 );
498// (
499// this->vertex(0)->featureEdgePoint()
500// || this->vertex(1)->featureEdgePoint()
501// || this->vertex(2)->featureEdgePoint()
502// || this->vertex(3)->featureEdgePoint()
503// )
504// && (
505// (
506// this->vertex(0)->featureEdgePoint()
507// || this->vertex(0)->featurePoint()
508// )
509// && (
510// this->vertex(1)->featureEdgePoint()
511// || this->vertex(1)->featurePoint()
512// )
513// && (
514// this->vertex(2)->featureEdgePoint()
515// || this->vertex(2)->featurePoint()
516// )
517// && (
518// this->vertex(3)->featureEdgePoint()
519// || this->vertex(3)->featurePoint()
520// )
521// )
522// );
523}
524
525
526template<class Gt, class Cb>
528{
529 return
530 (
531 this->vertex(0)->featurePoint()
532 && this->vertex(1)->featurePoint()
533 && this->vertex(2)->featurePoint()
534 && this->vertex(3)->featurePoint()
535 );
536}
537
538
539template<class Gt, class Cb>
541{
542 return
543 (
544 this->vertex(0)->nearProcBoundary()
545 || this->vertex(1)->nearProcBoundary()
546 || this->vertex(2)->nearProcBoundary()
547 || this->vertex(3)->nearProcBoundary()
548 );
549}
550
551
552template<class Gt, class Cb>
554{
555 Foam::label nMasters = 0;
556 Foam::label nSlaves = 0;
557
558 Vertex_handle vM[2];
559 Vertex_handle vS[2];
560
561 for (Foam::label i = 0; i < 4; ++i)
562 {
563 Vertex_handle v = this->vertex(i);
564
565 if (v->internalBoundaryPoint())
566 {
567 vM[nMasters] = v;
568 nMasters++;
569 }
570
571 if (v->externalBoundaryPoint())
572 {
573 vS[nSlaves] = v;
574 nSlaves++;
575 }
576 }
577
578 Foam::label nPairs = 0;
579
580 if (nMasters == 2 && nSlaves == 2)
581 {
584
585 if
586 (
587 vM[0]->type() == vS[0]->index()
588 && vM[0]->index() == vS[0]->type()
589 )
590 {
591 vp0 = reinterpret_cast<const Foam::point&>(vM[0]->point())
592 - reinterpret_cast<const Foam::point&>(vS[0]->point());
593 nPairs++;
594 }
595 else if
596 (
597 vM[0]->type() == vS[1]->index()
598 && vM[0]->index() == vS[1]->type()
599 )
600 {
601 vp0 = reinterpret_cast<const Foam::point&>(vM[0]->point())
602 - reinterpret_cast<const Foam::point&>(vS[1]->point());
603 nPairs++;
604 }
605
606 if
607 (
608 vM[1]->type() == vS[0]->index()
609 && vM[1]->index() == vS[0]->type()
610 )
611 {
612 vp1 = reinterpret_cast<const Foam::point&>(vM[1]->point())
613 - reinterpret_cast<const Foam::point&>(vS[0]->point());
614 nPairs++;
615 }
616 else if
617 (
618 vM[1]->type() == vS[1]->index()
619 && vM[1]->index() == vS[1]->type()
620 )
621 {
622 vp1 = reinterpret_cast<const Foam::point&>(vM[1]->point())
623 - reinterpret_cast<const Foam::point&>(vS[1]->point());
624 nPairs++;
625 }
626
627 if (nPairs == 2)
628 {
629 if (Foam::vectorTools::areParallel(vp0, vp1))
630 {
631 Foam::Pout<< "PARALLEL" << Foam::endl;
632
633 return true;
634 }
635 }
636 }
637
638 return false;
639}
640
641
642template<class Gt, class Cb>
644{
645 int featureVertex = -1;
646 for (int i = 0; i < 4; ++i)
647 {
648 if (this->vertex(i)->constrained())
649 {
650 featureVertex = i;
651 }
652 }
653
654 // Pick cell with a face attached to an infinite cell
655
656 if (featureVertex != -1)
657 {
658 Vertex_handle v1 =
659 this->vertex(Tds::vertex_triple_index(featureVertex, 0));
660 Vertex_handle v2 =
661 this->vertex(Tds::vertex_triple_index(featureVertex, 1));
662 Vertex_handle v3 =
663 this->vertex(Tds::vertex_triple_index(featureVertex, 2));
664
665 if (v1->internalBoundaryPoint())
666 {
667 if
668 (
669 v2->externalBoundaryPoint()
670 && v3->externalBoundaryPoint()
671 )
672 {
673 return true;
674 }
675 }
676 else if (v2->internalBoundaryPoint())
677 {
678 if
679 (
680 v1->externalBoundaryPoint()
681 && v3->externalBoundaryPoint()
682 )
683 {
684 return true;
685 }
686 }
687 else if (v3->internalBoundaryPoint())
688 {
689 if
690 (
691 v1->externalBoundaryPoint()
692 && v2->externalBoundaryPoint()
693 )
694 {
695 return true;
696 }
697 }
698 }
699
700 return false;
701}
702
703
704template<class Gt, class Cb>
706{
707 int featureVertex = -1;
708 for (int i = 0; i < 4; ++i)
709 {
710 if (this->vertex(i)->constrained())
711 {
712 featureVertex = i;
713 }
714 }
715
716 // Pick cell with a face attached to an infinite cell
717
718 if (featureVertex != -1)
719 {
720 Vertex_handle v1 =
721 this->vertex(Tds::vertex_triple_index(featureVertex, 0));
722 Vertex_handle v2 =
723 this->vertex(Tds::vertex_triple_index(featureVertex, 1));
724 Vertex_handle v3 =
725 this->vertex(Tds::vertex_triple_index(featureVertex, 2));
726
727 if (v1->externalBoundaryPoint())
728 {
729 if
730 (
731 v2->internalBoundaryPoint()
732 && v3->internalBoundaryPoint()
733 )
734 {
735 return true;
736 }
737 }
738 else if (v2->externalBoundaryPoint())
739 {
740 if
741 (
742 v1->internalBoundaryPoint()
743 && v3->internalBoundaryPoint()
744 )
745 {
746 return true;
747 }
748 }
749 else if (v3->externalBoundaryPoint())
750 {
751 if
752 (
753 v1->internalBoundaryPoint()
754 && v2->internalBoundaryPoint()
755 )
756 {
757 return true;
758 }
759 }
760 }
761
762 return false;
763}
764
765
766// * * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * //
An indexed form of CGAL::Triangulation_cell_base_3<K> used to keep track of the Delaunay cells (tets)...
Definition: indexedCell.H:95
bool featurePointExternalCell() const
Definition: indexedCellI.H:643
Foam::FixedList< Foam::label, 4 > globallyOrderedCellVertices(const Foam::globalIndex &globalDelaunayVertexIndices) const
Using the globalIndex object, return a list of four vertices with.
Definition: indexedCellI.H:341
Foam::label vertexLowestProc() const
Definition: indexedCellI.H:296
bool internalOrBoundaryDualVertex() const
Is the Delaunay cell part of the final dual mesh, i.e. any vertex.
Definition: indexedCellI.H:374
bool hasSeedPoint() const
Does the Delaunay cell have a seed point.
Definition: indexedCellI.H:220
bool hasFeaturePoint() const
Does the Delaunay cell have a feature point.
Definition: indexedCellI.H:207
bool parallelDualVertex() const
Does the Dual vertex form part of a processor patch.
Definition: indexedCellI.H:272
Foam::tetCell vertexGlobalIndices(const Foam::globalIndex &globalDelaunayVertexIndices) const
Using the globalIndex object, return a list of four (sorted) global.
Definition: indexedCellI.H:314
bool baffleEdgeDualVertex() const
Definition: indexedCellI.H:468
bool hasFarPoint() const
Does the Delaunay cell have a far point.
Definition: indexedCellI.H:181
bool hasBoundaryPoint() const
Definition: indexedCellI.H:246
bool hasInternalPoint() const
Definition: indexedCellI.H:233
bool anyInternalOrBoundaryDualVertex() const
Is the Delaunay cell real or referred (or mixed), i.e. all vertices.
Definition: indexedCellI.H:387
bool featurePointInternalCell() const
Definition: indexedCellI.H:705
const Foam::point dual()
Definition: indexedCellI.H:122
bool potentialCoplanarCell() const
Definition: indexedCellI.H:553
Cb::Cell_handle Cell_handle
Definition: indexedCell.H:126
bool baffleSurfaceDualVertex() const
Definition: indexedCellI.H:447
Foam::label & cellIndex()
Definition: indexedCellI.H:98
Cb::Vertex_handle Vertex_handle
Definition: indexedCell.H:125
bool featurePointDualVertex() const
A dual vertex on a feature point will result from this Delaunay cell.
Definition: indexedCellI.H:527
bool hasReferredPoint() const
Does the Delaunay cell have a referred point.
Definition: indexedCellI.H:194
bool unassigned() const
Definition: indexedCellI.H:138
bool boundaryDualVertex() const
A dual vertex on the boundary will result from a Delaunay cell with.
Definition: indexedCellI.H:404
bool hasConstrainedPoint() const
Definition: indexedCellI.H:259
bool nearProcBoundary() const
Definition: indexedCellI.H:540
bool featureEdgeDualVertex() const
A dual vertex on a feature edge will result from this Delaunay cell.
Definition: indexedCellI.H:489
bool real() const
Is the Delaunay cell real, i.e. any real vertex.
Definition: indexedCellI.H:159
A 1D vector of objects of type <T> with a fixed length <N>.
Definition: FixedList.H:81
static constexpr label size() noexcept
Return the number of elements in the FixedList.
Definition: FixedList.H:416
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:68
label toGlobal(const label i) const
From local to global index.
Definition: globalIndexI.H:260
int myProcNo() const noexcept
Return processor number.
A tetrahedral cell primitive.
Definition: tetCell.H:67
bool areParallel(const Vector< T > &a, const Vector< T > &b, const T &tolerance=SMALL)
Test if a and b are parallel: a^b = 0.
Definition: vectorTools.H:56
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i)
Definition: labelList.C:38
vector point
Point is a vector.
Definition: point.H:43
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.