Loading...
Searching...
No Matches
motionSmootherAlgoCheck.C
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) 2011-2014 OpenFOAM Foundation
9 Copyright (C) 2020,2022 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 "motionSmootherAlgo.H"
30#include "polyMeshGeometry.H"
31#include "IOmanip.H"
32#include "pointSet.H"
33
34// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
35
36bool Foam::motionSmootherAlgo::checkMesh
37(
38 const bool report,
39 const polyMesh& mesh,
40 const dictionary& dict,
41 const labelList& checkFaces,
42 labelHashSet& wrongFaces,
43 const bool dryRun
44)
45{
46 List<labelPair> emptyBaffles;
47 return checkMesh
48 (
49 report,
50 mesh,
51 dict,
52 checkFaces,
53 emptyBaffles,
54 wrongFaces,
55 dryRun
56 );
57}
58
59bool Foam::motionSmootherAlgo::checkMesh
60(
61 const bool report,
62 const polyMesh& mesh,
63 const dictionary& dict,
64 const labelList& checkFaces,
65 const List<labelPair>& baffles,
66 labelHashSet& wrongFaces,
67 const bool dryRun
68)
69{
70 const scalar maxNonOrtho
71 (
72 get<scalar>(dict, "maxNonOrtho", dryRun, keyType::REGEX_RECURSIVE)
73 );
74 const scalar minVol
75 (
76 get<scalar>(dict, "minVol", dryRun, keyType::REGEX_RECURSIVE)
77 );
78 const scalar minTetQuality
79 (
80 get<scalar>(dict, "minTetQuality", dryRun, keyType::REGEX_RECURSIVE)
81 );
82 const scalar maxConcave
83 (
84 get<scalar>(dict, "maxConcave", dryRun, keyType::REGEX_RECURSIVE)
85 );
86 const scalar minArea
87 (
88 get<scalar>(dict, "minArea", dryRun, keyType::REGEX_RECURSIVE)
89 );
90 const scalar maxIntSkew
91 (
92 get<scalar>
93 (
94 dict, "maxInternalSkewness", dryRun, keyType::REGEX_RECURSIVE
95 )
96 );
97 const scalar maxBounSkew
98 (
99 get<scalar>
100 (
101 dict, "maxBoundarySkewness", dryRun, keyType::REGEX_RECURSIVE
102 )
103 );
104 const scalar minWeight
105 (
106 get<scalar>(dict, "minFaceWeight", dryRun, keyType::REGEX_RECURSIVE)
107 );
108 const scalar minVolRatio
109 (
110 get<scalar>(dict, "minVolRatio", dryRun, keyType::REGEX_RECURSIVE)
111 );
112 const scalar minTwist
113 (
114 get<scalar>(dict, "minTwist", dryRun, keyType::REGEX_RECURSIVE)
115 );
116 const scalar minTriangleTwist
117 (
118 get<scalar>(dict, "minTriangleTwist", dryRun, keyType::REGEX_RECURSIVE)
119 );
120 const scalar minFaceFlatness
121 (
122 dict.getOrDefault<scalar>
123 (
124 "minFaceFlatness", -1, keyType::REGEX_RECURSIVE
125 )
126 );
127 const scalar minDet
128 (
129 get<scalar>(dict, "minDeterminant", dryRun, keyType::REGEX_RECURSIVE)
130 );
131 const scalar minEdgeLength
132 (
133 dict.getOrDefault<scalar>
134 (
135 "minEdgeLength", -1, keyType::REGEX_RECURSIVE
136 )
137 );
138
139
140 if (dryRun)
141 {
142 string errorMsg(FatalError.message());
143 string IOerrorMsg(FatalIOError.message());
144
145 if (errorMsg.size() || IOerrorMsg.size())
146 {
147 //errorMsg = "[dryRun] " + errorMsg;
148 //errorMsg.replaceAll("\n", "\n[dryRun] ");
149 //IOerrorMsg = "[dryRun] " + IOerrorMsg;
150 //IOerrorMsg.replaceAll("\n", "\n[dryRun] ");
151
152 IOWarningInFunction(dict)
153 << nl
154 << "Missing/incorrect required dictionary entries:" << nl
155 << nl
156 << IOerrorMsg.c_str() << nl
157 << errorMsg.c_str() << nl
158 //<< nl << "Exiting dry-run" << nl
159 << endl;
160
161 FatalError.clear();
162 FatalIOError.clear();
163 }
164 return false;
165 }
166
167
168 label nWrongFaces = 0;
169
170 Info<< "Checking faces in error :" << endl;
171 //Pout.setf(ios_base::left);
172
173 if (maxNonOrtho < 180.0-SMALL)
174 {
175 polyMeshGeometry::checkFaceDotProduct
176 (
177 report,
178 maxNonOrtho,
179 mesh,
180 mesh.cellCentres(),
181 mesh.faceAreas(),
182 checkFaces,
183 baffles,
184 &wrongFaces
185 );
186
187 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
188
189 Info<< " non-orthogonality > "
190 << setw(3) << maxNonOrtho
191 << " degrees : "
192 << nNewWrongFaces-nWrongFaces << endl;
193
194 nWrongFaces = nNewWrongFaces;
195 }
196
197 if (minVol > -GREAT)
198 {
199 polyMeshGeometry::checkFacePyramids
200 (
201 report,
202 minVol,
203 mesh,
204 mesh.cellCentres(),
205 mesh.points(),
206 checkFaces,
207 baffles,
208 &wrongFaces
209 );
210
211 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
212
213 Info<< " faces with face pyramid volume < "
214 << setw(5) << minVol << " : "
215 << nNewWrongFaces-nWrongFaces << endl;
216
217 nWrongFaces = nNewWrongFaces;
218 }
219
220 if (minTetQuality > -GREAT)
221 {
222 polyMeshGeometry::checkFaceTets
223 (
224 report,
225 minTetQuality,
226 mesh,
227 mesh.cellCentres(),
228 mesh.faceCentres(),
229 mesh.points(),
230 checkFaces,
231 baffles,
232 &wrongFaces
233 );
234
235 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
236
237 Info<< " faces with face-decomposition tet quality < "
238 << setw(5) << minTetQuality << " : "
239 << nNewWrongFaces-nWrongFaces << endl;
240
241 nWrongFaces = nNewWrongFaces;
242 }
243
244 if (maxConcave < 180.0-SMALL)
245 {
246 polyMeshGeometry::checkFaceAngles
247 (
248 report,
249 maxConcave,
250 mesh,
251 mesh.faceAreas(),
252 mesh.points(),
253 checkFaces,
254 &wrongFaces
255 );
256
257 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
258
259 Info<< " faces with concavity > "
260 << setw(3) << maxConcave
261 << " degrees : "
262 << nNewWrongFaces-nWrongFaces << endl;
263
264 nWrongFaces = nNewWrongFaces;
265 }
266
267 if (minArea > -SMALL)
268 {
269 polyMeshGeometry::checkFaceArea
270 (
271 report,
272 minArea,
273 mesh,
274 mesh.faceAreas(),
275 checkFaces,
276 &wrongFaces
277 );
278
279 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
280
281 Info<< " faces with area < "
282 << setw(5) << minArea
283 << " m^2 : "
284 << nNewWrongFaces-nWrongFaces << endl;
285
286 nWrongFaces = nNewWrongFaces;
287 }
288
289 if (maxIntSkew > 0 || maxBounSkew > 0)
290 {
291 polyMeshGeometry::checkFaceSkewness
292 (
293 report,
294 maxIntSkew,
295 maxBounSkew,
296 mesh,
297 mesh.points(),
298 mesh.cellCentres(),
299 mesh.faceCentres(),
300 mesh.faceAreas(),
301 checkFaces,
302 baffles,
303 &wrongFaces
304 );
305
306 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
307
308 Info<< " faces with skewness > "
309 << setw(3) << maxIntSkew
310 << " (internal) or " << setw(3) << maxBounSkew
311 << " (boundary) : " << nNewWrongFaces-nWrongFaces << endl;
312
313 nWrongFaces = nNewWrongFaces;
314 }
315
316 if (minWeight >= 0 && minWeight < 1)
317 {
318 polyMeshGeometry::checkFaceWeights
319 (
320 report,
321 minWeight,
322 mesh,
323 mesh.cellCentres(),
324 mesh.faceCentres(),
325 mesh.faceAreas(),
326 checkFaces,
327 baffles,
328 &wrongFaces
329 );
330
331 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
332
333 Info<< " faces with interpolation weights (0..1) < "
334 << setw(5) << minWeight
335 << " : "
336 << nNewWrongFaces-nWrongFaces << endl;
337
338 nWrongFaces = nNewWrongFaces;
339 }
340
341 if (minVolRatio >= 0)
342 {
343 polyMeshGeometry::checkVolRatio
344 (
345 report,
346 minVolRatio,
347 mesh,
348 mesh.cellVolumes(),
349 checkFaces,
350 baffles,
351 &wrongFaces
352 );
353
354 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
355
356 Info<< " faces with volume ratio of neighbour cells < "
357 << setw(5) << minVolRatio
358 << " : "
359 << nNewWrongFaces-nWrongFaces << endl;
360
361 nWrongFaces = nNewWrongFaces;
362 }
363
364 if (minTwist > -1)
365 {
366 //Pout<< "Checking face twist: dot product of face normal "
367 // << "with face triangle normals" << endl;
368 polyMeshGeometry::checkFaceTwist
369 (
370 report,
371 minTwist,
372 mesh,
373 mesh.cellCentres(),
374 mesh.faceAreas(),
375 mesh.faceCentres(),
376 mesh.points(),
377 checkFaces,
378 &wrongFaces
379 );
380
381 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
382
383 Info<< " faces with face twist < "
384 << setw(5) << minTwist
385 << " : "
386 << nNewWrongFaces-nWrongFaces << endl;
387
388 nWrongFaces = nNewWrongFaces;
389 }
390
391 if (minTriangleTwist > -1)
392 {
393 //Pout<< "Checking triangle twist: dot product of consecutive triangle"
394 // << " normals resulting from face-centre decomposition" << endl;
395 polyMeshGeometry::checkTriangleTwist
396 (
397 report,
398 minTriangleTwist,
399 mesh,
400 mesh.faceAreas(),
401 mesh.faceCentres(),
402 mesh.points(),
403 checkFaces,
404 &wrongFaces
405 );
406
407 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
408
409 Info<< " faces with triangle twist < "
410 << setw(5) << minTriangleTwist
411 << " : "
412 << nNewWrongFaces-nWrongFaces << endl;
413
414 nWrongFaces = nNewWrongFaces;
415 }
416
417 if (minFaceFlatness > -SMALL)
418 {
419 polyMeshGeometry::checkFaceFlatness
420 (
421 report,
422 minFaceFlatness,
423 mesh,
424 mesh.faceAreas(),
425 mesh.faceCentres(),
426 mesh.points(),
427 checkFaces,
428 &wrongFaces
429 );
430
431 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
432
433 Info<< " faces with flatness < "
434 << setw(5) << minFaceFlatness
435 << " : "
436 << nNewWrongFaces-nWrongFaces << endl;
437
438 nWrongFaces = nNewWrongFaces;
439 }
440
441 if (minDet > -1)
442 {
443 polyMeshGeometry::checkCellDeterminant
444 (
445 report,
446 minDet,
447 mesh,
448 mesh.faceAreas(),
449 checkFaces,
450 polyMeshGeometry::affectedCells(mesh, checkFaces),
451 &wrongFaces
452 );
453
454 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
455
456 Info<< " faces on cells with determinant < "
457 << setw(5) << minDet << " : "
458 << nNewWrongFaces-nWrongFaces << endl;
459
460 nWrongFaces = nNewWrongFaces;
461 }
462
463 if (minEdgeLength >= 0)
464 {
465 pointSet edgePoints(mesh, "smallEdgePoints", mesh.nPoints()/1000);
466 polyMeshGeometry::checkEdgeLength
467 (
468 report,
469 minEdgeLength,
470 mesh,
471 checkFaces,
472 &edgePoints
473 );
474
475 const auto& pointFaces = mesh.pointFaces();
476
477 label nNewWrongFaces = 0;
478 for (const label pointi : edgePoints)
479 {
480 const auto& pFaces = pointFaces[pointi];
481 for (const label facei : pFaces)
482 {
483 if (wrongFaces.insert(facei))
484 {
485 nNewWrongFaces++;
486 }
487 }
488 }
489
490 Info<< " faces with edge length < "
491 << setw(5) << minEdgeLength << " : "
492 << returnReduce(nNewWrongFaces, sumOp<label>()) << endl;
493
494 nWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
495 }
496
497 //Pout.setf(ios_base::right);
498
499 return nWrongFaces > 0;
500}
501
502
503bool Foam::motionSmootherAlgo::checkMesh
504(
505 const bool report,
506 const polyMesh& mesh,
507 const dictionary& dict,
508 labelHashSet& wrongFaces,
509 const bool dryRun
510)
511{
512 return checkMesh
513 (
514 report,
515 mesh,
516 dict,
517 identity(mesh.nFaces()),
518 wrongFaces,
519 dryRun
520 );
521}
522
523bool Foam::motionSmootherAlgo::checkMesh
524(
525 const bool report,
526 const dictionary& dict,
527 const polyMeshGeometry& meshGeom,
528 const pointField& points,
529 const labelList& checkFaces,
530 labelHashSet& wrongFaces,
531 const bool dryRun
532)
533{
534 List<labelPair> emptyBaffles;
535
536 return checkMesh
537 (
538 report,
539 dict,
540 meshGeom,
541 points,
542 checkFaces,
543 emptyBaffles,
544 wrongFaces,
545 dryRun
546 );
547}
548
549
550bool Foam::motionSmootherAlgo::checkMesh
551(
552 const bool report,
553 const dictionary& dict,
554 const polyMeshGeometry& meshGeom,
555 const pointField& points,
556 const labelList& checkFaces,
557 const List<labelPair>& baffles,
558 labelHashSet& wrongFaces,
559 const bool dryRun
560)
561{
562 const scalar maxNonOrtho
563 (
564 get<scalar>(dict, "maxNonOrtho", dryRun, keyType::REGEX_RECURSIVE)
565 );
566 const scalar minVol
567 (
568 get<scalar>(dict, "minVol", dryRun, keyType::REGEX_RECURSIVE)
569 );
570 const scalar minTetQuality
571 (
572 get<scalar>(dict, "minTetQuality", dryRun, keyType::REGEX_RECURSIVE)
573 );
574 const scalar maxConcave
575 (
576 get<scalar>(dict, "maxConcave", dryRun, keyType::REGEX_RECURSIVE)
577 );
578 const scalar minArea
579 (
580 get<scalar>(dict, "minArea", dryRun, keyType::REGEX_RECURSIVE)
581 );
582 const scalar maxIntSkew
583 (
584 get<scalar>
585 (
586 dict,
587 "maxInternalSkewness",
588 dryRun,
589 keyType::REGEX_RECURSIVE
590 )
591 );
592 const scalar maxBounSkew
593 (
594 get<scalar>
595 (
596 dict,
597 "maxBoundarySkewness",
598 dryRun,
599 keyType::REGEX_RECURSIVE
600 )
601 );
602 const scalar minWeight
603 (
604 get<scalar>(dict, "minFaceWeight", dryRun, keyType::REGEX_RECURSIVE)
605 );
606 const scalar minVolRatio
607 (
608 get<scalar>(dict, "minVolRatio", dryRun, keyType::REGEX_RECURSIVE)
609 );
610 const scalar minTwist
611 (
612 get<scalar>(dict, "minTwist", dryRun, keyType::REGEX_RECURSIVE)
613 );
614 const scalar minTriangleTwist
615 (
616 get<scalar>(dict, "minTriangleTwist", dryRun, keyType::REGEX_RECURSIVE)
617 );
618 scalar minFaceFlatness = -1.0;
619 dict.readIfPresent
620 (
621 "minFaceFlatness",
622 minFaceFlatness,
623 keyType::REGEX_RECURSIVE
624 );
625 const scalar minDet
626 (
627 get<scalar>(dict, "minDeterminant", dryRun, keyType::REGEX_RECURSIVE)
628 );
629 const scalar minEdgeLength
630 (
631 dict.getOrDefault<scalar>
632 (
633 "minEdgeLength", -1, keyType::REGEX_RECURSIVE
634 )
635 );
636
637 if (dryRun)
638 {
639 string errorMsg(FatalError.message());
640 string IOerrorMsg(FatalIOError.message());
641
642 if (errorMsg.size() || IOerrorMsg.size())
643 {
644 //errorMsg = "[dryRun] " + errorMsg;
645 //errorMsg.replaceAll("\n", "\n[dryRun] ");
646 //IOerrorMsg = "[dryRun] " + IOerrorMsg;
647 //IOerrorMsg.replaceAll("\n", "\n[dryRun] ");
648
649 Perr<< nl
650 << "Missing/incorrect required dictionary entries:" << nl
651 << nl
652 << IOerrorMsg.c_str() << nl
653 << errorMsg.c_str() << nl
654 //<< nl << "Exiting dry-run" << nl
655 << endl;
656
657 FatalError.clear();
658 FatalIOError.clear();
659 }
660 return false;
661 }
662
663
664 label nWrongFaces = 0;
665
666 Info<< "Checking faces in error :" << endl;
667 //Pout.setf(ios_base::left);
668
669 if (maxNonOrtho < 180.0-SMALL)
670 {
671 meshGeom.checkFaceDotProduct
672 (
673 report,
674 maxNonOrtho,
675 checkFaces,
676 baffles,
677 &wrongFaces
678 );
679
680 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
681
682 Info<< " non-orthogonality > "
683 << setw(3) << maxNonOrtho
684 << " degrees : "
685 << nNewWrongFaces-nWrongFaces << endl;
686
687 nWrongFaces = nNewWrongFaces;
688 }
689
690 if (minVol > -GREAT)
691 {
692 meshGeom.checkFacePyramids
693 (
694 report,
695 minTetQuality,
696 points,
697 checkFaces,
698 baffles,
699 &wrongFaces
700 );
701
702 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
703
704 Info<< " faces with face pyramid volume < "
705 << setw(5) << minVol << " : "
706 << nNewWrongFaces-nWrongFaces << endl;
707
708 nWrongFaces = nNewWrongFaces;
709 }
710
711 if (minTetQuality > -GREAT)
712 {
713 meshGeom.checkFaceTets
714 (
715 report,
716 minTetQuality,
717 points,
718 checkFaces,
719 baffles,
720 &wrongFaces
721 );
722
723 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
724
725 Info<< " faces with face-decomposition tet quality < "
726 << setw(5) << minTetQuality << " : "
727 << nNewWrongFaces-nWrongFaces << endl;
728
729 nWrongFaces = nNewWrongFaces;
730 }
731
732 if (maxConcave < 180.0-SMALL)
733 {
734 meshGeom.checkFaceAngles
735 (
736 report,
737 maxConcave,
738 points,
739 checkFaces,
740 &wrongFaces
741 );
742
743 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
744
745 Info<< " faces with concavity > "
746 << setw(3) << maxConcave
747 << " degrees : "
748 << nNewWrongFaces-nWrongFaces << endl;
749
750 nWrongFaces = nNewWrongFaces;
751 }
752
753 if (minArea > -SMALL)
754 {
755 meshGeom.checkFaceArea
756 (
757 report,
758 minArea,
759 checkFaces,
760 &wrongFaces
761 );
762
763 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
764
765 Info<< " faces with area < "
766 << setw(5) << minArea
767 << " m^2 : "
768 << nNewWrongFaces-nWrongFaces << endl;
769
770 nWrongFaces = nNewWrongFaces;
771 }
772
773 if (maxIntSkew > 0 || maxBounSkew > 0)
774 {
775 polyMeshGeometry::checkFaceSkewness
776 (
777 report,
778 maxIntSkew,
779 maxBounSkew,
780 meshGeom.mesh(),
781 points,
782 meshGeom.cellCentres(),
783 meshGeom.faceCentres(),
784 meshGeom.faceAreas(),
785 checkFaces,
786 baffles,
787 &wrongFaces
788 );
789
790 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
791
792 Info<< " faces with skewness > "
793 << setw(3) << maxIntSkew
794 << " (internal) or " << setw(3) << maxBounSkew
795 << " (boundary) : " << nNewWrongFaces-nWrongFaces << endl;
796
797 nWrongFaces = nNewWrongFaces;
798 }
799
800 if (minWeight >= 0 && minWeight < 1)
801 {
802 meshGeom.checkFaceWeights
803 (
804 report,
805 minWeight,
806 checkFaces,
807 baffles,
808 &wrongFaces
809 );
810
811 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
812
813 Info<< " faces with interpolation weights (0..1) < "
814 << setw(5) << minWeight
815 << " : "
816 << nNewWrongFaces-nWrongFaces << endl;
817
818 nWrongFaces = nNewWrongFaces;
819 }
820
821 if (minVolRatio >= 0)
822 {
823 meshGeom.checkVolRatio
824 (
825 report,
826 minVolRatio,
827 checkFaces,
828 baffles,
829 &wrongFaces
830 );
831
832 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
833
834 Info<< " faces with volume ratio of neighbour cells < "
835 << setw(5) << minVolRatio
836 << " : "
837 << nNewWrongFaces-nWrongFaces << endl;
838
839 nWrongFaces = nNewWrongFaces;
840 }
841
842 if (minTwist > -1)
843 {
844 //Pout<< "Checking face twist: dot product of face normal "
845 // << "with face triangle normals" << endl;
846 meshGeom.checkFaceTwist
847 (
848 report,
849 minTwist,
850 points,
851 checkFaces,
852 &wrongFaces
853 );
854
855 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
856
857 Info<< " faces with face twist < "
858 << setw(5) << minTwist
859 << " : "
860 << nNewWrongFaces-nWrongFaces << endl;
861
862 nWrongFaces = nNewWrongFaces;
863 }
864
865 if (minTriangleTwist > -1)
866 {
867 //Pout<< "Checking triangle twist: dot product of consecutive triangle"
868 // << " normals resulting from face-centre decomposition" << endl;
869 meshGeom.checkTriangleTwist
870 (
871 report,
872 minTriangleTwist,
873 points,
874 checkFaces,
875 &wrongFaces
876 );
877
878 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
879
880 Info<< " faces with triangle twist < "
881 << setw(5) << minTriangleTwist
882 << " : "
883 << nNewWrongFaces-nWrongFaces << endl;
884
885 nWrongFaces = nNewWrongFaces;
886 }
887
888 if (minFaceFlatness > -SMALL)
889 {
890 meshGeom.checkFaceFlatness
891 (
892 report,
893 minFaceFlatness,
894 points,
895 checkFaces,
896 &wrongFaces
897 );
898
899 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
900
901 Info<< " faces with flatness < "
902 << setw(5) << minFaceFlatness
903 << " : "
904 << nNewWrongFaces-nWrongFaces << endl;
905
906 nWrongFaces = nNewWrongFaces;
907 }
908
909 if (minDet > -1)
910 {
911 meshGeom.checkCellDeterminant
912 (
913 report,
914 minDet,
915 checkFaces,
916 polyMeshGeometry::affectedCells(meshGeom.mesh(), checkFaces),
917 &wrongFaces
918 );
919
920 label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
921
922 Info<< " faces on cells with determinant < "
923 << setw(5) << minDet << " : "
924 << nNewWrongFaces-nWrongFaces << endl;
925
926 nWrongFaces = nNewWrongFaces;
927 }
928
929 if (minEdgeLength >= 0)
930 {
931 pointSet edgePoints
932 (
933 meshGeom.mesh(),
934 "smallEdgePoints",
935 meshGeom.mesh().nPoints()/1000
936 );
937 meshGeom.checkEdgeLength
938 (
939 report,
940 minEdgeLength,
941 checkFaces,
942 &edgePoints
943 );
944
945 const auto& pointFaces = meshGeom.mesh().pointFaces();
946
947 label nNewWrongFaces = 0;
948 for (const label pointi : edgePoints)
949 {
950 const auto& pFaces = pointFaces[pointi];
951 for (const label facei : pFaces)
952 {
953 if (wrongFaces.insert(facei))
954 {
955 nNewWrongFaces++;
956 }
957 }
958 }
959
960 Info<< " faces with edge length < "
961 << setw(5) << minEdgeLength << " : "
962 << returnReduce(nNewWrongFaces, sumOp<label>()) << endl;
963
964 nWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
965 }
966
967 //Pout.setf(ios_base::right);
968
969 return nWrongFaces > 0;
970}
971
972
973// ************************************************************************* //
IOmanip.H
Istream and Ostream manipulators taking arguments.
Foam::HashSet::insert
bool insert(const Key &key)
Insert a new entry, not overwriting existing entries.
Definition HashSet.H:229
Foam::HashTable::size
label size() const noexcept
The number of elements in table.
Definition HashTable.H:358
Foam::List
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition List.H:72
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition dictionary.H:133
Foam::keyType::REGEX_RECURSIVE
@ REGEX_RECURSIVE
Definition keyType.H:88
Foam::motionSmootherAlgo::get
static Type get(const dictionary &dict, const word &keyword, const bool noExit, enum keyType::option matchOpt, const Type &defaultValue=Zero)
Wrapper around dictionary::get which does not exit.
Definition motionSmootherAlgoTemplates.C:288
Foam::motionSmootherAlgo::mesh
const polyMesh & mesh() const
Reference to mesh.
Definition motionSmootherAlgo.C:354
Foam::motionSmootherAlgo::checkMesh
static bool checkMesh(const bool report, const polyMesh &mesh, const dictionary &dict, labelHashSet &wrongFaces, const bool dryRun=false)
Check mesh with mesh settings in dict. Collects incorrect faces.
Definition motionSmootherAlgoCheck.C:497
Foam::pointSet
A set of point labels.
Definition pointSet.H:50
Foam::polyMeshGeometry
Updateable mesh geometry and checking routines.
Definition polyMeshGeometry.H:52
Foam::polyMeshGeometry::checkFacePyramids
static bool checkFacePyramids(const bool report, const scalar minPyrVol, const polyMesh &, const vectorField &cellCentres, const pointField &p, const labelList &checkFaces, const List< labelPair > &baffles, labelHashSet *)
See primitiveMesh.
Definition polyMeshGeometry.C:565
Foam::polyMeshGeometry::checkTriangleTwist
static bool checkTriangleTwist(const bool report, const scalar minTwist, const polyMesh &, const vectorField &faceAreas, const vectorField &faceCentres, const pointField &p, const labelList &checkFaces, labelHashSet *setPtr)
Consecutive triangle (from face-centre decomposition) normals.
Definition polyMeshGeometry.C:1718
Foam::polyMeshGeometry::checkFaceTwist
static bool checkFaceTwist(const bool report, const scalar minTwist, const polyMesh &, const vectorField &cellCentres, const vectorField &faceAreas, const vectorField &faceCentres, const pointField &p, const labelList &checkFaces, labelHashSet *setPtr)
Triangle (from face-centre decomposition) normal v.s.
Definition polyMeshGeometry.C:1572
Foam::polyMeshGeometry::checkFaceSkewness
static bool checkFaceSkewness(const bool report, const scalar internalSkew, const scalar boundarySkew, const polyMesh &mesh, const pointField &points, const vectorField &cellCentres, const vectorField &faceCentres, const vectorField &faceAreas, const labelList &checkFaces, const List< labelPair > &baffles, labelHashSet *setPtr)
See primitiveMesh.
Definition polyMeshGeometry.C:957
Foam::polyMeshGeometry::faceAreas
const vectorField & faceAreas() const
Definition polyMeshGeometry.H:165
Foam::polyMeshGeometry::mesh
const polyMesh & mesh() const
Definition polyMeshGeometry.H:160
Foam::polyMeshGeometry::checkFaceArea
static bool checkFaceArea(const bool report, const scalar minArea, const polyMesh &, const vectorField &faceAreas, const labelList &checkFaces, labelHashSet *setPtr)
Small faces.
Definition polyMeshGeometry.C:1957
Foam::polyMeshGeometry::checkFaceDotProduct
static bool checkFaceDotProduct(const bool report, const scalar orthWarn, const polyMesh &, const vectorField &cellCentres, const vectorField &faceAreas, const labelList &checkFaces, const List< labelPair > &baffles, labelHashSet *setPtr)
See primitiveMesh.
Definition polyMeshGeometry.C:387
Foam::polyMeshGeometry::checkFaceAngles
static bool checkFaceAngles(const bool report, const scalar maxDeg, const polyMesh &mesh, const vectorField &faceAreas, const pointField &p, const labelList &checkFaces, labelHashSet *setPtr)
See primitiveMesh.
Definition polyMeshGeometry.C:1447
Foam::polyMeshGeometry::checkFaceWeights
static bool checkFaceWeights(const bool report, const scalar warnWeight, const polyMesh &mesh, const vectorField &cellCentres, const vectorField &faceCentres, const vectorField &faceAreas, const labelList &checkFaces, const List< labelPair > &baffles, labelHashSet *setPtr)
Interpolation weights (0.5 for regular mesh).
Definition polyMeshGeometry.C:1163
Foam::polyMeshGeometry::faceCentres
const vectorField & faceCentres() const
Definition polyMeshGeometry.H:169
Foam::polyMeshGeometry::checkVolRatio
static bool checkVolRatio(const bool report, const scalar warnRatio, const polyMesh &mesh, const scalarField &cellVolumes, const labelList &checkFaces, const List< labelPair > &baffles, labelHashSet *setPtr)
Cell volume ratio of neighbouring cells (1 for regular mesh).
Definition polyMeshGeometry.C:1311
Foam::polyMeshGeometry::checkEdgeLength
static bool checkEdgeLength(const bool report, const scalar minEdgeLength, const polyMesh &mesh, const labelList &checkFaces, labelHashSet *pointSetPtr)
Small edges. Optionally collects points of small edges.
Definition polyMeshGeometry.C:2106
Foam::polyMeshGeometry::checkFaceTets
static bool checkFaceTets(const bool report, const scalar minPyrVol, const polyMesh &, const vectorField &cellCentres, const vectorField &faceCentres, const pointField &p, const labelList &checkFaces, const List< labelPair > &baffles, labelHashSet *)
See primitiveMesh.
Definition polyMeshGeometry.C:736
Foam::polyMeshGeometry::checkFaceFlatness
static bool checkFaceFlatness(const bool report, const scalar minFlatness, const polyMesh &, const vectorField &faceAreas, const vectorField &faceCentres, const pointField &p, const labelList &checkFaces, labelHashSet *setPtr)
Area of faces v.s. sum of triangle areas.
Definition polyMeshGeometry.C:1864
Foam::polyMeshGeometry::cellCentres
const vectorField & cellCentres() const
Definition polyMeshGeometry.H:173
Foam::polyMeshGeometry::checkCellDeterminant
static bool checkCellDeterminant(const bool report, const scalar minDet, const polyMesh &, const vectorField &faceAreas, const labelList &checkFaces, const labelList &affectedCells, labelHashSet *setPtr)
Area of internal faces v.s. boundary faces.
Definition polyMeshGeometry.C:2014
Foam::polyMeshGeometry::affectedCells
static labelList affectedCells(const polyMesh &, const labelList &changedFaces)
Helper function: get affected cells from faces.
Definition polyMeshGeometry.C:210
Foam::polyMesh
Mesh consisting of general polyhedral cells.
Definition polyMesh.H:79
Foam::polyMesh::points
virtual const pointField & points() const
Return raw points.
Definition polyMesh.C:1063
Foam::primitiveMesh::faceCentres
const vectorField & faceCentres() const
Definition primitiveMeshFaceCentresAndAreas.C:71
Foam::primitiveMesh::cellVolumes
const scalarField & cellVolumes() const
Definition primitiveMeshCellCentresAndVols.C:90
Foam::primitiveMesh::cellCentres
const vectorField & cellCentres() const
Definition primitiveMeshCellCentresAndVols.C:78
Foam::primitiveMesh::nPoints
label nPoints() const noexcept
Number of mesh points.
Definition primitiveMeshI.H:30
Foam::primitiveMesh::pointFaces
const labelListList & pointFaces() const
Definition primitiveMeshPointFaces.C:27
Foam::primitiveMesh::faceAreas
const vectorField & faceAreas() const
Definition primitiveMeshFaceCentresAndAreas.C:83
mesh
dynamicFvMesh & mesh
Definition createDynamicFvMesh.H:6
points
const pointField & points
Definition gmvOutputHeader.H:1
IOWarningInFunction
#define IOWarningInFunction(ios)
Report an IO warning using Foam::Warning.
Definition messageStream.H:481
Foam::labelList
List< label > labelList
A List of labels.
Definition List.H:62
Foam::labelHashSet
HashSet< label, Hash< label > > labelHashSet
A HashSet of labels, uses label hasher.
Definition HashSet.H:85
Foam::Perr
prefixOSstream Perr
OSstream wrapped stderr (std::cerr) with parallel prefix.
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere).
Foam::returnReduce
T returnReduce(const T &value, BinaryOp bop, const int tag=UPstream::msgType(), const int communicator=UPstream::worldComm)
Perform reduction on a copy, using specified binary operation.
Definition PstreamReduceOps.H:333
Foam::setw
Omanip< int > setw(const int i)
Definition IOmanip.H:199
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition Ostream.H:519
Foam::identity
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i), works like std::iota() but returning a...
Foam::FatalIOError
IOerror FatalIOError
Error stream (stdout output on all processes), with additional 'FOAM FATAL IO ERROR' header text and ...
Foam::FatalError
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
Foam::pointField
vectorField pointField
pointField is a vectorField.
Definition pointFieldFwd.H:38
Foam::nl
constexpr char nl
The newline '\n' character (0x0a).
Definition Ostream.H:50
pFaces
Info<< "Finished reading KIVA file"<< endl;cellShapeList cellShapes(nPoints);labelList cellZoning(nPoints, -1);const cellModel &hex=cellModel::ref(cellModel::HEX);labelList hexLabels(8);label activeCells=0;labelList pointMap(nPoints);forAll(pointMap, i){ pointMap[i]=i;}for(label i=0;i< nPoints;i++){ if(f[i] > 0.0) { hexLabels[0]=i;hexLabels[1]=i1tab[i];hexLabels[2]=i3tab[i1tab[i]];hexLabels[3]=i3tab[i];hexLabels[4]=i8tab[i];hexLabels[5]=i1tab[i8tab[i]];hexLabels[6]=i3tab[i1tab[i8tab[i]]];hexLabels[7]=i3tab[i8tab[i]];cellShapes[activeCells].reset(hex, hexLabels);edgeList edges=cellShapes[activeCells].edges();forAll(edges, ei) { if(edges[ei].mag(points)< SMALL) { label start=pointMap[edges[ei].start()];while(start !=pointMap[start]) { start=pointMap[start];} label end=pointMap[edges[ei].end()];while(end !=pointMap[end]) { end=pointMap[end];} pointMap[start]=pointMap[end]=Foam::min(start, end);} } cellZoning[activeCells]=idreg[i];activeCells++;}}cellShapes.setSize(activeCells);cellZoning.setSize(activeCells);forAll(cellShapes, celli){ cellShape &cs=cellShapes[celli];forAll(cs, i) { cs[i]=pointMap[cs[i]];} cs.collapse();}label bcIDs[11]={-1, 0, 2, 4, -1, 5, -1, 6, 7, 8, 9};constexpr label nBCs=12;const word *kivaPatchTypes[nBCs]={ &wallPolyPatch::typeName, &wallPolyPatch::typeName, &wallPolyPatch::typeName, &wallPolyPatch::typeName, &symmetryPolyPatch::typeName, &wedgePolyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &symmetryPolyPatch::typeName, &oldCyclicPolyPatch::typeName};enum patchTypeNames{ PISTON, VALVE, LINER, CYLINDERHEAD, AXIS, WEDGE, INFLOW, OUTFLOW, PRESIN, PRESOUT, SYMMETRYPLANE, CYCLIC};const char *kivaPatchNames[nBCs]={ "piston", "valve", "liner", "cylinderHead", "axis", "wedge", "inflow", "outflow", "presin", "presout", "symmetryPlane", "cyclic"};List< DynamicList< face > > pFaces[nBCs]
Definition readKivaGrid.H:233
dict
dictionary dict
Definition searchingEngine.H:11