MPLIC Class Reference

Multicut Piecewise-Linear Interface Calculation (MPLIC) corrected scheme is a surface interpolation scheme for flux calculation in advection of a bounded variable, e.g. phase fraction and for interface capturing in the volume of fluid (VoF) method. More...

Inheritance diagram for MPLIC:
Collaboration diagram for MPLIC:

Public Member Functions

 TypeName ("MPLIC")
 Runtime type information. More...
 
 MPLIC (const fvMesh &mesh, const surfaceScalarField &faceFlux, Istream &is)
 Construct from faceFlux and Istream. More...
 
virtual tmp< surfaceScalarFieldweights (const VolField< scalar > &) const
 Return the interpolation weighting factors. More...
 
virtual tmp< surfaceScalarFieldinterpolate (const VolField< scalar > &vf) const
 Return the face-interpolate of the given cell field. More...
 
void operator= (const MPLIC &)=delete
 Disallow default bitwise assignment. More...
 
- Public Member Functions inherited from surfaceInterpolationScheme< scalar >
 TypeName ("surfaceInterpolationScheme")
 Runtime type information. More...
 
 declareRunTimeSelectionTable (tmp, surfaceInterpolationScheme, Mesh,(const fvMesh &mesh, Istream &schemeData),(mesh, schemeData))
 
 declareRunTimeSelectionTable (tmp, surfaceInterpolationScheme, MeshFlux,(const fvMesh &mesh, const surfaceScalarField &faceFlux, Istream &schemeData),(mesh, faceFlux, schemeData))
 
 surfaceInterpolationScheme (const fvMesh &mesh)
 Construct from mesh. More...
 
 surfaceInterpolationScheme (const surfaceInterpolationScheme &)=delete
 Disallow default bitwise copy construction. More...
 
virtual ~surfaceInterpolationScheme ()
 Destructor. More...
 
const fvMeshmesh () const
 Return mesh reference. More...
 
virtual tmp< SurfaceField< scalar > > interpolate (const VolField< scalar > &) const
 Return the face-interpolate of the given cell field. More...
 
tmp< SurfaceField< scalar > > interpolate (const tmp< VolField< scalar >> &) const
 Return the face-interpolate of the given tmp cell field. More...
 
virtual tmp< SurfaceField< typename innerProduct< vector, scalar >::type > > dotInterpolate (const surfaceVectorField &Sf, const VolField< scalar > &vf) const
 Return the face-interpolate of the given cell field. More...
 
tmp< SurfaceField< typename innerProduct< vector, scalar >::type > > dotInterpolate (const surfaceVectorField &Sf, const tmp< VolField< scalar >> &) const
 Return the face-interpolate of the given tmp cell field. More...
 
tmp< SurfaceField< typename innerProduct< vector, scalar >::type > > dotInterpolate (const surfaceVectorField &Sf, const VolField< scalar > &) const
 
virtual bool corrected () const
 Return true if this scheme uses an explicit correction. More...
 
virtual tmp< SurfaceField< scalar > > correction (const VolField< scalar > &) const
 Return the explicit correction to the face-interpolate. More...
 
void operator= (const surfaceInterpolationScheme &)=delete
 Disallow default bitwise assignment. More...
 
- Public Member Functions inherited from refCount
int count () const
 Return the current reference count. More...
 
bool unique () const
 Return true if the reference count is zero. More...
 
void operator++ ()
 Increment the reference count. More...
 
void operator++ (int)
 Increment the reference count. More...
 
void operator-- ()
 Decrement the reference count. More...
 
void operator-- (int)
 Decrement the reference count. More...
 
void operator= (const refCount &)=delete
 Disallow bitwise assignment. More...
 

Protected Member Functions

void setCellAlphaf (const label celli, const scalarField &phi, scalarField &alphaf, boolList &correctedFaces, const DynamicList< scalar > &cellAlphaf, const fvMesh &mesh) const
 Set alphaPhi for the faces of the given cell. More...
 
tmp< surfaceScalarFieldsurfaceAlpha (const volScalarField &alpha, const surfaceScalarField &phi, scalarField &splicedTvff, const bool unweighted, const scalar tol, const bool isMPLIC=true) const
 Return alpha interpolation. More...
 
