cutFaceAdvect Class Reference

Calculates the face fluxes. More...

Inheritance diagram for cutFaceAdvect:
[legend]
Collaboration diagram for cutFaceAdvect:
[legend]

Public Member Functions

 cutFaceAdvect (const fvMesh &mesh, const volScalarField &alpha1)
 Construct from fvMesh and a scalarField. More...
 
label calcSubFace (const label faceI, const vector &normal, const vector &base)
 Calculates cut centre and cut area (plicReconstruction) More...
 
label calcSubFace (const face &f, const pointField &points, const scalarField &val, const scalar cutValue)
 Calculates cut centre and cut area (iso reconstruction) More...
 
label calcSubFace (const label faceI, const label sign, const scalar cutValue)
 Calculates cut centre and cut area (iso reconstruction) More...
 
scalar timeIntegratedFaceFlux (const label faceI, const vector &x0, const vector &n0, const scalar Un0, const scalar dt, const scalar phi, const scalar magSf)
 Calculate time integrated flux for a face. More...
 
scalar timeIntegratedFaceFlux (const label faceI, const scalarField &times, const scalar Un0, const scalar dt, const scalar phi, const scalar magSf)
 Calculate time integrated flux for a face. More...
 
scalar timeIntegratedArea (const label faceI, const scalar dt, const scalar magSf, const scalar Un0)
 Calculate time integrated area for a face. More...
 
scalar timeIntegratedArea (const pointField &fPts, const scalarField &pTimes, const scalar dt, const scalar magSf, const scalar Un0)
 Calculate time integrated area for a face. More...
 
void quadAreaCoeffs (const DynamicPointList &pf0, const DynamicPointList &pf1, scalar &quadArea, scalar &intQuadArea) const
 For face with vertices p calculate its area and integrated area. More...
 
void cutPoints (const label faceI, const scalar f0, DynamicList< point > &cutPoints)
 Get cutPoints of face. More...
 
void cutPoints (const pointField &pts, const scalarField &f, const scalar f0, DynamicList< point > &cutPoints)
 Get cutPoints of face. More...
 
const pointsubFaceCentre () const noexcept
 Returns centre of cutted face. More...
 
const vectorsubFaceArea () const noexcept
 Returns area vector of cutted face. More...
 
const DynamicList< point > & subFacePoints () const noexcept
 Returns the cut edge of the cutted face. More...
 
const DynamicList< point > & surfacePoints () const noexcept
 Returns point of the face in sorted of cutted face. More...
 
void clearStorage ()
 Resets internal variables. More...
 
- Public Member Functions inherited from cutFace
 cutFace (const fvMesh &mesh)
 Construct from fvMesh. More...
 

Additional Inherited Members

- Static Public Attributes inherited from cutFace
static int debug = 0
 
- Protected Member Functions inherited from cutFace
void calcSubFace (const label faceI, const scalarList &pointStatus, label firstFullySubmergedPoint, DynamicList< point > &subFacePoints, DynamicList< point > &surfacePoints, label &faceStatus, vector &subFaceCentre, vector &subFaceArea)
 
void calcSubFace (const label faceI, const scalarList &pointStatus, const scalarList &weights, label firstFullySubmergedPoint, DynamicList< point > &subFacePoints, DynamicList< point > &surfacePoints, label &faceStatus, vector &subFaceCentre, vector &subFaceArea)
 
void calcSubFaceCentreAndArea (DynamicList< point > &subFacePoints, vector &subFaceCentre, vector &subFaceArea)
 Calculates centre and normal of the face. More...
 
void calcSubFace (const face &f, const pointField &points, const scalarList &pointStatus, label firstFullySubmergedPoint, DynamicList< point > &subFacePoints, DynamicList< point > &surfacePoints, label &faceStatus, vector &subFaceCentre, vector &subFaceArea)
 

Detailed Description

Calculates the face fluxes.

Reference:

    Roenby, J., Bredmose, H. and Jasak, H. (2016).
    A computational method for sharp interface advection
    Royal Society Open Science, 3
    doi 10.1098/rsos.160405

Original code supplied by Johan Roenby, DHI (2016)

Source files

Definition at line 66 of file cutFaceAdvect.H.

Constructor & Destructor Documentation

◆ cutFaceAdvect()

cutFaceAdvect ( const fvMesh mesh,
const volScalarField alpha1 
)

Construct from fvMesh and a scalarField.

Length of scalarField should equal number of mesh points

Definition at line 35 of file cutFaceAdvect.C.

References cutFaceAdvect::clearStorage().

Here is the call graph for this function:

Member Function Documentation

◆ calcSubFace() [1/3]

Foam::label calcSubFace ( const label  faceI,
const vector normal,
const vector base 
)

Calculates cut centre and cut area (plicReconstruction)

Returns
face status

Definition at line 59 of file cutFaceAdvect.C.

References cutFace::calcSubFace(), f(), forAll, Foam::mag(), UList< T >::size(), and Foam::Zero.

Here is the call graph for this function:

◆ calcSubFace() [2/3]

Foam::label calcSubFace ( const face f,
const pointField points,
const scalarField val,
const scalar  cutValue 
)

Calculates cut centre and cut area (iso reconstruction)

Returns
face status

Definition at line 123 of file cutFaceAdvect.C.

