Go to the documentation of this file.
47 void Foam::polyMeshGeometry::updateFaceCentresAndAreas
55 for (
const label facei : changedFaces)
64 faceCentres_[facei] = (1.0/3.0)*(
p[
f[0]] +
p[
f[1]] +
p[
f[2]]);
65 faceAreas_[facei] = 0.5*((
p[
f[1]] -
p[
f[0]])^(
p[
f[2]] -
p[
f[0]]));
94 faceCentres_[facei] = (1.0/3.0)*sumAc/(sumA + VSMALL);
95 faceAreas_[facei] = 0.5*sumN;
101 void Foam::polyMeshGeometry::updateCellCentresAndVols
111 UIndirectList<vector>(cellCentres_, changedCells) =
Zero;
112 UIndirectList<scalar>(cellVolumes_, changedCells) =
Zero;
116 for (
const label celli : changedCells)
124 for (
const label facei : cFaces)
126 const point& fc = faceCentres_[facei];
130 cEst /= cFaces.
size();
134 for (
const label facei : cFaces)
137 scalar pyr3Vol = faceAreas_[facei] & (faceCentres_[facei] - cEst);
139 if (own[facei] != celli)
145 cellVolumes_[celli] += pyr3Vol;
148 const vector pCtr = (3.0/4.0)*faceCentres_[facei] + (1.0/4.0)*cEst;
151 cellCentres_[celli] += pyr3Vol*pCtr;
156 if (
mag(cellVolumes_[celli]) > VSMALL)
158 point cc = cellCentres_[celli] / cellVolumes_[celli];
165 cellCentres_[celli] = cc;
169 cellCentres_[celli] = cEst;
174 cellCentres_[celli] = cEst;
177 cellVolumes_[celli] *= (1.0/3.0);
193 for (
const label facei : changedFaces)
195 affectedCells.insert(own[facei]);
199 affectedCells.insert(nei[facei]);
202 return affectedCells.toc();
206 Foam::scalar Foam::polyMeshGeometry::checkNonOrtho
210 const scalar severeNonorthogonalityThreshold,
215 label& severeNonOrth,
220 scalar dDotS = (d &
s)/(
mag(d)*
mag(
s) + VSMALL);
222 if (dDotS < severeNonorthogonalityThreshold)
236 Pout<<
"Severe non-orthogonality for face " << facei
252 <<
"Severe non-orthogonality detected for face "
274 bool Foam::polyMeshGeometry::checkFaceTet
276 const polyMesh&
mesh,
278 const scalar minTetQuality,
299 if (tetQual < minTetQuality)
303 Pout<<
"bool polyMeshGeometry::checkFaceTets("
304 <<
"const bool, const scalar, const pointField&"
305 <<
", const pointField&"
306 <<
", const labelList&, labelHashSet*) : "
308 <<
" has a triangle that points the wrong way." <<
nl
309 <<
"Tet quality: " << tetQual
315 setPtr->insert(facei);
340 faceAreas_ = mesh_.faceAreas();
341 faceCentres_ = mesh_.faceCentres();
342 cellCentres_ = mesh_.cellCentres();
343 cellVolumes_ = mesh_.cellVolumes();
354 updateFaceCentresAndAreas(
p, changedFaces);
356 updateCellCentresAndVols(affectedCells(mesh_, changedFaces), changedFaces);
363 const scalar orthWarn,
380 const scalar severeNonorthogonalityThreshold =
::cos(
degToRad(orthWarn));
392 scalar minDDotS = GREAT;
397 label severeNonOrth = 0;
399 label errorNonOrth = 0;
401 for (
const label facei : checkFaces)
403 const point& ownCc = cellCentres[own[facei]];
407 scalar dDotS = checkNonOrtho
411 severeNonorthogonalityThreshold,
414 cellCentres[nei[facei]] - ownCc,
421 if (dDotS < minDDotS)
435 scalar dDotS = checkNonOrtho
439 severeNonorthogonalityThreshold,
449 if (dDotS < minDDotS)
462 const label face0 = baffle.first();
463 const label face1 = baffle.second();
465 const point& ownCc = cellCentres[own[face0]];
467 scalar dDotS = checkNonOrtho
471 severeNonorthogonalityThreshold,
474 cellCentres[own[face1]] - ownCc,
481 if (dDotS < minDDotS)
499 if (report && minDDotS < severeNonorthogonalityThreshold)
501 Info<<
"Number of non-orthogonality errors: " << errorNonOrth
502 <<
". Number of severely non-orthogonal faces: "
503 << severeNonOrth <<
"." <<
endl;
511 Info<<
"Mesh non-orthogonality Max: "
518 if (errorNonOrth > 0)
523 <<
"Error in non-orthogonality detected" <<
endl;
531 Info<<
"Non-orthogonality check OK.\n" <<
endl;
541 const scalar minPyrVol,
556 label nErrorPyrs = 0;
558 for (
const label facei : checkFaces)
564 cellCentres[own[facei]]
567 if (pyrVol > -minPyrVol)
573 Pout<<
"bool polyMeshGeometry::checkFacePyramids("
574 <<
"const bool, const scalar, const pointField&"
575 <<
", const labelList&, labelHashSet*): "
576 <<
"face " << facei <<
" points the wrong way. " <<
endl
577 <<
"Pyramid volume: " << -pyrVol
578 <<
" Face " <<
f[facei] <<
" area: " <<
f[facei].mag(
p)
579 <<
" Owner cell: " << own[facei] <<
endl
580 <<
"Owner cell vertex labels: "
596 if (pyrVol < minPyrVol)
602 Pout<<
"bool polyMeshGeometry::checkFacePyramids("
603 <<
"const bool, const scalar, const pointField&"
604 <<
", const labelList&, labelHashSet*): "
605 <<
"face " << facei <<
" points the wrong way. " <<
endl
606 <<
"Pyramid volume: " << -pyrVol
607 <<
" Face " <<
f[facei] <<
" area: " <<
f[facei].mag(
p)
608 <<
" Neighbour cell: " << nei[facei] <<
endl
609 <<
"Neighbour cell vertex labels: "
623 const label face0 = baffle.first();
624 const label face1 = baffle.second();
626 const point& ownCc = cellCentres[own[face0]];
635 if (pyrVolOwn > -minPyrVol)
641 Pout<<
"bool polyMeshGeometry::checkFacePyramids("
642 <<
"const bool, const scalar, const pointField&"
643 <<
", const labelList&, labelHashSet*): "
644 <<
"face " << face0 <<
" points the wrong way. " <<
endl
645 <<
"Pyramid volume: " << -pyrVolOwn
646 <<
" Face " <<
f[face0] <<
" area: " <<
f[face0].mag(
p)
647 <<
" Owner cell: " << own[face0] <<
endl
648 <<
"Owner cell vertex labels: "
662 if (pyrVolNbr < minPyrVol)
668 Pout<<
"bool polyMeshGeometry::checkFacePyramids("
669 <<
"const bool, const scalar, const pointField&"
670 <<
", const labelList&, labelHashSet*): "
671 <<
"face " << face0 <<
" points the wrong way. " <<
endl
672 <<
"Pyramid volume: " << -pyrVolNbr
673 <<
" Face " <<
f[face0] <<
" area: " <<
f[face0].mag(
p)
674 <<
" Neighbour cell: " << own[face1] <<
endl
675 <<
"Neighbour cell vertex labels: "
693 <<
"Error in face pyramids: faces pointing the wrong way."
702 Info<<
"Face pyramids OK.\n" <<
endl;
712 const scalar minTetQuality,
738 label nErrorTets = 0;
740 for (
const label facei : checkFaces)
744 bool tetError = checkFaceTet
751 cellCentres[own[facei]],
764 bool tetError = checkFaceTet
772 cellCentres[nei[facei]],
849 const label face0 = baffle.first();
850 const label face1 = baffle.second();
852 bool tetError = checkFaceTet
859 cellCentres[own[face0]],
870 tetError = checkFaceTet
878 cellCentres[own[face1]],
893 cellCentres[own[face1]],
914 <<
"Error in face decomposition: negative tets."
933 const scalar internalSkew,
934 const scalar boundarySkew,
960 for (
const label facei : checkFaces)
972 cellCentres[own[facei]],
973 cellCentres[nei[facei]]
979 if (skewness > internalSkew)
985 Pout<<
"Severe skewness for face " << facei
986 <<
" skewness = " << skewness <<
endl;
994 maxSkew =
max(maxSkew, skewness);
1006 cellCentres[own[facei]],
1013 if (skewness > internalSkew)
1019 Pout<<
"Severe skewness for coupled face " << facei
1020 <<
" skewness = " << skewness <<
endl;
1028 maxSkew =
max(maxSkew, skewness);
1040 cellCentres[own[facei]]
1047 if (skewness > boundarySkew)
1053 Pout<<
"Severe skewness for boundary face " << facei
1054 <<
" skewness = " << skewness <<
endl;
1062 maxSkew =
max(maxSkew, skewness);
1068 const label face0 = baffle.first();
1069 const label face1 = baffle.second();
1071 const point& ownCc = cellCentres[own[face0]];
1072 const point& neiCc = cellCentres[own[face1]];
1089 if (skewness > internalSkew)
1095 Pout<<
"Severe skewness for face " << face0
1096 <<
" skewness = " << skewness <<
endl;
1104 maxSkew =
max(maxSkew, skewness);
1117 <<
" percent.\nThis may impair the quality of the result." <<
nl
1118 << nWarnSkew <<
" highly skew faces detected."
1127 Info<<
"Max skewness = " << 100*maxSkew
1128 <<
" percent. Face skewness OK.\n" <<
endl;
1138 const scalar warnWeight,
1164 scalar minWeight = GREAT;
1166 label nWarnWeight = 0;
1168 for (
const label facei : checkFaces)
1170 const point& fc = faceCentres[facei];
1171 const vector& fa = faceAreas[facei];
1173 scalar dOwn =
mag(fa & (fc-cellCentres[own[facei]]));
1177 scalar dNei =
mag(fa & (cellCentres[nei[facei]]-fc));
1178 scalar weight =
min(dNei,dOwn)/(dNei+dOwn+VSMALL);
1180 if (weight < warnWeight)
1186 Pout<<
"Small weighting factor for face " << facei
1187 <<
" weight = " << weight <<
endl;
1195 minWeight =
min(minWeight, weight);
1201 if (
patches[patchi].coupled())
1204 scalar weight =
min(dNei,dOwn)/(dNei+dOwn+VSMALL);
1206 if (weight < warnWeight)
1212 Pout<<
"Small weighting factor for face " << facei
1213 <<
" weight = " << weight <<
endl;
1221 minWeight =
min(minWeight, weight);
1228 const label face0 = baffle.first();
1229 const label face1 = baffle.second();
1231 const point& ownCc = cellCentres[own[face0]];
1232 const point& fc = faceCentres[face0];
1233 const vector& fa = faceAreas[face0];
1235 scalar dOwn =
mag(fa & (fc-ownCc));
1236 scalar dNei =
mag(fa & (cellCentres[own[face1]]-fc));
1237 scalar weight =
min(dNei,dOwn)/(dNei+dOwn+VSMALL);
1239 if (weight < warnWeight)
1245 Pout<<
"Small weighting factor for face " << face0
1246 <<
" weight = " << weight <<
endl;
1254 minWeight =
min(minWeight, weight);
1260 if (minWeight < warnWeight)
1265 << minWeight <<
'.' <<
nl
1266 << nWarnWeight <<
" faces with small weights detected."
1275 Info<<
"Min weight = " << minWeight
1276 <<
". Weights OK.\n" <<
endl;
1286 const scalar warnRatio,
1310 scalar minRatio = GREAT;
1312 label nWarnRatio = 0;
1314 for (
const label facei : checkFaces)
1316 scalar ownVol =
mag(cellVolumes[own[facei]]);
1318 scalar neiVol = -GREAT;
1322 neiVol =
mag(cellVolumes[nei[facei]]);
1328 if (
patches[patchi].coupled())
1336 scalar ratio =
min(ownVol, neiVol) / (
max(ownVol, neiVol) + VSMALL);
1338 if (ratio < warnRatio)
1344 Pout<<
"Small ratio for face " << facei
1345 <<
" ratio = " << ratio <<
endl;
1353 minRatio =
min(minRatio, ratio);
1359 const label face0 = baffle.first();
1360 const label face1 = baffle.second();
1362 scalar ownVol =
mag(cellVolumes[own[face0]]);
1364 scalar neiVol =
mag(cellVolumes[own[face1]]);
1368 scalar ratio =
min(ownVol, neiVol) / (
max(ownVol, neiVol) + VSMALL);
1370 if (ratio < warnRatio)
1376 Pout<<
"Small ratio for face " << face0
1377 <<
" ratio = " << ratio <<
endl;
1385 minRatio =
min(minRatio, ratio);
1392 if (minRatio < warnRatio)
1397 << minRatio <<
'.' <<
nl
1398 << nWarnRatio <<
" faces with small ratios detected."
1407 Info<<
"Min ratio = " << minRatio
1408 <<
". Ratios OK.\n" <<
endl;
1422 const scalar maxDeg,
1430 if (maxDeg < -SMALL || maxDeg > 180+SMALL)
1433 <<
"maxDeg should be [0..180] but is now " << maxDeg
1441 scalar maxEdgeSin = 0.0;
1445 label errorFacei = -1;
1447 for (
const label facei : checkFaces)
1449 const face&
f = fcs[facei];
1455 scalar magEPrev =
mag(ePrev);
1456 ePrev /= magEPrev + VSMALL;
1461 label fp1 =
f.fcIndex(fp0);
1465 scalar magE10 =
mag(e10);
1466 e10 /= magE10 + VSMALL;
1468 if (magEPrev > SMALL && magE10 > SMALL)
1470 vector edgeNormal = ePrev ^ e10;
1471 scalar magEdgeNormal =
mag(edgeNormal);
1473 if (magEdgeNormal < maxSin)
1480 edgeNormal /= magEdgeNormal;
1482 if ((edgeNormal & faceNormal) < SMALL)
1484 if (facei != errorFacei)
1491 maxEdgeSin =
max(maxEdgeSin, magEdgeNormal);
1511 if (maxEdgeSin > SMALL)
1513 scalar maxConcaveDegr =
1516 Info<<
"There are " << nConcave
1517 <<
" faces with concave angles between consecutive"
1518 <<
" edges. Max concave angle = " << maxConcaveDegr
1519 <<
" degrees.\n" <<
endl;
1523 Info<<
"All angles in faces are convex or less than " << maxDeg
1524 <<
" degrees concave.\n" <<
endl;
1533 << nConcave <<
" face points with severe concave angle (> "
1534 << maxDeg <<
" deg) found.\n"
1550 const scalar minTwist,
1560 if (minTwist < -1-SMALL || minTwist > 1+SMALL)
1563 <<
"minTwist should be [-1..1] but is now " << minTwist
1585 for (
const label facei : checkFaces)
1587 const face&
f = fcs[facei];
1598 cellCentres[nei[facei]]
1599 - cellCentres[own[facei]]
1608 - cellCentres[own[facei]]
1617 - cellCentres[own[facei]]
1623 const point& fc = faceCentres[facei];
1632 p[
f.nextLabel(fpI)],
1637 scalar magTri =
mag(triArea);
1639 if (magTri > VSMALL && ((nf & triArea/magTri) < minTwist))
1661 Info<<
"There are " << nWarped
1662 <<
" faces with cosine of the angle"
1663 <<
" between triangle normal and face normal less than "
1664 << minTwist <<
nl <<
endl;
1668 Info<<
"All faces are flat in that the cosine of the angle"
1669 <<
" between triangle normal and face normal less than "
1670 << minTwist <<
nl <<
endl;
1679 << nWarped <<
" faces with severe warpage "
1680 <<
"(cosine of the angle between triangle normal and "
1681 <<
"face normal < " << minTwist <<
") found.\n"
1696 const scalar minTwist,
1705 if (minTwist < -1-SMALL || minTwist > 1+SMALL)
1708 <<
"minTwist should be [-1..1] but is now " << minTwist
1716 for (
const label facei : checkFaces)
1718 const face&
f = fcs[facei];
1722 const point& fc = faceCentres[facei];
1737 scalar magTri =
mag(prevN);
1739 if (magTri > VSMALL)
1764 scalar magTri =
mag(triN);
1766 if (magTri > VSMALL)
1770 if ((prevN & triN) < minTwist)
1784 else if (minTwist > 0)
1796 while (fp != startFp);
1808 Info<<
"There are " << nWarped
1809 <<
" faces with cosine of the angle"
1810 <<
" between consecutive triangle normals less than "
1811 << minTwist <<
nl <<
endl;
1815 Info<<
"All faces are flat in that the cosine of the angle"
1816 <<
" between consecutive triangle normals is less than "
1817 << minTwist <<
nl <<
endl;
1826 << nWarped <<
" faces with severe warpage "
1827 <<
"(cosine of the angle between consecutive triangle normals"
1828 <<
" < " << minTwist <<
") found.\n"
1842 const scalar minFlatness,
1851 if (minFlatness < -SMALL || minFlatness > 1+SMALL)
1854 <<
"minFlatness should be [0..1] but is now " << minFlatness
1862 for (
const label facei : checkFaces)
1864 const face&
f = fcs[facei];
1868 const point& fc = faceCentres[facei];
1871 scalar sumArea = 0.0;
1883 if (sumArea/
mag(faceAreas[facei]) < minFlatness)
1900 Info<<
"There are " << nWarped
1901 <<
" faces with area of invidual triangles"
1902 <<
" compared to overall area less than "
1903 << minFlatness <<
nl <<
endl;
1907 Info<<
"All faces are flat in that the area of invidual triangles"
1908 <<
" compared to overall area is less than "
1909 << minFlatness <<
nl <<
endl;
1918 << nWarped <<
" non-flat faces "
1919 <<
"(area of invidual triangles"
1920 <<
" compared to overall area"
1921 <<
" < " << minFlatness <<
") found.\n"
1935 const scalar minArea,
1942 label nZeroArea = 0;
1944 for (
const label facei : checkFaces)
1946 if (
mag(faceAreas[facei]) < minArea)
1963 Info<<
"There are " << nZeroArea
1964 <<
" faces with area < " << minArea <<
'.' <<
nl <<
endl;
1968 Info<<
"All faces have area > " << minArea <<
'.' <<
nl <<
endl;
1977 << nZeroArea <<
" faces with area < " << minArea
1992 const scalar warnDet,
2002 scalar minDet = GREAT;
2003 scalar sumDet = 0.0;
2007 for (
const label celli : affectedCells)
2012 scalar magAreaSum = 0;
2014 for (
const label facei : cFaces)
2016 scalar magArea =
mag(faceAreas[facei]);
2018 magAreaSum += magArea;
2019 areaSum += faceAreas[facei]*(faceAreas[facei]/(magArea+VSMALL));
2022 scalar scaledDet =
det(areaSum/(magAreaSum+VSMALL))/0.037037037037037;
2024 minDet =
min(minDet, scaledDet);
2025 sumDet += scaledDet;
2028 if (scaledDet < warnDet)
2047 Info<<
"Cell determinant (1 = uniform cube) : average = "
2048 << sumDet / nSumDet <<
" min = " << minDet <<
endl;
2053 Info<<
"There are " << nWarnDet
2054 <<
" cells with determinant < " << warnDet <<
'.' <<
nl
2059 Info<<
"All faces have determinant > " << warnDet <<
'.' <<
nl
2069 << nWarnDet <<
" cells with determinant < " << warnDet
2084 const scalar orthWarn,
2090 return checkFaceDotProduct
2107 const scalar minPyrVol,
2114 return checkFacePyramids
2131 const scalar minTetQuality,
2138 return checkFaceTets
2182 const scalar warnWeight,
2188 return checkFaceWeights
2206 const scalar warnRatio,
2212 return checkVolRatio
2228 const scalar maxDeg,
2234 return checkFaceAngles
2250 const scalar minTwist,
2256 return checkFaceTwist
2274 const scalar minTwist,
2280 return checkTriangleTwist
2297 const scalar minFlatness,
2303 return checkFaceFlatness
2320 const scalar minArea,
2325 return checkFaceArea
2340 const scalar warnDet,
2346 return checkCellDeterminant
List< label > labelList
A List of labels.
vectorField pointField
pointField is a vectorField.
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.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
static constexpr const zero Zero
Global zero.
label nInternalFaces() const
Number of internal faces.
dimensionedScalar sin(const dimensionedScalar &ds)
label nFaces() const
Number of mesh faces.
static const boundBox invertedBox
A large inverted boundBox: min/max == +/- ROOTVGREAT.
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)
static constexpr direction size()
Return the number of elements in the VectorSpace = Ncmpts.
const cellList & cells() const
Unit conversion functions.
void correct()
Take over properties from mesh.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
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.
prefixOSstream Pout
An Ostream wrapper for parallel output to std::cout.
polyMeshGeometry(const polyMesh &)
Construct from mesh.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Mesh consisting of general polyhedral cells.
tetrahedron< point, const point & > tetPointRef
#define forAll(list, i)
Loop across all elements in list.
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.
static label findSharedBasePoint(const polyMesh &mesh, label fI, const point &nCc, scalar tol, bool report=false)
Find the first suitable base point to use for a minimum.
A triangle primitive used to calculate face normals and swept volumes.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
static label findBasePoint(const polyMesh &mesh, label fI, scalar tol, bool report=false)
Find the base point to use for a minimum triangle.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
static labelList affectedCells(const polyMesh &, const labelList &changedFaces)
Helper function: get affected cells from faces.
static bool checkFaceArea(const bool report, const scalar minArea, const polyMesh &, const vectorField &faceAreas, const labelList &checkFaces, labelHashSet *setPtr)
Small faces.
messageStream Info
Information stream (uses stdout - output is on the master only)
virtual const labelList & faceOwner() const
Return face owner.
#define SeriousErrorInFunction
Report an error message using Foam::SeriousError.
List< cell > cellList
A List of cells.
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
constexpr scalar radToDeg(const scalar rad) noexcept
Conversion from radians to degrees.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
label nBoundaryFaces() const
Number of boundary faces (== nFaces - nInternalFaces)
errorManip< error > abort(error &err)
Vector< scalar > vector
A scalar version of the templated Vector.
constexpr scalar degToRad(const scalar deg) noexcept
Conversion from degrees to radians.
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.
VectorSpace< Form, Cmpt, Ncmpts > normalised(const VectorSpace< Form, Cmpt, Ncmpts > &vs)
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.
virtual const faceList & faces() const
Return raw faces.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
constexpr scalar pi(M_PI)
An ordered pair of two objects of type <T> with first() and second() elements.
List< face > faceList
A List of faces.
bool isInternalFace(const label faceIndex) const
Return true if given face label is internal to the mesh.
dimensionedScalar acos(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
bool insert(const Key &key)
Insert a new entry, not overwriting existing entries.
const polyBoundaryMesh & patches
const dimensionedScalar c
Speed of light in a vacuum.
static const Vector< scalar > zero
A face is a list of labels corresponding to mesh vertices.
dimensionedScalar det(const dimensionedSphericalTensor &dt)
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.
vector point
Point is a vector.
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.
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys and label hasher.
defineTypeNameAndDebug(combustionModel, 0)
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.
pyramid< point, const point &, const face & > pyramidPointFaceRef
#define WarningInFunction
Report a warning using Foam::Warning.
A cell is defined as a list of faces with extra functionality.
dimensionedScalar asin(const dimensionedScalar &ds)
virtual const labelList & faceNeighbour() const
Return face neighbour.
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)
triangle< point, const point & > triPointRef
dimensionedScalar cos(const dimensionedScalar &ds)
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.