Go to the documentation of this file.
43 stabilisedFvGeometryScheme,
70 fCtrs[facei] = (1.0/3.0)*(
p[
f[0]] +
p[
f[1]] +
p[
f[2]]);
71 fAreas[facei] = 0.5*((
p[
f[1]] -
p[
f[0]])^(
p[
f[2]] -
p[
f[0]]));
80 solveVector fCentre =
p[
f[0]];
83 fCentre += solveVector(
p[
f[
pi]]);
90 solveVector sumA =
Zero;
94 const solveVector nextPoint(
p[
f[nextPi]]);
95 const solveVector thisPoint(
p[
f[
pi]]);
98 (nextPoint - thisPoint)^(fCentre - thisPoint);
108 solveScalar sumAn = 0.0;
109 solveVector sumAnc =
Zero;
113 const solveVector nextPoint(
p[
f[nextPi]]);
114 const solveVector thisPoint(
p[
f[
pi]]);
116 const solveVector
c = thisPoint + nextPoint + fCentre;
117 const solveVector a =
118 (nextPoint - thisPoint)^(fCentre - thisPoint);
120 const scalar an = a & sumAHat;
129 if (sumAn > ROOTVSMALL)
131 fCtrs[facei] = (1.0/3.0)*sumAnc/sumAn;
135 fCtrs[facei] = fCentre;
137 fAreas[facei] = 0.5*sumA;
145 Foam::stabilisedFvGeometryScheme::stabilisedFvGeometryScheme
164 Pout<<
"stabilisedFvGeometryScheme::movePoints() : "
165 <<
"recalculating primitiveMesh centres" <<
endl;
201 std::move(faceCentres),
202 std::move(faceAreas),
203 std::move(cellCentres),
204 std::move(cellVolumes)
int debug
Static debugging option.
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
virtual const pointField & points() const
Return raw points.
static constexpr const zero Zero
Global zero (0)
label nFaces() const
Number of mesh faces.
bool hasCellVolumes() const
Default geometry calculation scheme. Slight stabilisation for bad meshes.
Ostream & endl(Ostream &os)
Add newline and flush stream.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
Mesh consisting of general polyhedral cells.
#define forAll(list, i)
Loop across all elements in list.
label nCells() const
Number of mesh cells.
bool hasCellCentres() const
virtual void movePoints()
Do what is necessary if the mesh has moved.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Macros for easy insertion into run-time selection tables.
Mesh data needed to do the Finite Volume discretisation.
bool hasFaceAreas() const
bool hasFaceCentres() const
VectorSpace< Form, Cmpt, Ncmpts > normalised(const VectorSpace< Form, Cmpt, Ncmpts > &vs)
static void makeFaceCentresAndAreas(const polyMesh &mesh, const pointField &p, vectorField &fCtrs, vectorField &fAreas)
Calculate face area and centre weighted using pyramid height.
constexpr scalar pi(M_PI)
Templated 3D Vector derived from VectorSpace adding construction from 3 components,...
void resetGeometry(pointField &&faceCentres, pointField &&faceAreas, pointField &&cellCentres, scalarField &&cellVolumes)
Reset the local geometry.
const dimensionedScalar c
Speed of light in a vacuum.
defineTypeNameAndDebug(combustionModel, 0)
const fvMesh & mesh_
Hold reference to mesh.