verticalDamping Class Reference

This fvModel applies an explicit forcing force to components of the vector field in the direction of gravity. Its intended purpose is to damp the vertical motions of an interface in the region approaching an outlet so that no reflections are generated. More...

Inheritance diagram for verticalDamping:
Collaboration diagram for verticalDamping:

Public Member Functions

 TypeName ("verticalDamping")
 Runtime type information. More...
 
 verticalDamping (const word &name, const word &modelType, const fvMesh &mesh, const dictionary &dict)
 Construct from components. More...
 
virtual ~verticalDamping ()
 Destructor. More...
 
virtual bool read (const dictionary &dict)
 Read dictionary. More...
 
virtual wordList addSupFields () const
 Return the list of fields for which the fvModel adds source term. More...
 
virtual void addSup (const volVectorField &U, fvMatrix< vector > &eqn) const
 Source term to momentum equation. More...
 
virtual void addSup (const volScalarField &rho, const volVectorField &U, fvMatrix< vector > &eqn) const
 Source term to compressible momentum equation. More...
 
virtual void addSup (const volScalarField &alpha, const volScalarField &rho, const volVectorField &U, fvMatrix< vector > &eqn) const
 Source term to phase momentum equation. More...
 
virtual bool movePoints ()
 Update for mesh motion. More...
 
virtual void topoChange (const polyTopoChangeMap &)
 Update topology using the given map. More...
 
virtual void mapMesh (const polyMeshMap &)
 Update from another mesh using the given map. More...
 
virtual void distribute (const polyDistributionMap &)
 Redistribute or update using the given distribution map. More...
 
- Public Member Functions inherited from forcing
 TypeName ("forcing")
 Runtime type information. More...
 
 forcing (const word &name, const word &modelType, const fvMesh &mesh, const dictionary &dict)
 Construct from components. More...
 
virtual ~forcing ()
 Destructor. More...
 
- Public Member Functions inherited from fvModel
 TypeName ("fvModel")
 Runtime type information. More...
 
 declareRunTimeSelectionTable (autoPtr, fvModel, dictionary,(const word &name, const word &modelType, const fvMesh &mesh, const dictionary &dict),(name, modelType, mesh, dict))
 
 fvModel (const word &name, const word &modelType, const fvMesh &mesh, const dictionary &dict)
 Construct from components. More...
 
autoPtr< fvModelclone () const
 Return clone. More...
 
virtual ~fvModel ()
 Destructor. More...
 
const wordname () const
 Return const access to the source name. More...
 
const fvMeshmesh () const
 Return const access to the mesh database. More...
 
const dictionarycoeffs () const
 Return dictionary. More...
 
virtual bool addsSupToField (const word &fieldName) const
 Return true if the fvModel adds a source term to the given. More...
 
virtual scalar maxDeltaT () const
 Return the maximum time-step for stable operation. More...
 
virtual void addSup (fvMatrix< scalar > &eqn) const
 Add a source term to a field-less proxy equation. More...
 
template<class Type >
tmp< fvMatrix< Type > > sourceProxy (const VolField< Type > &eqnField) const
 Add a source term to an equation. More...
 
template<class Type >
tmp< fvMatrix< Type > > source (const VolField< Type > &field) const
 Return source for an equation. More...
 
template<class Type >
tmp< fvMatrix< Type > > sourceProxy (const VolField< Type > &field, const VolField< Type > &eqnField) const
 Return source for an equation. More...
 
template<class Type >
tmp< fvMatrix< Type > > source (const volScalarField &rho, const VolField< Type > &field) const
 Return source for a compressible equation. More...
 
template<class Type >
tmp< fvMatrix< Type > > sourceProxy (const volScalarField &rho, const VolField< Type > &field, const VolField< Type > &eqnField) const
 Return source for a compressible equation. More...
 
