19 const wedgePolyPatch& wpp
22 const polyBoundaryMesh& patches = mesh.boundaryMesh();
24 scalar wppCosAngle = wpp.cosAngle();
32 && isA<wedgePolyPatch>(patches[
patchi])
35 const wedgePolyPatch& pp =
36 refCast<const wedgePolyPatch>(patches[
patchi]);
39 scalar ppCosAngle = wpp.centreNormal() & pp.n();
43 pp.size() == wpp.size()
44 &&
mag(pp.axis() & wpp.axis()) >= (1-1
e-3)
45 &&
mag(ppCosAngle - wppCosAngle) >= 1
e-3
60 const Vector<label>& directions,
65 EdgeMap<label> edgesInError;
71 const polyBoundaryMesh& patches = mesh.boundaryMesh();
74 if (patches[patchi].size() && isA<wedgePolyPatch>(patches[patchi]))
76 const wedgePolyPatch& pp =
77 refCast<const wedgePolyPatch>(patches[
patchi]);
79 scalar wedgeAngle =
acos(pp.cosAngle());
83 Info<<
" Wedge " << pp.name() <<
" with angle " 84 <<
radToDeg(wedgeAngle) <<
" degrees" 91 if (oppositePatchi == -1)
95 Info<<
" ***Cannot find opposite wedge for wedge " 101 const wedgePolyPatch& opp =
102 refCast<const wedgePolyPatch>(patches[oppositePatchi]);
105 if (
mag(opp.axis() & pp.axis()) < (1-1
e-3))
109 Info<<
" ***Wedges do not have the same axis." 110 <<
" Encountered " << pp.axis()
111 <<
" on patch " << pp.name()
112 <<
" which differs from " << opp.axis()
113 <<
" on opposite wedge patch" << opp.axis()
124 const face& f = pp[i];
128 label p1 = f.nextLabel(fp);
129 edgesInError.insert(edge(p0, p1), -1);
135 const point& p0 = p[pp.meshPoints()[0]];
136 forAll(pp.meshPoints(), i)
138 const point& pt = p[pp.meshPoints()[i]];
139 scalar d =
mag((pt - p0) & pp.n());
145 Info<<
" ***Wedge patch " << pp.name() <<
" not planar." 146 <<
" Point " << pt <<
" is not in patch plane by " 159 label nEdgesInError = 0;
163 const face& f = fcs[facei];
168 label p1 = f.nextLabel(fp);
172 scalar magD =
mag(d);
174 if (magD > ROOTVSMALL)
179 label nEmptyDirs = 0;
180 label nNonEmptyDirs = 0;
181 for (
direction cmpt=0; cmpt<vector::nComponents; cmpt++)
183 if (
mag(d[cmpt]) > 1
e-6)
185 if (directions[cmpt] == 0)
200 else if (nEmptyDirs == 1)
203 if (nNonEmptyDirs > 0)
205 if (edgesInError.insert(edge(p0, p1), facei))
211 else if (nEmptyDirs > 1)
214 if (edgesInError.insert(edge(p0, p1), facei))
230 Info<<
" ***Number of edges not aligned with or perpendicular to " 231 <<
"non-empty directions: " << nErrorEdges <<
endl;
236 setPtr->resize(2*nEdgesInError);
241 setPtr->insert(iter.key()[0]);
242 setPtr->insert(iter.key()[1]);
253 Info<<
" All edges aligned with or perpendicular to " 254 <<
"non-empty directions." <<
endl;
264 class transformPositionList
271 const coupledPolyPatch& cpp,
272 List<pointField>& pts
278 List<pointField> newPts(pts.size());
281 newPts[facei].setSize(pts[facei].size());
294 if (facePts.size() > index)
296 ptsAtIndex[facei] = facePts[index];
308 cpp.transformPosition(ptsAtIndex);
314 if (facePts.size() > index)
316 facePts[index] = ptsAtIndex[facei];
331 const polyMesh& mesh,
338 const polyBoundaryMesh& patches = mesh.boundaryMesh();
342 List<pointField> nbrPoints(fcs.size() - mesh.nInternalFaces());
347 if (patches[patchi].coupled())
349 const coupledPolyPatch& cpp = refCast<const coupledPolyPatch>
356 label bFacei = cpp.start() + i - mesh.nInternalFaces();
357 const face& f = cpp[i];
358 nbrPoints[bFacei].setSize(f.size());
361 const point& p0 = p[f[fp]];
362 nbrPoints[bFacei][fp] = p0;
372 transformPositionList()
376 label nErrorFaces = 0;
377 scalar avgMismatch = 0;
378 label nCoupledPoints = 0;
382 if (patches[patchi].coupled())
384 const coupledPolyPatch& cpp =
385 refCast<const coupledPolyPatch>(patches[
patchi]);
402 label bFacei = cpp.start() + i - mesh.nInternalFaces();
403 const face& f = cpp[i];
405 if (f.size() != nbrPoints[bFacei].size())
408 <<
"Local face size : " << f.size()
409 <<
" does not equal neighbour face size : " 410 << nbrPoints[bFacei].size()
417 const point& p0 = p[f[fp]];
418 scalar d =
mag(p0 - nbrPoints[bFacei][j]);
420 if (d > smallDist[i])
424 setPtr->insert(cpp.start()+i);
441 reduce(nErrorFaces, sumOp<label>());
442 reduce(avgMismatch, maxOp<scalar>());
443 reduce(nCoupledPoints, sumOp<label>());
445 if (nCoupledPoints > 0)
447 avgMismatch /= nCoupledPoints;
454 Info<<
" **Error in coupled point location: " 456 <<
" faces have their 0th or consecutive vertex not opposite" 457 <<
" their coupled equivalent. Average mismatch " 458 << avgMismatch <<
"." 468 Info<<
" Coupled point location match (average " 469 << avgMismatch <<
") OK." <<
endl;
479 const polyMesh& mesh,
480 const bool allGeometry,
481 const autoPtr<surfaceWriter>& writer
484 label noFailedChecks = 0;
486 Info<<
"\nChecking geometry..." <<
endl;
489 const boundBox& globalBb = mesh.bounds();
491 Info<<
" Overall domain bounding box " 492 << globalBb.min() <<
" " << globalBb.max() <<
endl;
496 scalar minDistSqr =
magSqr(1
e-6 * globalBb.span());
500 Info<<
" Mesh has " << mesh.nGeometricD()
501 <<
" geometric (non-empty/wedge) directions " << validDirs <<
endl;
505 Info<<
" Mesh has " << mesh.nSolutionD()
506 <<
" solution (non-empty) directions " << solDirs <<
endl;
508 if (mesh.nGeometricD() < 3)
510 pointSet nonAlignedPoints(mesh,
"nonAlignedEdges", mesh.nPoints()/100);
516 &&
checkWedges(mesh,
true, validDirs, &nonAlignedPoints)
520 && mesh.checkEdgeAlignment(
true, validDirs, &nonAlignedPoints)
527 nonAlignedPoints.size(),
533 Info<<
" <<Writing " << nNonAligned
534 <<
" points on non-aligned edges to set " 535 << nonAlignedPoints.name() <<
endl;
536 nonAlignedPoints.instance() = mesh.pointsInstance();
537 nonAlignedPoints.
write();
542 if (mesh.checkClosedBoundary(
true)) noFailedChecks++;
545 cellSet
cells(mesh,
"nonClosedCells", mesh.nCells()/100+1);
546 cellSet aspectCells(mesh,
"highAspectRatioCells", mesh.nCells()/100+1);
549 mesh.checkClosedCells
564 Info<<
" <<Writing " << nNonClosed
565 <<
" non closed cells to set " <<
cells.name() <<
endl;
566 cells.instance() = mesh.pointsInstance();
579 Info<<
" <<Writing " << nHighAspect
580 <<
" cells with high aspect ratio to set " 581 << aspectCells.name() <<
endl;
582 aspectCells.instance() = mesh.pointsInstance();
592 faceSet faces(mesh,
"zeroAreaFaces", mesh.nFaces()/100+1);
593 if (mesh.checkFaceAreas(
true, &faces))
601 Info<<
" <<Writing " << nFaces
602 <<
" zero area faces to set " << faces.name() <<
endl;
603 faces.instance() = mesh.pointsInstance();
614 cellSet
cells(mesh,
"zeroVolumeCells", mesh.nCells()/100+1);
615 if (mesh.checkCellVolumes(
true, &
cells))
623 Info<<
" <<Writing " << nCells
624 <<
" zero volume cells to set " <<
cells.name() <<
endl;
625 cells.instance() = mesh.pointsInstance();
636 faceSet faces(mesh,
"nonOrthoFaces", mesh.nFaces()/100+1);
637 if (mesh.checkFaceOrthogonality(
true, &faces))
646 Info<<
" <<Writing " << nFaces
647 <<
" non-orthogonal faces to set " << faces.name() <<
endl;
648 faces.instance() = mesh.pointsInstance();
658 faceSet faces(mesh,
"wrongOrientedFaces", mesh.nFaces()/100 + 1);
659 if (mesh.checkFacePyramids(
true, -SMALL, &faces))
667 Info<<
" <<Writing " << nFaces
668 <<
" faces with incorrect orientation to set " 669 << faces.name() <<
endl;
670 faces.instance() = mesh.pointsInstance();
681 faceSet faces(mesh,
"skewFaces", mesh.nFaces()/100+1);
682 if (mesh.checkFaceSkewness(
true, &faces))
690 Info<<
" <<Writing " << nFaces
691 <<
" skew faces to set " << faces.name() <<
endl;
692 faces.instance() = mesh.pointsInstance();
703 faceSet faces(mesh,
"coupledFaces", mesh.nFaces()/100 + 1);
712 Info<<
" <<Writing " << nFaces
713 <<
" faces with incorrectly matched 0th (or consecutive)" 715 << faces.name() <<
endl;
716 faces.instance() = mesh.pointsInstance();
728 faceSet faces(mesh,
"lowQualityTetFaces", mesh.nFaces()/100+1);
746 Info<<
" <<Writing " << nFaces
747 <<
" faces with low quality or negative volume " 748 <<
"decomposition tets to set " << faces.name() <<
endl;
749 faces.instance() = mesh.pointsInstance();
762 pointSet
points(mesh,
"shortEdges", mesh.nPoints()/1000 + 1);
763 if (mesh.checkEdgeLength(
true, minDistSqr, &
points))
771 Info<<
" <<Writing " << nPoints
772 <<
" points on short edges to set " <<
points.name()
774 points.instance() = mesh.pointsInstance();
781 if (mesh.checkPointNearness(
false, minDistSqr, &
points))
787 if (nPoints > nEdgeClose)
789 pointSet nearPoints(mesh,
"nearPoints",
points);
790 Info<<
" <<Writing " << nPoints
791 <<
" near (closer than " <<
Foam::sqrt(minDistSqr)
792 <<
" apart) points to set " << nearPoints.
name() <<
endl;
793 nearPoints.instance() = mesh.pointsInstance();
801 faceSet faces(mesh,
"concaveFaces", mesh.nFaces()/100 + 1);
802 if (mesh.checkFaceAngles(
true, 10, &faces))
810 Info<<
" <<Writing " << nFaces
811 <<
" faces with concave angles to set " << faces.name()
813 faces.instance() = mesh.pointsInstance();
825 faceSet faces(mesh,
"warpedFaces", mesh.nFaces()/100 + 1);
826 if (mesh.checkFaceFlatness(
true, 0.8, &faces))
834 Info<<
" <<Writing " << nFaces
835 <<
" warped faces to set " << faces.name() <<
endl;
836 faces.instance() = mesh.pointsInstance();
848 cellSet
cells(mesh,
"underdeterminedCells", mesh.nCells()/100);
849 if (mesh.checkCellDeterminant(
true, &
cells))
855 Info<<
" <<Writing " << nCells
856 <<
" under-determined cells to set " <<
cells.name() <<
endl;
857 cells.instance() = mesh.pointsInstance();
868 cellSet
cells(mesh,
"concaveCells", mesh.nCells()/100);
869 if (mesh.checkConcaveCells(
true, &
cells))
875 Info<<
" <<Writing " << nCells
876 <<
" concave cells to set " <<
cells.name() <<
endl;
877 cells.instance() = mesh.pointsInstance();
888 faceSet faces(mesh,
"lowWeightFaces", mesh.nFaces()/100);
889 if (mesh.checkFaceWeight(
true, 0.05, &faces))
895 Info<<
" <<Writing " << nFaces
896 <<
" faces with low interpolation weights to set " 897 << faces.name() <<
endl;
898 faces.instance() = mesh.pointsInstance();
909 faceSet faces(mesh,
"lowVolRatioFaces", mesh.nFaces()/100);
910 if (mesh.checkVolRatio(
true, 0.01, &faces))
916 Info<<
" <<Writing " << nFaces
917 <<
" faces with low volume ratio cells to set " 918 << faces.name() <<
endl;
919 faces.instance() = mesh.pointsInstance();
928 return noFailedChecks;
dimensionedScalar acos(const dimensionedScalar &ds)
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
scalar radToDeg(const scalar rad)
Conversion from radians to degrees.
const double e
Elementary charge.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
static bool checkFaceTets(const polyMesh &mesh, scalar tol=minTetQuality, const bool report=false, labelHashSet *setPtr=NULL)
Check face-decomposition tet volume.
Unit conversion functions.
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
static const Vector< label > one
Vector< scalar > vector
A scalar version of the templated Vector.
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys.
const word & name() const
Return const reference to name.
vectorField pointField
pointField is a vectorField.
const dimensionedScalar e
Elementary charge.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
label checkGeometry(const polyMesh &mesh, const bool allGeometry, const autoPtr< surfaceWriter > &)
bool checkWedges(const polyMesh &, const bool report, const Vector< label > &, labelHashSet *)
Check wedge orientation.
void mergeAndWrite(const polyMesh &mesh, const surfaceWriter &writer, const word &name, const indirectPrimitivePatch setPatch, const fileName &outputDir)
Generate merged surface on master and write. Needs input patch.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
errorManip< error > abort(error &err)
label findOppositeWedge(const polyMesh &, const wedgePolyPatch &)
dimensioned< scalar > magSqr(const dimensioned< Type > &)
bool checkCoupledPoints(const polyMesh &, const bool report, labelHashSet *)
Check 0th vertex on coupled faces.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
vector point
Point is a vector.
dimensioned< scalar > mag(const dimensioned< Type > &)
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
virtual Ostream & write(const token &)=0
Write next token to stream.
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
static const scalar minTetQuality
Minimum tetrahedron quality.