- Protected Member Functions inherited from refCount
 refCount ()
 Construct null initialising count to 0. More...
 
 refCount (const refCount &)=delete
 Disallow copy. More...
 

Protected Attributes

const surfaceScalarFieldphi_
 

Additional Inherited Members

- Static Public Member Functions inherited from surfaceInterpolationScheme< scalar >
static tmp< surfaceInterpolationScheme< scalar > > New (const fvMesh &mesh, Istream &schemeData)
 Return new tmp interpolation scheme. More...
 
static tmp< surfaceInterpolationScheme< scalar > > New (const fvMesh &mesh, const surfaceScalarField &faceFlux, Istream &schemeData)
 Return new tmp interpolation scheme. More...
 
static tmp< SurfaceField< scalar > > interpolate (const VolField< scalar > &, const tmp< surfaceScalarField > &, const tmp< surfaceScalarField > &)
 Return the face-interpolate of the given cell field. More...
 
static tmp< SurfaceField< scalar > > interpolate (const VolField< scalar > &, const tmp< surfaceScalarField > &)
 Return the face-interpolate of the given cell field. More...
 
static tmp< SurfaceField< typename innerProduct< typename SFType::value_type, scalar >::type > > dotInterpolate (const SFType &Sf, const VolField< scalar > &vf, const tmp< surfaceScalarField > &tlambdas)
 Return the face-interpolate of the given cell field. More...
 

Detailed Description

Multicut Piecewise-Linear Interface Calculation (MPLIC) corrected scheme is a surface interpolation scheme for flux calculation in advection of a bounded variable, e.g. phase fraction and for interface capturing in the volume of fluid (VoF) method.

The interface is represented by multiple cuts which split each cell to match the volume fraction of the phase in the cell. The cut planes are oriented according to the point field of the local phase fraction. The phase fraction at each cell face - the interpolated value - is then calculated from the face area on either side of the cuts.

Three progressively more complex algorithms are used to ensure the cell volume fraction is accurately reproduced:

  1. single cut: cuts all the cell faces regardless the order
  2. multi cut: topological face-edge-face walk which can split cell into multiple sub-volumes
  3. tetrahedron cut: decomposes cell into tetrahedrons which are cut

Example:

divSchemes
{
    .
    .
    div(phi,alpha)      Gauss MPLIC;
    .
    .
}
See also
Foam::MPLICU Foam::PLIC Foam::PLICU Foam::interfaceCompression
Source files

Definition at line 83 of file MPLIC.H.

Constructor & Destructor Documentation

◆ MPLIC()

MPLIC ( const fvMesh mesh,
const surfaceScalarField faceFlux,
Istream is 
)
inline

Construct from faceFlux and Istream.

Definition at line 128 of file MPLIC.H.

Member Function Documentation

◆ setCellAlphaf()

void setCellAlphaf ( const label  celli,
const scalarField phi,
scalarField alphaf,
boolList correctedFaces,
const DynamicList< scalar > &  cellAlphaf,
const fvMesh mesh 
) const
protected

Set alphaPhi for the faces of the given cell.

Definition at line 46 of file MPLIC.C.

References primitiveMesh::cells(), polyMesh::faceOwner(), forAll, surfaceInterpolationScheme< scalar >::mesh(), and Foam::sign().

Here is the call graph for this function:

◆ surfaceAlpha()

◆ TypeName()

TypeName ( "MPLIC"  )

Runtime type information.

◆ weights()

virtual tmp<surfaceScalarField> weights ( const VolField< scalar > &  ) const
inlinevirtual

Return the interpolation weighting factors.

Implements surfaceInterpolationScheme< scalar >.

Reimplemented in MPLICU.

Definition at line 143 of file MPLIC.H.

References NotImplemented.

◆ interpolate()

Foam::tmp< Foam::surfaceScalarField > interpolate ( const VolField< scalar > &  vf) const
virtual

Return the face-interpolate of the given cell field.

Reimplemented in MPLICU.

Definition at line 236 of file MPLIC.C.

References Foam::e, Foam::fvc::interpolate(), and Foam::name().

Here is the call graph for this function:

◆ operator=()

void operator= ( const MPLIC )
delete

Disallow default bitwise assignment.

Member Data Documentation

◆ phi_

const surfaceScalarField& phi_
protected

Definition at line 91 of file MPLIC.H.

Referenced by MPLICU::interpolate().


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