template<class Type >
tmp< fvMatrix< Type > > source (const volScalarField &alpha, const volScalarField &rho, const VolField< Type > &field) const
 Return source for a phase equation. More...
 
template<class Type >
tmp< fvMatrix< Type > > sourceProxy (const volScalarField &alpha, const volScalarField &rho, const VolField< Type > &field, const VolField< Type > &eqnField) const
 Return source for a phase equation. More...
 
template<class Type >
tmp< fvMatrix< Type > > source (const volScalarField &alpha, const geometricOneField &rho, const VolField< Type > &field) const
 Return source for a phase equation. More...
 
template<class Type >
tmp< fvMatrix< Type > > source (const geometricOneField &alpha, const volScalarField &rho, const VolField< Type > &field) const
 Return source for a phase equation. More...
 
template<class Type >
tmp< fvMatrix< Type > > source (const geometricOneField &alpha, const geometricOneField &rho, const VolField< Type > &field) const
 Return source for a phase equation. More...
 
template<class Type >
tmp< fvMatrix< Type > > d2dt2 (const VolField< Type > &field) const
 Return source for an equation with a second time derivative. More...
 
virtual void preUpdateMesh ()
 Prepare for mesh update. More...
 
virtual void correct ()
 Correct the fvModel. More...
 
virtual bool write (const bool write=true) const
 Write fvModel data. More...
 
template<class AlphaRhoFieldType , class ... AlphaRhoFieldTypes>
Foam::dimensionSet sourceDims (const dimensionSet &ds, const AlphaRhoFieldType &alphaRhoField, const AlphaRhoFieldTypes &... alphaRhoFields)
 
template<class AlphaRhoFieldType , class ... AlphaRhoFieldTypes>
const Foam::wordfieldName (const AlphaRhoFieldType &alphaRhoField, const AlphaRhoFieldTypes &... alphaRhoFields)
 
template<class AlphaRhoFieldType >
const Foam::wordfieldName (const AlphaRhoFieldType &alphaRhoField)
 
template<class Type , class ... AlphaRhoFieldTypes>
Foam::tmp< Foam::fvMatrix< Type > > sourceTerm (const VolField< Type > &eqnField, const dimensionSet &ds, const AlphaRhoFieldTypes &... alphaRhoFields) const
 
template<class Type >
Foam::tmp< Foam::fvMatrix< Type > > sourceProxy (const VolField< Type > &eqnField) const
 
template<class Type >
Foam::tmp< Foam::fvMatrix< Type > > source (const VolField< Type > &field) const
 
template<class Type >
Foam::tmp< Foam::fvMatrix< Type > > sourceProxy (const VolField< Type > &field, const VolField< Type > &eqnField) const
 
template<class Type >
Foam::tmp< Foam::fvMatrix< Type > > source (const volScalarField &rho, const VolField< Type > &field) const
 
template<class Type >
Foam::tmp< Foam::fvMatrix< Type > > sourceProxy (const volScalarField &rho, const VolField< Type > &field, const VolField< Type > &eqnField) const
 
template<class Type >
Foam::tmp< Foam::fvMatrix< Type > > source (const volScalarField &alpha, const volScalarField &rho, const VolField< Type > &field) const
 
template<class Type >
Foam::tmp< Foam::fvMatrix< Type > > sourceProxy (const volScalarField &alpha, const volScalarField &rho, const VolField< Type > &field, const VolField< Type > &eqnField) const
 
template<class Type >
Foam::tmp< Foam::fvMatrix< Type > > source (const geometricOneField &alpha, const geometricOneField &rho, const VolField< Type > &field) const
 
template<class Type >
Foam::tmp< Foam::fvMatrix< Type > > source (const volScalarField &alpha, const geometricOneField &rho, const VolField< Type > &field) const
 
template<class Type >
Foam::tmp< Foam::fvMatrix< Type > > source (const geometricOneField &alpha, const volScalarField &rho, const VolField< Type > &field) const
 
