Go to the documentation of this file.
50 vector& subCellCentre, scalar& subCellVolume
63 forAll(cutFaceCentres, facei)
67 mag(cutFaceAreas[facei] & (cutFaceCentres[facei] - cEst)), VSMALL);
70 vector pc = 0.75 * cutFaceCentres[facei] + 0.25 * cEst;
73 subCellCentre += pyr3Vol * pc;
76 subCellVolume += pyr3Vol;
79 subCellCentre /= subCellVolume;
87 const vector& subCellCentre,
97 for (
const point&
p : edgePoints)
105 fCentre /= nEdgePoints;
116 const label
nPoints = edgePoints.size();
119 const point& nextPoint = edgePoints[
pi + 1];
121 vector c = edgePoints[
pi] + nextPoint + fCentre;
123 (nextPoint - edgePoints[
pi]) ^ (fCentre - edgePoints[
pi]);
135 if (sumA < ROOTVSMALL)
137 faceCentre = fCentre;
142 faceCentre = (1.0/3.0)*sumAc/sumA;
147 if ((faceArea & (faceCentre - subCellCentre)) >= 0)
163 vector xhat = faceEdges[0][0] - faceCentre;
164 xhat = (xhat - (xhat & zhat)*zhat);
173 for (
const point&
p : edgePoints)
176 unsortedFacePointAngles.append
180 ((
p - faceCentre) & yhat),
181 ((
p - faceCentre) & xhat)
189 facePoints.
append(unsortedFacePoints[order[0]]);
190 for (label
pi = 1;
pi < order.size(); ++
pi)
196 unsortedFacePointAngles[order[
pi]]
197 - unsortedFacePointAngles[order[
pi - 1]]
200 facePoints.
append(unsortedFacePoints[order[
pi]]);
static constexpr const zero Zero
Global zero (0)
dimensionedScalar atan2(const dimensionedScalar &x, const dimensionedScalar &y)
static void calcCellData(const DynamicList< point > &cutFaceCentres, const DynamicList< vector > &cutFaceAreas, vector &subCellCentre, scalar &subCellVolume)
Calculates volume and centre of the cutted cell.
Vector< Cmpt > & normalise()
Normalise the vector by its magnitude.
dimensionedScalar sign(const dimensionedScalar &ds)
#define forAll(list, i)
Loop across all elements in list.
cutCell(const fvMesh &unused)
Construct from fvMesh.
DynamicList< T, SizeMin > & append(const T &val)
Append an element to the end of this list.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
static void calcIsoFacePointsFromEdges(const vector &faceArea, const vector &faceCentre, const DynamicList< DynamicList< point >> &faceEdges, DynamicList< point > &facePoints)
Calculates the point of the cutting face.
Mesh data needed to do the Finite Volume discretisation.
VectorSpace< Form, Cmpt, Ncmpts > normalised(const VectorSpace< Form, Cmpt, Ncmpts > &vs)
static void calcGeomDataCutFace(const DynamicList< DynamicList< point >> &faceEdges, const vector &subCellCentre, vector &faceArea, vector &faceCentre)
Calculates area and centre of the cutting face.
constexpr scalar pi(M_PI)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const dimensionedScalar e
Elementary charge.
const dimensionedScalar c
Speed of light in a vacuum.
labelList sortedOrder(const UList< T > &input)
Return the (stable) sort order for the list.
dimensioned< Type > average(const DimensionedField< Type, GeoMesh > &df)