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);
339 faceAreas_ = mesh_.faceAreas();
340 faceCentres_ = mesh_.faceCentres();
341 cellCentres_ = mesh_.cellCentres();
342 cellVolumes_ = mesh_.cellVolumes();
353 updateFaceCentresAndAreas(
p, changedFaces);
355 updateCellCentresAndVols(affectedCells(mesh_, changedFaces), changedFaces);
362 const scalar orthWarn,
379 const scalar severeNonorthogonalityThreshold =
::cos(
degToRad(orthWarn));
391 scalar minDDotS = GREAT;
396 label severeNonOrth = 0;
398 label errorNonOrth = 0;
400 for (
const label facei : checkFaces)
402 const point& ownCc = cellCentres[own[facei]];
406 scalar dDotS = checkNonOrtho
410 severeNonorthogonalityThreshold,
413 cellCentres[nei[facei]] - ownCc,
420 if (dDotS < minDDotS)
434 scalar dDotS = checkNonOrtho
438 severeNonorthogonalityThreshold,
448 if (dDotS < minDDotS)
461 const label face0 = baffle.first();
462 const label face1 = baffle.second();
464 const point& ownCc = cellCentres[own[face0]];
466 scalar dDotS = checkNonOrtho
470 severeNonorthogonalityThreshold,
473 cellCentres[own[face1]] - ownCc,
480 if (dDotS < minDDotS)
498 if (report && minDDotS < severeNonorthogonalityThreshold)
500 Info<<
"Number of non-orthogonality errors: " << errorNonOrth
501 <<
". Number of severely non-orthogonal faces: "
502 << severeNonOrth <<
"." <<
endl;
510 Info<<
"Mesh non-orthogonality Max: "
517 if (errorNonOrth > 0)
522 <<
"Error in non-orthogonality detected" <<
endl;
530 Info<<
"Non-orthogonality check OK.\n" <<
endl;
540 const scalar minPyrVol,
555 label nErrorPyrs = 0;
557 for (
const label facei : checkFaces)
563 cellCentres[own[facei]]
566 if (pyrVol > -minPyrVol)
572 Pout<<
"bool polyMeshGeometry::checkFacePyramids("
573 <<
"const bool, const scalar, const pointField&"
574 <<
", const labelList&, labelHashSet*): "
575 <<
"face " << facei <<
" points the wrong way. " <<
endl
576 <<
"Pyramid volume: " << -pyrVol
577 <<
" Face " <<
f[facei] <<
" area: " <<
f[facei].mag(
p)
578 <<
" Owner cell: " << own[facei] <<
endl
579 <<
"Owner cell vertex labels: "
595 if (pyrVol < minPyrVol)
601 Pout<<
"bool polyMeshGeometry::checkFacePyramids("
602 <<
"const bool, const scalar, const pointField&"
603 <<
", const labelList&, labelHashSet*): "
604 <<
"face " << facei <<
" points the wrong way. " <<
endl
605 <<
"Pyramid volume: " << -pyrVol
606 <<
" Face " <<
f[facei] <<
" area: " <<
f[facei].mag(
p)
607 <<
" Neighbour cell: " << nei[facei] <<
endl
608 <<
"Neighbour cell vertex labels: "
622 const label face0 = baffle.first();
623 const label face1 = baffle.second();
625 const point& ownCc = cellCentres[own[face0]];
634 if (pyrVolOwn > -minPyrVol)
640 Pout<<
"bool polyMeshGeometry::checkFacePyramids("
641 <<
"const bool, const scalar, const pointField&"
642 <<
", const labelList&, labelHashSet*): "
643 <<
"face " << face0 <<
" points the wrong way. " <<
endl
644 <<
"Pyramid volume: " << -pyrVolOwn
645 <<
" Face " <<
f[face0] <<
" area: " <<
f[face0].mag(
p)
646 <<
" Owner cell: " << own[face0] <<
endl
647 <<
"Owner cell vertex labels: "
661 if (pyrVolNbr < minPyrVol)
667 Pout<<
"bool polyMeshGeometry::checkFacePyramids("
668 <<
"const bool, const scalar, const pointField&"
669 <<
", const labelList&, labelHashSet*): "
670 <<
"face " << face0 <<
" points the wrong way. " <<
endl
671 <<
"Pyramid volume: " << -pyrVolNbr
672 <<
" Face " <<
f[face0] <<
" area: " <<
f[face0].mag(
p)
673 <<
" Neighbour cell: " << own[face1] <<
endl
674 <<
"Neighbour cell vertex labels: "
692 <<
"Error in face pyramids: faces pointing the wrong way."
701 Info<<
"Face pyramids OK.\n" <<
endl;
711 const scalar minTetQuality,
737 label nErrorTets = 0;
739 for (
const label facei : checkFaces)
743 bool tetError = checkFaceTet
750 cellCentres[own[facei]],
763 bool tetError = checkFaceTet
771 cellCentres[nei[facei]],
848 const label face0 = baffle.first();
849 const label face1 = baffle.second();
851 bool tetError = checkFaceTet
858 cellCentres[own[face0]],
869 tetError = checkFaceTet
877 cellCentres[own[face1]],
892 cellCentres[own[face1]],
913 <<
"Error in face decomposition: negative tets."
932 const scalar internalSkew,
933 const scalar boundarySkew,
959 for (
const label facei : checkFaces)
971 cellCentres[own[facei]],
972 cellCentres[nei[facei]]
978 if (skewness > internalSkew)
984 Pout<<
"Severe skewness for face " << facei
985 <<
" skewness = " << skewness <<
endl;
993 maxSkew =
max(maxSkew, skewness);
1005 cellCentres[own[facei]],
1012 if (skewness > internalSkew)
1018 Pout<<
"Severe skewness for coupled face " << facei
1019 <<
" skewness = " << skewness <<
endl;
1027 maxSkew =
max(maxSkew, skewness);
1039 cellCentres[own[facei]]
1046 if (skewness > boundarySkew)
1052 Pout<<
"Severe skewness for boundary face " << facei
1053 <<
" skewness = " << skewness <<
endl;
1061 maxSkew =
max(maxSkew, skewness);
1067 const label face0 = baffle.first();
1068 const label face1 = baffle.second();
1070 const point& ownCc = cellCentres[own[face0]];
1071 const point& neiCc = cellCentres[own[face1]];
1088 if (skewness > internalSkew)
1094 Pout<<
"Severe skewness for face " << face0
1095 <<
" skewness = " << skewness <<
endl;
1103 maxSkew =
max(maxSkew, skewness);
1116 <<
" percent.\nThis may impair the quality of the result." <<
nl
1117 << nWarnSkew <<
" highly skew faces detected."
1126 Info<<
"Max skewness = " << 100*maxSkew
1127 <<
" percent. Face skewness OK.\n" <<
endl;
1137 const scalar warnWeight,
1163 scalar minWeight = GREAT;
1165 label nWarnWeight = 0;
1167 for (
const label facei : checkFaces)
1169 const point& fc = faceCentres[facei];
1170 const vector& fa = faceAreas[facei];
1172 scalar dOwn =
mag(fa & (fc-cellCentres[own[facei]]));
1176 scalar dNei =
mag(fa & (cellCentres[nei[facei]]-fc));
1177 scalar weight =
min(dNei,dOwn)/(dNei+dOwn+VSMALL);
1179 if (weight < warnWeight)
1185 Pout<<
"Small weighting factor for face " << facei
1186 <<
" weight = " << weight <<
endl;
1194 minWeight =
min(minWeight, weight);
1203 scalar weight =
min(dNei,dOwn)/(dNei+dOwn+VSMALL);
1205 if (weight < warnWeight)
1211 Pout<<
"Small weighting factor for face " << facei
1212 <<
" weight = " << weight <<
endl;
1220 minWeight =
min(minWeight, weight);
1227 const label face0 = baffle.first();
1228 const label face1 = baffle.second();
1230 const point& ownCc = cellCentres[own[face0]];
1231 const point& fc = faceCentres[face0];
1232 const vector& fa = faceAreas[face0];
1234 scalar dOwn =
mag(fa & (fc-ownCc));
1235 scalar dNei =
mag(fa & (cellCentres[own[face1]]-fc));
1236 scalar weight =
min(dNei,dOwn)/(dNei+dOwn+VSMALL);
1238 if (weight < warnWeight)
1244 Pout<<
"Small weighting factor for face " << face0
1245 <<
" weight = " << weight <<
endl;
1253 minWeight =
min(minWeight, weight);
1259 if (minWeight < warnWeight)
1264 << minWeight <<
'.' <<
nl
1265 << nWarnWeight <<
" faces with small weights detected."
1274 Info<<
"Min weight = " << minWeight
1275 <<
". Weights OK.\n" <<
endl;
1285 const scalar warnRatio,
1309 scalar minRatio = GREAT;
1311 label nWarnRatio = 0;
1313 for (
const label facei : checkFaces)
1315 scalar ownVol =
mag(cellVolumes[own[facei]]);
1317 scalar neiVol = -GREAT;
1321 neiVol =
mag(cellVolumes[nei[facei]]);
1335 scalar ratio =
min(ownVol, neiVol) / (
max(ownVol, neiVol) + VSMALL);
1337 if (ratio < warnRatio)
1343 Pout<<
"Small ratio for face " << facei
1344 <<
" ratio = " << ratio <<
endl;
1352 minRatio =
min(minRatio, ratio);
1358 const label face0 = baffle.first();
1359 const label face1 = baffle.second();
1361 scalar ownVol =
mag(cellVolumes[own[face0]]);
1363 scalar neiVol =
mag(cellVolumes[own[face1]]);
1367 scalar ratio =
min(ownVol, neiVol) / (
max(ownVol, neiVol) + VSMALL);
1369 if (ratio < warnRatio)
1375 Pout<<
"Small ratio for face " << face0
1376 <<
" ratio = " << ratio <<
endl;
1384 minRatio =
min(minRatio, ratio);
1391 if (minRatio < warnRatio)
1396 << minRatio <<
'.' <<
nl
1397 << nWarnRatio <<
" faces with small ratios detected."
1406 Info<<
"Min ratio = " << minRatio
1407 <<
". Ratios OK.\n" <<
endl;
1421 const scalar maxDeg,
1429 if (maxDeg < -SMALL || maxDeg > 180+SMALL)
1432 <<
"maxDeg should be [0..180] but is now " << maxDeg
1440 scalar maxEdgeSin = 0.0;
1444 label errorFacei = -1;
1446 for (
const label facei : checkFaces)
1448 const face&
f = fcs[facei];
1454 scalar magEPrev =
mag(ePrev);
1455 ePrev /= magEPrev + VSMALL;
1460 label fp1 =
f.fcIndex(fp0);
1464 scalar magE10 =
mag(e10);
1465 e10 /= magE10 + VSMALL;
1467 if (magEPrev > SMALL && magE10 > SMALL)
1469 vector edgeNormal = ePrev ^ e10;
1470 scalar magEdgeNormal =
mag(edgeNormal);
1472 if (magEdgeNormal < maxSin)
1479 edgeNormal /= magEdgeNormal;
1481 if ((edgeNormal & faceNormal) < SMALL)
1483 if (facei != errorFacei)
1490 maxEdgeSin =
max(maxEdgeSin, magEdgeNormal);
1510 if (maxEdgeSin > SMALL)
1512 scalar maxConcaveDegr =
1515 Info<<
"There are " << nConcave
1516 <<
" faces with concave angles between consecutive"
1517 <<
" edges. Max concave angle = " << maxConcaveDegr
1518 <<
" degrees.\n" <<
endl;
1522 Info<<
"All angles in faces are convex or less than " << maxDeg
1523 <<
" degrees concave.\n" <<
endl;
1532 << nConcave <<
" face points with severe concave angle (> "
1533 << maxDeg <<
" deg) found.\n"
1549 const scalar minTwist,
1559 if (minTwist < -1-SMALL || minTwist > 1+SMALL)
1562 <<
"minTwist should be [-1..1] but is now " << minTwist
1584 for (
const label facei : checkFaces)
1586 const face&
f = fcs[facei];
1597 cellCentres[nei[facei]]
1598 - cellCentres[own[facei]]
1607 - cellCentres[own[facei]]
1616 - cellCentres[own[facei]]
1622 const point& fc = faceCentres[facei];
1631 p[
f.nextLabel(fpI)],
1636 scalar magTri =
mag(triArea);
1638 if (magTri > VSMALL && ((nf & triArea/magTri) < minTwist))
1660 Info<<
"There are " << nWarped
1661 <<
" faces with cosine of the angle"
1662 <<
" between triangle normal and face normal less than "
1663 << minTwist <<
nl <<
endl;
1667 Info<<
"All faces are flat in that the cosine of the angle"
1668 <<
" between triangle normal and face normal less than "
1669 << minTwist <<
nl <<
endl;
1678 << nWarped <<
" faces with severe warpage "
1679 <<
"(cosine of the angle between triangle normal and "
1680 <<
"face normal < " << minTwist <<
") found.\n"
1695 const scalar minTwist,
1704 if (minTwist < -1-SMALL || minTwist > 1+SMALL)
1707 <<
"minTwist should be [-1..1] but is now " << minTwist
1715 for (
const label facei : checkFaces)
1717 const face&
f = fcs[facei];
1721 const point& fc = faceCentres[facei];
1736 scalar magTri =
mag(prevN);
1738 if (magTri > VSMALL)
1763 scalar magTri =
mag(triN);
1765 if (magTri > VSMALL)
1769 if ((prevN & triN) < minTwist)
1783 else if (minTwist > 0)
1795 while (fp != startFp);
1807 Info<<
"There are " << nWarped
1808 <<
" faces with cosine of the angle"
1809 <<
" between consecutive triangle normals less than "
1810 << minTwist <<
nl <<
endl;
1814 Info<<
"All faces are flat in that the cosine of the angle"
1815 <<
" between consecutive triangle normals is less than "
1816 << minTwist <<
nl <<
endl;
1825 << nWarped <<
" faces with severe warpage "
1826 <<
"(cosine of the angle between consecutive triangle normals"
1827 <<
" < " << minTwist <<
") found.\n"
1841 const scalar minFlatness,
1850 if (minFlatness < -SMALL || minFlatness > 1+SMALL)
1853 <<
"minFlatness should be [0..1] but is now " << minFlatness
1861 for (
const label facei : checkFaces)
1863 const face&
f = fcs[facei];
1867 const point& fc = faceCentres[facei];
1870 scalar sumArea = 0.0;
1882 if (sumArea/
mag(faceAreas[facei]) < minFlatness)
1899 Info<<
"There are " << nWarped
1900 <<
" faces with area of individual triangles"
1901 <<
" compared to overall area less than "
1902 << minFlatness <<
nl <<
endl;
1906 Info<<
"All faces are flat in that the area of individual triangles"
1907 <<
" compared to overall area is less than "
1908 << minFlatness <<
nl <<
endl;
1917 << nWarped <<
" non-flat faces "
1918 <<
"(area of individual triangles"
1919 <<
" compared to overall area"
1920 <<
" < " << minFlatness <<
") found.\n"
1934 const scalar minArea,
1941 label nZeroArea = 0;
1943 for (
const label facei : checkFaces)
1945 if (
mag(faceAreas[facei]) < minArea)
1962 Info<<
"There are " << nZeroArea
1963 <<
" faces with area < " << minArea <<
'.' <<
nl <<
endl;
1967 Info<<
"All faces have area > " << minArea <<
'.' <<
nl <<
endl;
1976 << nZeroArea <<
" faces with area < " << minArea
1991 const scalar warnDet,
2001 scalar minDet = GREAT;
2002 scalar sumDet = 0.0;
2006 for (
const label celli : affectedCells)
2011 scalar magAreaSum = 0;
2013 for (
const label facei : cFaces)
2015 scalar magArea =
mag(faceAreas[facei]);
2017 magAreaSum += magArea;
2018 areaSum += faceAreas[facei]*(faceAreas[facei]/(magArea+VSMALL));
2021 scalar scaledDet =
det(areaSum/(magAreaSum+VSMALL))/0.037037037037037;
2023 minDet =
min(minDet, scaledDet);
2024 sumDet += scaledDet;
2027 if (scaledDet < warnDet)
2046 Info<<
"Cell determinant (1 = uniform cube) : average = "
2047 << sumDet / nSumDet <<
" min = " << minDet <<
endl;
2052 Info<<
"There are " << nWarnDet
2053 <<
" cells with determinant < " << warnDet <<
'.' <<
nl
2058 Info<<
"All faces have determinant > " << warnDet <<
'.' <<
nl
2068 << nWarnDet <<
" cells with determinant < " << warnDet
2083 const scalar orthWarn,
2089 return checkFaceDotProduct
2106 const scalar minPyrVol,
2113 return checkFacePyramids
2130 const scalar minTetQuality,
2137 return checkFaceTets
2181 const scalar warnWeight,
2187 return checkFaceWeights
2205 const scalar warnRatio,
2211 return checkVolRatio
2227 const scalar maxDeg,
2233 return checkFaceAngles
2249 const scalar minTwist,
2255 return checkFaceTwist
2273 const scalar minTwist,
2279 return checkTriangleTwist
2296 const scalar minFlatness,
2302 return checkFaceFlatness
2319 const scalar minArea,
2324 return checkFaceArea
2339 const scalar warnDet,
2345 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 (0)
dimensionedScalar sin(const dimensionedScalar &ds)
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)
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
OSstream wrapped stdout (std::cout) with parallel prefix.
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
A tetrahedron using referred points.
#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.
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 (stdout output on master, null elsewhere)
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 nBoundaryFaces() const noexcept
Number of boundary faces (== nFaces - nInternalFaces)
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.
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.
bool isInternalFace(const label faceIndex) const noexcept
Return true if given face label is internal to the mesh.
virtual const faceList & faces() const
Return raw faces.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
constexpr scalar pi(M_PI)
List< face > faceList
A List of faces.
label nInternalFaces() const noexcept
Number of internal faces.
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
label nFaces() const noexcept
Number of mesh faces.
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.
bool coupled(solutionDict.getOrDefault("coupledEnergyField", false))
static constexpr direction size() noexcept
The number of elements in the VectorSpace = Ncmpts.
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 of labels, uses 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
A pyramid using referred points and faces.
#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
A triangle using referred points.
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.