template<class Type >
Foam::tmp< Foam::fvMatrix< Type > > d2dt2 (const VolField< Type > &field) const
 

Additional Inherited Members

- Static Public Member Functions inherited from fvModel
template<class AlphaRhoFieldType , class ... AlphaRhoFieldTypes>
static dimensionSet sourceDims (const dimensionSet &ds, const AlphaRhoFieldType &alphaRhoField, const AlphaRhoFieldTypes &... alphaRhoFields)
 Return the dimensions of the matrix of a source term. More...
 
static const dimensionSetsourceDims (const dimensionSet &ds)
 Return the dimensions of the matrix of a source term (base. More...
 
template<class AlphaRhoFieldType , class ... AlphaRhoFieldTypes>
static const wordfieldName (const AlphaRhoFieldType &alphaRhoField, const AlphaRhoFieldTypes &... alphaRhoFields)
 Return the name of the field associated with a source term. More...
 
template<class AlphaRhoFieldType >
static const wordfieldName (const AlphaRhoFieldType &alphaRhoField)
 Return the name of the field associated with a source term (base. More...
 
static const wordfieldName ()
 Return the name of the field associated with a source term. Special. More...
 
static autoPtr< fvModelNew (const word &name, const fvMesh &mesh, const dictionary &dict)
 Return a reference to the selected fvModel. More...
 
- Protected Member Functions inherited from forcing
void readCoeffs ()
 Non-virtual read. More...
 
void readLambda ()
 Read the forcing coefficients. More...
 
dimensionedScalar regionLength () const
 Calculate and return the volume average forcing region length. More...
 
virtual tmp< volScalarField::Internalscale () const
 Return the scale distribution. More...
 
virtual tmp< volScalarField::InternalforceCoeff () const
 Return the force coefficient. More...
 
void writeForceFields () const
 Optionally write the forcing fields: More...
 
- Protected Member Functions inherited from fvModel
template<class Type >
void addSupType (const VolField< Type > &field, fvMatrix< Type > &eqn) const
 Add a source term to an equation. More...
 
template<class Type >
void addSupType (const volScalarField &rho, const VolField< Type > &field, fvMatrix< Type > &eqn) const
 Add a source term to a compressible equation. More...
 
template<class Type >
void addSupType (const volScalarField &alpha, const volScalarField &rho, const VolField< Type > &field, fvMatrix< Type > &eqn) const
 Add a source term to a phase equation. More...
 
template<class Type , class ... AlphaRhoFieldTypes>
tmp< fvMatrix< Type > > sourceTerm (const VolField< Type > &eqnField, const dimensionSet &ds, const AlphaRhoFieldTypes &... alphaRhoFields) const
 Return a source for an equation. More...
 
- Protected Attributes inherited from forcing
bool writeForceFields_
 Optionally write the forcing fields. More...
 
dimensionedScalar lambda_
 Forcing coefficient [1/s]. More...
 
dimensionedScalar lambdaBoundary_
 Optional boundary forcing coefficient [1/s]. More...
 
autoPtr< Function1< scalar > > scale_
 The scaling function. More...
 
vectorField origins_
 Origins of the scaling coordinate. More...
 
vectorField directions_
 Directions of increasing scaling coordinate. More...
 

Detailed Description

This fvModel applies an explicit forcing force to components of the vector field in the direction of gravity. Its intended purpose is to damp the vertical motions of an interface in the region approaching an outlet so that no reflections are generated.

Damping is achieved by applying a force to the momentum equation proportional to the momentum of the flow in the direction of gravity. The constant of proportionality is given by a coefficient $\lambda$ which has units of inverse-time. In the absence of any other forces this would generate an exponential decay of the vertical velocity.

\[ \frac{d (m u_z)}{d t} = - \lambda m u_z \]

\[ u_z = u_{z0} e^{- \lambda t} \]

The coefficient $\lambda$ should be set based on the desired level of forcing and the residence time of a perturbation through the forcing zone. For example, if waves moving at 2 [m/s] are travelling through a forcing zone 8 [m] in length, then the residence time is 4 [s]. If it is deemed necessary to damp for 5 time-scales, then $\lambda$ should be set to equal 5/(4 [s]) = 1.2 [1/s].

Usage
Example usage:
verticalDamping1
{
    type            verticalDamping;

    libs            ("libwaves.so");

    // Define the line along which to apply the graduation
    origin          (1200 0 0);
    direction       (1 0 0);

    // Or, define multiple lines
    // origins         ((1200 0 0) (1200 -300 0) (1200 300 0));
    // directions      ((1 0 0) (0 -1 0) (0 1 0));

    scale
    {
        type        halfCosineRamp;
        start       0;
        duration    600;
    }

    lambda          [0 0 -1 0 0 0 0] 1; // Damping coefficient

    timeStart       0;
    duration        1e6;

    // Write the forcing fields: forcing:scale, forcing:forceCoeff
    writeForceFields true;
}
See also
Foam::fv::forcing Foam::fv::isotropicDamping
Source files

Definition at line 111 of file verticalDamping.H.

Constructor & Destructor Documentation

◆ verticalDamping()

verticalDamping ( const word name,
const word modelType,
const fvMesh mesh,
const dictionary dict 
)

Construct from components.

Definition at line 68 of file verticalDamping.C.

References forcing::writeForceFields().

Here is the call graph for this function:

◆ ~verticalDamping()

virtual ~verticalDamping ( )
inlinevirtual

Destructor.

Definition at line 153 of file verticalDamping.H.

Member Function Documentation

◆ TypeName()

TypeName ( "verticalDamping"  )

Runtime type information.

◆ read()

bool read ( const dictionary dict)
virtual

Read dictionary.

Reimplemented from forcing.

Definition at line 142 of file verticalDamping.C.

References dict, and forcing::read().

Here is the call graph for this function:

◆ addSupFields()

Foam::wordList addSupFields ( ) const
virtual

Return the list of fields for which the fvModel adds source term.

to the transport equation

Reimplemented from fvModel.

Definition at line 85 of file verticalDamping.C.

◆ addSup() [1/3]

void addSup ( const volVectorField U,
fvMatrix< vector > &  eqn 
) const
virtual

Source term to momentum equation.

Definition at line 91 of file verticalDamping.C.

References Foam::add(), and U.

Here is the call graph for this function:

◆ addSup() [2/3]

void addSup ( const volScalarField rho,
const volVectorField U,
fvMatrix< vector > &  eqn 
) const
virtual

Source term to compressible momentum equation.

Definition at line 101 of file verticalDamping.C.

References Foam::add(), rho, and U.

Here is the call graph for this function:

◆ addSup() [3/3]

void addSup ( const volScalarField alpha,
const volScalarField rho,
const volVectorField U,
fvMatrix< vector > &  eqn 
) const
virtual

Source term to phase momentum equation.

Definition at line 112 of file verticalDamping.C.

References Foam::add(), alpha(), rho, and U.

Here is the call graph for this function:

◆ movePoints()

bool movePoints ( )
virtual

Update for mesh motion.

Implements fvModel.

Definition at line 124 of file verticalDamping.C.

◆ topoChange()

void topoChange ( const polyTopoChangeMap )
virtual

Update topology using the given map.

Implements fvModel.

Definition at line 130 of file verticalDamping.C.

◆ mapMesh()

void mapMesh ( const polyMeshMap map)
virtual

Update from another mesh using the given map.

Implements fvModel.

Definition at line 134 of file verticalDamping.C.

◆ distribute()

void distribute ( const polyDistributionMap )
virtual

Redistribute or update using the given distribution map.

Implements fvModel.

Definition at line 138 of file verticalDamping.C.


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