58 const labelList& nei = mesh_.faceNeighbour();
62 vector d = centres[nei[facei]] - centres[own[facei]];
69 result[own[facei]] =
max(cosDDotS, result[own[facei]]);
71 result[nei[facei]] =
max(cosDDotS, result[nei[facei]]);
77 mesh_.boundaryMesh()[
patchi].faceCells();
80 mesh_.boundaryMesh()[
patchi].faceCentres();
83 mesh_.boundaryMesh()[
patchi].faceAreas();
87 vector d = faceCentres[facei] - centres[faceCells[facei]];
94 result[faceCells[facei]] =
max(cosDDotS, result[faceCells[facei]]);
119 const labelList& own = mesh_.faceOwner();
120 const labelList& nei = mesh_.faceNeighbour();
126 (faceCtrs[facei] - cellCtrs[own[facei]]) & areas[facei]
131 (cellCtrs[nei[facei]] - faceCtrs[facei]) & areas[facei]
134 point faceIntersection =
136 + (dOwn/(dOwn+dNei))*(cellCtrs[nei[facei]] - cellCtrs[own[facei]]);
139 mag(faceCtrs[facei] - faceIntersection)
140 /(
mag(cellCtrs[nei[facei]] - cellCtrs[own[facei]]) + vSmall);
142 result[own[facei]] =
max(skewness, result[own[facei]]);
144 result[nei[facei]] =
max(skewness, result[nei[facei]]);
150 mesh_.boundaryMesh()[
patchi].faceCells();
153 mesh_.boundaryMesh()[
patchi].faceCentres();
156 mesh_.boundaryMesh()[
patchi].faceAreas();
158 forAll(faceCentres, facei)
160 vector n = faceAreas[facei]/
mag(faceAreas[facei]);
162 point faceIntersection =
163 cellCtrs[faceCells[facei]]
164 + ((faceCentres[facei] - cellCtrs[faceCells[facei]])&
n)*
n;
167 mag(faceCentres[facei] - faceIntersection)
169 mag(faceCentres[facei] - cellCtrs[faceCells[facei]])
173 result[faceCells[facei]] =
max(skewness, result[faceCells[facei]]);
196 const labelList& own = mesh_.faceOwner();
197 const labelList& nei = mesh_.faceNeighbour();
201 vector d = centres[nei[facei]] - centres[own[facei]];
203 scalar magS =
mag(
s);
208 result[facei] = cosDDotS;
211 label globalFacei = mesh_.nInternalFaces();
216 mesh_.boundaryMesh()[
patchi].faceCells();
219 mesh_.boundaryMesh()[
patchi].faceCentres();
222 mesh_.boundaryMesh()[
patchi].faceAreas();
224 forAll(faceCentres, facei)
226 vector d = faceCentres[facei] - centres[faceCells[facei]];
228 scalar magS =
mag(
s);
233 result[globalFacei++] = cosDDotS;
257 const labelList& own = mesh_.faceOwner();
258 const labelList& nei = mesh_.faceNeighbour();
264 (faceCtrs[facei] - cellCtrs[own[facei]]) & areas[facei]
269 (cellCtrs[nei[facei]] - faceCtrs[facei]) & areas[facei]
272 point faceIntersection =
274 + (dOwn/(dOwn+dNei))*(cellCtrs[nei[facei]] - cellCtrs[own[facei]]);
277 mag(faceCtrs[facei] - faceIntersection)
278 /(
mag(cellCtrs[nei[facei]] - cellCtrs[own[facei]]) + vSmall);
282 label globalFacei = mesh_.nInternalFaces();
287 mesh_.boundaryMesh()[
patchi].faceCells();
290 mesh_.boundaryMesh()[
patchi].faceCentres();
293 mesh_.boundaryMesh()[
patchi].faceAreas();
295 forAll(faceCentres, facei)
297 vector n = faceAreas[facei]/
mag(faceAreas[facei]);
299 point faceIntersection =
300 cellCtrs[faceCells[facei]]
301 + ((faceCentres[facei] - cellCtrs[faceCells[facei]])&
n)*
n;
303 result[globalFacei++] =
304 mag(faceCentres[facei] - faceIntersection)
306 mag(faceCentres[facei] - cellCtrs[faceCells[facei]])
#define forAll(list, i)
Loop across all elements in list.
Pre-declare related SubField type.
tmp< scalarField > nonOrthogonality() const
Return cell non-orthogonality.
tmp< scalarField > skewness() const
Return cell skewness.
tmp< scalarField > faceSkewness() const
Return face skewness.
tmp< scalarField > faceNonOrthogonality() const
Return face non-orthogonality.
cellQuality(const polyMesh &mesh)
Construct from mesh.
Mesh consisting of general polyhedral cells.
A class for managing temporary objects.
T & ref() const
Return non-const reference or generate a fatal error.
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.name(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
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.
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
dimensioned< scalar > mag(const dimensioned< Type > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar acos(const dimensionedScalar &ds)
Unit conversion functions.