MomentumTransferPhaseSystem< BasePhaseSystem > Class Template Reference

Class which models interfacial momentum transfer between a number of phases. Drag, virtual mass, lift, wall lubrication and turbulent dispersion are all modelled. The explicit contribution from the drag is omitted from the transfer matrices, as this forms part of the solution of the pressure equation. More...

Inheritance diagram for MomentumTransferPhaseSystem< BasePhaseSystem >:
Collaboration diagram for MomentumTransferPhaseSystem< BasePhaseSystem >:

Public Member Functions

 MomentumTransferPhaseSystem (const fvMesh &)
 Construct from fvMesh. More...
 
virtual ~MomentumTransferPhaseSystem ()
 Destructor. More...
 
virtual autoPtr< phaseSystem::momentumTransferTablemomentumTransfer ()
 Return the momentum transfer matrices for the cell-based algorithm. More...
 
virtual autoPtr< phaseSystem::momentumTransferTablemomentumTransferf ()
 As momentumTransfer, but for the face-based algorithm. More...
 
virtual PtrList< surfaceScalarFieldFs () const
 Return the explicit force fluxes for the cell-based algorithm, that. More...
 
virtual PtrList< surfaceScalarFieldFfs () const
 As Fs, but for the face-based algorithm. More...
 
virtual void invADVs (const PtrList< volScalarField > &As, PtrList< volVectorField > &HVms, PtrList< PtrList< volScalarField >> &invADVs, PtrList< PtrList< surfaceScalarField >> &invADVfs) const
 Return the inverse of the central + drag + virtual mass. More...
 
virtual PtrList< PtrList< surfaceScalarField > > invADVfs (const PtrList< surfaceScalarField > &Afs, PtrList< surfaceScalarField > &HVmfs) const
 Return the inverse of the central + drag + virtual mass. More...
 
virtual bool implicitPhasePressure (const phaseModel &phase) const
 Returns true if the phase pressure is treated implicitly. More...
 
virtual bool implicitPhasePressure () const
 Returns true if the phase pressure is treated implicitly. More...
 
virtual tmp< surfaceScalarFieldalphaDByAf (const PtrList< volScalarField > &rAs) const
 Return the phase diffusivity. More...
 
virtual PtrList< surfaceScalarFieldddtCorrs () const
 Return the flux corrections for the cell-based algorithm. These. More...
 
virtual void dragCorrs (PtrList< volVectorField > &dragCorrs, PtrList< surfaceScalarField > &dragCorrf) const
 Set the cell and faces drag correction fields. More...
 
virtual bool read ()
 Read base phaseProperties dictionary. More...
 

Protected Member Functions

void addDmdtUfs (const phaseSystem::dmdtfTable &dmdtfs, phaseSystem::momentumTransferTable &eqns)
 Add momentum transfer terms which result from bulk mass transfers. More...
 
void addTmpField (tmp< surfaceScalarField > &result, const tmp< surfaceScalarField > &field) const
 
void invADVs (List< UPtrList< scalarField >> &ADVs) const
 Invert the ADVs matrix inplace. More...
 
template<template< class > class PatchField, class GeoMesh >
void invADVs (PtrList< PtrList< GeometricField< scalar, PatchField, GeoMesh >>> &ADVs) const
 Invert the ADVs matrix inplace. More...
 

Detailed Description

template<class BasePhaseSystem>
class Foam::MomentumTransferPhaseSystem< BasePhaseSystem >

Class which models interfacial momentum transfer between a number of phases. Drag, virtual mass, lift, wall lubrication and turbulent dispersion are all modelled. The explicit contribution from the drag is omitted from the transfer matrices, as this forms part of the solution of the pressure equation.

Source files

Definition at line 61 of file MomentumTransferPhaseSystem.H.

Constructor & Destructor Documentation

◆ MomentumTransferPhaseSystem()

Construct from fvMesh.

Definition at line 112 of file MomentumTransferPhaseSystem.C.

◆ ~MomentumTransferPhaseSystem()

Destructor.

Definition at line 131 of file MomentumTransferPhaseSystem.C.

Member Function Documentation

◆ addDmdtUfs()

void addDmdtUfs ( const phaseSystem::dmdtfTable dmdtfs,
phaseSystem::momentumTransferTable eqns 
)
protected

Add momentum transfer terms which result from bulk mass transfers.

Definition at line 52 of file MomentumTransferPhaseSystem.C.

References forAllConstIter, phaseModel::name(), Foam::negPart(), phaseInterface::phase1(), phaseInterface::phase2(), Foam::posPart(), Foam::fvc::Sp(), phaseModel::stationary(), phaseModel::U(), and phaseModel::URef().

Here is the call graph for this function:

◆ addTmpField()

void addTmpField ( tmp< surfaceScalarField > &  result,
const tmp< surfaceScalarField > &  field 
) const
protected

Definition at line 85 of file MomentumTransferPhaseSystem.C.

References tmp< T >::isTmp(), tmp< T >::ref(), and tmp< T >::valid().

Here is the call graph for this function:

◆ invADVs() [1/3]

void invADVs ( List< UPtrList< scalarField >> &  ADVs) const
protected

Invert the ADVs matrix inplace.

Definition at line 433 of file MomentumTransferPhaseSystem.C.

References forAll, Foam::LUBacksubstitute(), Foam::LUDecompose(), n, and Foam::Zero.

Referenced by MomentumTransferPhaseSystem< BasePhaseSystem >::invADVs().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ invADVs() [2/3]

void invADVs ( PtrList< PtrList< GeometricField< scalar, PatchField, GeoMesh >>> &  ADVs) const
protected

Invert the ADVs matrix inplace.

Definition at line 476 of file MomentumTransferPhaseSystem.C.

References Foam::dimensions(), Foam::dimless, forAll, n, patchi, and List< T >::setSize().

Here is the call graph for this function:

◆ momentumTransfer()

Return the momentum transfer matrices for the cell-based algorithm.

This includes implicit and explicit forces that add into the cell UEqn in the normal way.

Definition at line 140 of file MomentumTransferPhaseSystem.C.

References Foam::dimMass, Foam::dimTime, Foam::dimVelocity, forAll, HashTable< T, Key, Hash >::insert(), phaseModel::name(), and phaseModel::U().

Here is the call graph for this function:

◆ momentumTransferf()

As momentumTransfer, but for the face-based algorithm.

Definition at line 167 of file MomentumTransferPhaseSystem.C.

References Foam::dimMass, Foam::dimTime, Foam::dimVelocity, forAll, HashTable< T, Key, Hash >::insert(), phaseModel::name(), and phaseModel::U().

Here is the call graph for this function:

◆ Fs()

Return the explicit force fluxes for the cell-based algorithm, that.

do not depend on phase mass/volume fluxes, and can therefore be evaluated outside the corrector loop. This includes things like lift, turbulent dispersion, and wall lubrication.

Definition at line 194 of file MomentumTransferPhaseSystem.C.

References Foam::addField(), D, Foam::constant::physicoChemical::F, Foam::fvc::flux(), forAll, forAllConstIter, Foam::fvc::interpolate(), Foam::max(), phaseModel::pPrimef(), and Foam::fvc::snGrad().

Here is the call graph for this function:

◆ Ffs()

As Fs, but for the face-based algorithm.

Definition at line 314 of file MomentumTransferPhaseSystem.C.

References Foam::addField(), D, forAll, forAllConstIter, Foam::fvc::interpolate(), Foam::max(), phaseModel::pPrimef(), and Foam::fvc::snGrad().

Here is the call graph for this function:

◆ invADVs() [3/3]

◆ invADVfs()

◆ implicitPhasePressure() [1/2]

bool implicitPhasePressure ( const phaseModel phase) const
virtual

Returns true if the phase pressure is treated implicitly.

in the phase fraction equation

Definition at line 942 of file MomentumTransferPhaseSystem.C.

◆ implicitPhasePressure() [2/2]

bool implicitPhasePressure
virtual

Returns true if the phase pressure is treated implicitly.

in the phase fraction equation for any phase

Definition at line 956 of file MomentumTransferPhaseSystem.C.

References forAll.

◆ alphaDByAf()

Foam::tmp< Foam::surfaceScalarField > alphaDByAf ( const PtrList< volScalarField > &  rAs) const
virtual

Return the phase diffusivity.

divided by the momentum central coefficient

Definition at line 974 of file MomentumTransferPhaseSystem.C.

References D, forAll, forAllConstIter, phaseModel::index(), Foam::fvc::interpolate(), Foam::max(), and phaseModel::pPrimef().

Here is the call graph for this function:

◆ ddtCorrs()

Return the flux corrections for the cell-based algorithm. These.

depend on phase mass/volume fluxes, and must therefore be evaluated inside the corrector loop.

Definition at line 1039 of file MomentumTransferPhaseSystem.C.

References Foam::addField(), Foam::fvc::ddtCorr(), singleRegionSolutionControl::dict(), forAll, forAllConstIter, phaseModel::index(), Foam::fvc::interpolate(), K, dictionary::lookupOrDefault(), Foam::max(), phaseModel::phi(), pimple(), phaseModel::residualAlpha(), phaseModel::rho(), PtrList< T >::set(), phaseModel::U(), and phaseModel::Uf().

Here is the call graph for this function:

◆ dragCorrs()

void dragCorrs ( PtrList< volVectorField > &  dragCorrs,
PtrList< surfaceScalarField > &  dragCorrf 
) const
virtual

Set the cell and faces drag correction fields.

Definition at line 1134 of file MomentumTransferPhaseSystem.C.

References Foam::addField(), forAll, forAllConstIter, phaseModel::index(), Foam::fvc::interpolate(), K, K1, Foam::max(), phaseModel::phi(), Foam::fvc::reconstruct(), phaseModel::residualAlpha(), and PtrList< T >::set().

Here is the call graph for this function:

◆ read()

bool read
virtual

Read base phaseProperties dictionary.

Definition at line 1197 of file MomentumTransferPhaseSystem.C.

References Foam::blockMeshTools::read().

Here is the call graph for this function:

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