References cutFace::calcSubFace(), f(), forAll, Foam::mag(), points, UList< T >::size(), and Foam::Zero.

Here is the call graph for this function:

◆ calcSubFace() [3/3]

Foam::label calcSubFace ( const label  faceI,
const label  sign,
const scalar  cutValue 
)

Calculates cut centre and cut area (iso reconstruction)

Returns
face status

Definition at line 576 of file cutFaceAdvect.C.

References cutFace::calcSubFace(), f(), forAll, Foam::mag(), Foam::sign(), UList< T >::size(), and Foam::Zero.

Here is the call graph for this function:

◆ timeIntegratedFaceFlux() [1/2]

Foam::scalar timeIntegratedFaceFlux ( const label  faceI,
const vector x0,
const vector n0,
const scalar  Un0,
const scalar  dt,
const scalar  phi,
const scalar  magSf 
)

Calculate time integrated flux for a face.

Definition at line 187 of file cutFaceAdvect.C.

References e, Foam::endl(), f(), forAll, Foam::mag(), nPoints, phi, Foam::sign(), UList< T >::size(), and WarningInFunction.

Here is the call graph for this function:

◆ timeIntegratedFaceFlux() [2/2]

Foam::scalar timeIntegratedFaceFlux ( const label  faceI,
const scalarField times,
const scalar  Un0,
const scalar  dt,
const scalar  phi,
const scalar  magSf 
)

Calculate time integrated flux for a face.

Definition at line 334 of file cutFaceAdvect.C.

References forAll, nPoints, phi, Foam::sign(), and UList< T >::size().

Here is the call graph for this function:

◆ timeIntegratedArea() [1/2]

Foam::scalar timeIntegratedArea ( const label  faceI,
const scalar  dt,
const scalar  magSf,
const scalar  Un0 
)

Calculate time integrated area for a face.

Definition at line 639 of file cutFaceAdvect.C.

References alpha, DynamicList< T, SizeMin >::append(), beta(), e, UList< T >::first(), UList< T >::last(), Foam::mag(), Foam::max(), Foam::pos0(), Foam::sign(), and Foam::sortedOrder().

Here is the call graph for this function:

◆ timeIntegratedArea() [2/2]

Foam::scalar timeIntegratedArea ( const pointField fPts,
const scalarField pTimes,
const scalar  dt,
const scalar  magSf,
const scalar  Un0 
)

Calculate time integrated area for a face.

Definition at line 381 of file cutFaceAdvect.C.

References alpha, DynamicList< T, SizeMin >::append(), beta(), e, UList< T >::first(), Foam::identity(), UList< T >::last(), Foam::mag(), Foam::max(), Foam::pos0(), Foam::sign(), UList< T >::size(), and Foam::sortedOrder().

Here is the call graph for this function:

◆ quadAreaCoeffs()

void quadAreaCoeffs ( const DynamicPointList pf0,
const DynamicPointList pf1,
scalar &  quadArea,
scalar &  intQuadArea 
) const

For face with vertices p calculate its area and integrated area.

between isocutting lines with isovalues f0 and f1 given vertex values f

Definition at line 779 of file cutFaceAdvect.C.

References A, alpha, B, beta(), D, Foam::endl(), Foam::mag(), UList< T >::size(), WarningInFunction, and Foam::Zero.

Here is the call graph for this function:

◆ cutPoints() [1/2]

void cutPoints ( const label  faceI,
const scalar  f0,
DynamicList< point > &  cutPoints 
)

Get cutPoints of face.

Definition at line 881 of file cutFaceAdvect.C.

References DynamicList< T, SizeMin >::append(), Foam::endl(), f(), forAll, Foam::mag(), nPoints, s(), UList< T >::size(), and WarningInFunction.

Here is the call graph for this function:

◆ cutPoints() [2/2]

void cutPoints ( const pointField pts,
const scalarField f,
const scalar  f0,
DynamicList< point > &  cutPoints 
)

Get cutPoints of face.

Definition at line 936 of file cutFaceAdvect.C.

References DynamicList< T, SizeMin >::append(), Foam::endl(), f(), forAll, Foam::mag(), nPoints, s(), UList< T >::size(), and WarningInFunction.

Here is the call graph for this function:

◆ subFaceCentre()

const point & subFaceCentre ( ) const
inlinenoexcept

Returns centre of cutted face.

Definition at line 234 of file cutFaceAdvect.H.

◆ subFaceArea()

const vector & subFaceArea ( ) const
inlinenoexcept

Returns area vector of cutted face.

Definition at line 240 of file cutFaceAdvect.H.

◆ subFacePoints()

const DynamicList< point > & subFacePoints ( ) const
inlinenoexcept

Returns the cut edge of the cutted face.

Definition at line 246 of file cutFaceAdvect.H.

◆ surfacePoints()

const DynamicList< point > & surfacePoints ( ) const
inlinenoexcept

Returns point of the face in sorted of cutted face.

Definition at line 252 of file cutFaceAdvect.H.

◆ clearStorage()

void clearStorage ( )

Resets internal variables.

Definition at line 986 of file cutFaceAdvect.C.

References Foam::Zero.

Referenced by cutFaceAdvect::cutFaceAdvect().

Here is the caller graph for this function:

The documentation for this class was generated from the following files: