CrankNicolsonDdtScheme< Type > Class Template Reference

Second-oder Crank-Nicolson implicit ddt using the current and previous time-step fields as well as the previous time-step ddt. More...

Inheritance diagram for CrankNicolsonDdtScheme< Type >:
[legend]
Collaboration diagram for CrankNicolsonDdtScheme< Type >:
[legend]

Public Types

typedef ddtScheme< Type >::fluxFieldType fluxFieldType
 
- Public Types inherited from ddtScheme< Type >
typedef GeometricField< typename flux< Type >::type, fvsPatchField, surfaceMeshfluxFieldType
 

Public Member Functions

 TypeName ("CrankNicolson")
 Runtime type information. More...
 
 CrankNicolsonDdtScheme (const fvMesh &mesh)
 Construct from mesh. More...
 
 CrankNicolsonDdtScheme (const fvMesh &mesh, Istream &is)
 Construct from mesh and Istream. More...
 
const fvMeshmesh () const
 Return mesh reference. More...
 
scalar ocCoeff () const
 Return the current off-centreing coefficient. More...
 
tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt (const dimensioned< Type > &)
 
tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt (const GeometricField< Type, fvPatchField, volMesh > &)
 
tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt (const dimensionedScalar &, const GeometricField< Type, fvPatchField, volMesh > &)
 
tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt (const volScalarField &, const GeometricField< Type, fvPatchField, volMesh > &)
 
tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt (const volScalarField &alpha, const volScalarField &rho, const GeometricField< Type, fvPatchField, volMesh > &psi)
 
tmp< fvMatrix< Type > > fvmDdt (const GeometricField< Type, fvPatchField, volMesh > &)
 
tmp< fvMatrix< Type > > fvmDdt (const dimensionedScalar &, const GeometricField< Type, fvPatchField, volMesh > &)
 
tmp< fvMatrix< Type > > fvmDdt (const volScalarField &, const GeometricField< Type, fvPatchField, volMesh > &)
 
tmp< fvMatrix< Type > > fvmDdt (const volScalarField &alpha, const volScalarField &rho, const GeometricField< Type, fvPatchField, volMesh > &psi)
 
tmp< fluxFieldTypefvcDdtUfCorr (const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)
 
tmp< fluxFieldTypefvcDdtPhiCorr (const GeometricField< Type, fvPatchField, volMesh > &U, const fluxFieldType &phi)
 
tmp< fluxFieldTypefvcDdtUfCorr (const volScalarField &rho, const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)
 
tmp< fluxFieldTypefvcDdtPhiCorr (const volScalarField &rho, const GeometricField< Type, fvPatchField, volMesh > &U, const fluxFieldType &phi)
 
tmp< surfaceScalarFieldmeshPhi (const GeometricField< Type, fvPatchField, volMesh > &)
 
template<class GeoField >
CrankNicolsonDdtScheme< Type >::template DDt0Field< GeoField > & ddt0_ (const word &name, const dimensionSet &dims)
 
tmp< surfaceScalarFieldfvcDdtUfCorr (const GeometricField< scalar, fvPatchField, volMesh > &U, const GeometricField< scalar, fvsPatchField, surfaceMesh > &Uf)
 
tmp< surfaceScalarFieldfvcDdtPhiCorr (const volScalarField &U, const surfaceScalarField &phi)
 
tmp< surfaceScalarFieldfvcDdtUfCorr (const volScalarField &rho, const volScalarField &U, const surfaceScalarField &Uf)
 
tmp< surfaceScalarFieldfvcDdtPhiCorr (const volScalarField &rho, const volScalarField &U, const surfaceScalarField &phi)
 
- Public Member Functions inherited from ddtScheme< Type >
virtual const wordtype () const =0
 Runtime type information. More...
 
 declareRunTimeSelectionTable (tmp, ddtScheme, Istream,(const fvMesh &mesh, Istream &schemeData),(mesh, schemeData))
 
 ddtScheme (const fvMesh &mesh)
 Construct from mesh. More...
 
 ddtScheme (const fvMesh &mesh, Istream &is)
 Construct from mesh and Istream. More...
 
virtual ~ddtScheme ()=default
 Destructor. More...
 
const fvMeshmesh () const
 Return mesh reference. More...
 
virtual tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt (const dimensioned< Type > &)=0
 
virtual tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt (const GeometricField< Type, fvPatchField, volMesh > &)=0
 
virtual tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt (const dimensionedScalar &, const GeometricField< Type, fvPatchField, volMesh > &)=0
 
virtual tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt (const volScalarField &, const GeometricField< Type, fvPatchField, volMesh > &)=0
 
virtual tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt (const volScalarField &alpha, const volScalarField &rho, const GeometricField< Type, fvPatchField, volMesh > &)=0
 
virtual tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > fvcDdt (const GeometricField< Type, fvsPatchField, surfaceMesh > &)
 
virtual tmp< fvMatrix< Type > > fvmDdt (const GeometricField< Type, fvPatchField, volMesh > &)=0
 
virtual tmp< fvMatrix< Type > > fvmDdt (const dimensionedScalar &, const GeometricField< Type, fvPatchField, volMesh > &)=0
 
virtual tmp< fvMatrix< Type > > fvmDdt (const volScalarField &, const GeometricField< Type, fvPatchField, volMesh > &)=0
 
virtual tmp< fvMatrix< Type > > fvmDdt (const volScalarField &alpha, const volScalarField &rho, const GeometricField< Type, fvPatchField, volMesh > &vf)=0
 
tmp< surfaceScalarFieldfvcDdtPhiCoeff (const GeometricField< Type, fvPatchField, volMesh > &U, const fluxFieldType &phi, const fluxFieldType &phiCorr)
 
tmp< surfaceScalarFieldfvcDdtPhiCoeffExperimental (const GeometricField< Type, fvPatchField, volMesh > &U, const fluxFieldType &phi, const fluxFieldType &phiCorr)
 
tmp< surfaceScalarFieldfvcDdtPhiCoeff (const GeometricField< Type, fvPatchField, volMesh > &U, const fluxFieldType &phi, const fluxFieldType &phiCorr, const volScalarField &rho)
 
tmp< surfaceScalarFieldfvcDdtPhiCoeff (const GeometricField< Type, fvPatchField, volMesh > &U, const fluxFieldType &phi)
 
tmp< surfaceScalarFieldfvcDdtPhiCoeff (const GeometricField< Type, fvPatchField, volMesh > &rhoU, const fluxFieldType &phi, const volScalarField &rho)
 
virtual tmp< fluxFieldTypefvcDdtUfCorr (const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)=0
 
virtual tmp< fluxFieldTypefvcDdtPhiCorr (const GeometricField< Type, fvPatchField, volMesh > &U, const fluxFieldType &phi)=0
 
virtual tmp< fluxFieldTypefvcDdtUfCorr (const volScalarField &rho, const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)=0
 
virtual tmp< fluxFieldTypefvcDdtPhiCorr (const volScalarField &rho, const GeometricField< Type, fvPatchField, volMesh > &U, const fluxFieldType &phi)=0
 
virtual tmp< surfaceScalarFieldmeshPhi (const GeometricField< Type, fvPatchField, volMesh > &)=0
 
- Public Member Functions inherited from refCount
constexpr refCount () noexcept
 Default construct, initializing count to 0. More...
 
int count () const noexcept
 Return the current reference count. More...
 
bool unique () const noexcept
 Return true if the reference count is zero. More...
 
void operator++ () noexcept
 Increment the reference count. More...
 
void operator++ (int) noexcept
 Increment the reference count. More...
 
void operator-- () noexcept
 Decrement the reference count. More...
 
void operator-- (int) noexcept
 Decrement the reference count. More...
 
- Public Member Functions inherited from ddtSchemeBase
 ddtSchemeBase ()
 
 ddtSchemeBase ()
 

Additional Inherited Members

- Static Public Member Functions inherited from ddtScheme< Type >
static tmp< ddtScheme< Type > > New (const fvMesh &mesh, Istream &schemeData)
 Return a pointer to a new ddtScheme created on freestore. More...
 
- Static Public Attributes inherited from ddtSchemeBase
static bool experimentalDdtCorr
 
- Protected Member Functions inherited from ddtScheme< Type >
 ddtScheme (const ddtScheme &)=delete
 No copy construct. More...
 
void operator= (const ddtScheme &)=delete
 No copy assignment. More...
 
- Protected Attributes inherited from ddtScheme< Type >
const fvMeshmesh_
 
scalar ddtPhiCoeff_
 Input for fvcDdtPhiCoeff. More...
 

Detailed Description

template<class Type>
class Foam::fv::CrankNicolsonDdtScheme< Type >

Second-oder Crank-Nicolson implicit ddt using the current and previous time-step fields as well as the previous time-step ddt.

The Crank-Nicolson scheme is often unstable for complex flows in complex geometries and it is necessary to "off-centre" the scheme to stabilize it while retaining greater temporal accuracy than the first-order Euler-implicit scheme. Off-centering is specified via the mandatory coefficient ocCoeff in the range [0,1] following the scheme name e.g.

ddtSchemes
{
    default         CrankNicolson 0.9;
}

or with an optional "ramp" function to transition from the Euler scheme to Crank-Nicolson over a initial period to avoid start-up problems, e.g.

ddtSchemes
{
    default         CrankNicolson
    ocCoeff
    {
        type scale;
        scale linearRamp;
        duration 0.01;
        value 0.9;
    };
}

With a coefficient of 1 the scheme is fully centred and second-order, with a coefficient of 0 the scheme is equivalent to Euler-implicit. A coefficient of 0.9 has been found to be suitable for a range of cases for which higher-order accuracy in time is needed and provides similar accuracy and stability to the backward scheme. However, the backward scheme has been found to be more robust and provides formal second-order accuracy in time.

The advantage of the Crank-Nicolson scheme over backward is that only the new and old-time values are needed, the additional terms relating to the fluxes and sources are evaluated at the mid-point of the time-step which provides the opportunity to limit the fluxes in such a way as to ensure boundedness while maintaining greater accuracy in time compared to the Euler-implicit scheme. This approach is now used with MULES in the interFoam family of solvers. Boundedness cannot be guaranteed with the backward scheme.

Note
The Crank-Nicolson coefficient for the implicit part of the RHS is related to the off-centering coefficient by
    cnCoeff = 1.0/(1.0 + ocCoeff);
See also
Foam::Euler Foam::backward
Source files

Definition at line 116 of file CrankNicolsonDdtScheme.H.

Member Typedef Documentation

◆ fluxFieldType

typedef ddtScheme<Type>::fluxFieldType fluxFieldType

Definition at line 292 of file CrankNicolsonDdtScheme.H.

Constructor & Destructor Documentation

◆ CrankNicolsonDdtScheme() [1/2]

CrankNicolsonDdtScheme ( const fvMesh mesh)

Construct from mesh.

Definition at line 286 of file CrankNicolsonDdtScheme.C.

References CrankNicolsonDdtScheme< Type >::mesh(), polyMesh::moving(), and fvMesh::V00().

Here is the call graph for this function:

◆ CrankNicolsonDdtScheme() [2/2]

CrankNicolsonDdtScheme ( const fvMesh mesh,
Istream is 
)

Construct from mesh and Istream.

Definition at line 301 of file CrankNicolsonDdtScheme.C.

References dict, Foam::exit(), Foam::FatalIOError, FatalIOErrorInFunction, token::isNumber(), CrankNicolsonDdtScheme< Type >::mesh(), polyMesh::moving(), Time::New(), token::number(), CrankNicolsonDdtScheme< Type >::ocCoeff(), Istream::putBack(), and fvMesh::V00().

Here is the call graph for this function:

Member Function Documentation

◆ TypeName()

TypeName ( "CrankNicolson"  )

Runtime type information.

◆ mesh()

const fvMesh & mesh ( ) const
inline

Return mesh reference.

Definition at line 228 of file CrankNicolsonDdtScheme.H.

References ddtScheme< Type >::mesh().

Referenced by CrankNicolsonDdtScheme< Type >::CrankNicolsonDdtScheme(), and CrankNicolsonDdtScheme< Type >::ocCoeff().

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

◆ ocCoeff()

scalar ocCoeff ( ) const
inline

Return the current off-centreing coefficient.

Definition at line 234 of file CrankNicolsonDdtScheme.H.

References CrankNicolsonDdtScheme< Type >::mesh().

Referenced by CrankNicolsonDdtScheme< Type >::CrankNicolsonDdtScheme().

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

◆ fvcDdt() [1/5]

tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt ( const dimensioned< Type > &  dt)
virtual

Implements ddtScheme< Type >.

Definition at line 347 of file CrankNicolsonDdtScheme.C.

References dimensioned< Type >::dimensions(), Foam::dimTime, mesh, dimensioned< Type >::name(), tmp< T >::ref(), timeName, and Foam::Zero.

Here is the call graph for this function:

◆ fvcDdt() [2/5]

tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt ( const GeometricField< Type, fvPatchField, volMesh > &  vf)
virtual

Implements ddtScheme< Type >.

Definition at line 404 of file CrankNicolsonDdtScheme.C.

References GeometricField< Type, PatchField, GeoMesh >::boundaryField(), DimensionedField< Type, GeoMesh >::dimensions(), Foam::fv::ff(), mesh, IOobject::name(), GeometricField< Type, PatchField, GeoMesh >::oldTime(), timeName, and dimensioned< Type >::value().

Here is the call graph for this function:

◆ fvcDdt() [3/5]

◆ fvcDdt() [4/5]

◆ fvcDdt() [5/5]

tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt ( const volScalarField alpha,
const volScalarField rho,
const GeometricField< Type, fvPatchField, volMesh > &  psi 
)
virtual

◆ fvmDdt() [1/4]

tmp< fvMatrix< Type > > fvmDdt ( const GeometricField< Type, fvPatchField, volMesh > &  vf)
virtual

Implements ddtScheme< Type >.

Definition at line 811 of file CrankNicolsonDdtScheme.C.

References lduMatrix::diag(), DimensionedField< Type, GeoMesh >::dimensions(), Foam::dimTime, Foam::dimVol, Foam::fv::ff(), mesh, IOobject::name(), GeometricField< Type, PatchField, GeoMesh >::oldTime(), and fvMatrix< Type >::source().

Here is the call graph for this function:

◆ fvmDdt() [2/4]

tmp< fvMatrix< Type > > fvmDdt ( const dimensionedScalar rho,
const GeometricField< Type, fvPatchField, volMesh > &  vf 
)
virtual

Implements ddtScheme< Type >.

Definition at line 894 of file CrankNicolsonDdtScheme.C.

References lduMatrix::diag(), DimensionedField< Type, GeoMesh >::dimensions(), Foam::dimTime, Foam::dimVol, Foam::fv::ff(), mesh, IOobject::name(), GeometricField< Type, PatchField, GeoMesh >::oldTime(), rho, and fvMatrix< Type >::source().

Here is the call graph for this function:

◆ fvmDdt() [3/4]

tmp< fvMatrix< Type > > fvmDdt ( const volScalarField rho,
const GeometricField< Type, fvPatchField, volMesh > &  vf 
)
virtual

Implements ddtScheme< Type >.

Definition at line 976 of file CrankNicolsonDdtScheme.C.

References lduMatrix::diag(), DimensionedField< Type, GeoMesh >::dimensions(), Foam::dimTime, Foam::dimVol, Foam::fv::ff(), mesh, IOobject::name(), GeometricField< Type, PatchField, GeoMesh >::oldTime(), rho, and fvMatrix< Type >::source().

Here is the call graph for this function:

◆ fvmDdt() [4/4]

tmp< fvMatrix< Type > > fvmDdt ( const volScalarField alpha,
const volScalarField rho,
const GeometricField< Type, fvPatchField, volMesh > &  psi 
)
virtual

Implements ddtScheme< Type >.

Definition at line 1066 of file CrankNicolsonDdtScheme.C.

References alpha, lduMatrix::diag(), DimensionedField< Type, GeoMesh >::dimensions(), Foam::dimTime, Foam::dimVol, Foam::fv::ff(), mesh, IOobject::name(), GeometricField< Type, PatchField, GeoMesh >::oldTime(), rho, and fvMatrix< Type >::source().

Here is the call graph for this function:

◆ fvcDdtUfCorr() [1/4]

tmp< typename CrankNicolsonDdtScheme< Type >::fluxFieldType > fvcDdtUfCorr ( const GeometricField< Type, fvPatchField, volMesh > &  U,
const GeometricField< Type, fvsPatchField, surfaceMesh > &  Uf 
)
virtual

Implements ddtScheme< Type >.

Definition at line 1177 of file CrankNicolsonDdtScheme.C.

References Foam::fvc::interpolate(), mesh, timeName, U, and Uf.

Here is the call graph for this function:

◆ fvcDdtPhiCorr() [1/4]

tmp< typename CrankNicolsonDdtScheme< Type >::fluxFieldType > fvcDdtPhiCorr ( const GeometricField< Type, fvPatchField, volMesh > &  U,
const fluxFieldType phi 
)
virtual

Implements ddtScheme< Type >.

Definition at line 1238 of file CrankNicolsonDdtScheme.C.

References Foam::fvc::dotInterpolate(), mesh, phi, timeName, and U.

Here is the call graph for this function:

◆ fvcDdtUfCorr() [2/4]

tmp< typename CrankNicolsonDdtScheme< Type >::fluxFieldType > fvcDdtUfCorr ( const volScalarField rho,
const GeometricField< Type, fvPatchField, volMesh > &  U,
const GeometricField< Type, fvsPatchField, surfaceMesh > &  Uf 
)
virtual

Implements ddtScheme< Type >.

Definition at line 1301 of file CrankNicolsonDdtScheme.C.

References Foam::abort(), ddtCorr, Foam::dimVelocity, Foam::FatalError, FatalErrorInFunction, Foam::fvc::interpolate(), mesh, rho, timeName, U, and Uf.

Here is the call graph for this function:

◆ fvcDdtPhiCorr() [2/4]

tmp< typename CrankNicolsonDdtScheme< Type >::fluxFieldType > fvcDdtPhiCorr ( const volScalarField rho,
const GeometricField< Type, fvPatchField, volMesh > &  U,
const fluxFieldType phi 
)
virtual

Implements ddtScheme< Type >.

Definition at line 1458 of file CrankNicolsonDdtScheme.C.

References Foam::abort(), ddtCorr, Foam::dimArea, Foam::dimVelocity, Foam::fvc::dotInterpolate(), Foam::FatalError, FatalErrorInFunction, mesh, phi, rho, timeName, and U.

Here is the call graph for this function:

◆ meshPhi()

tmp< surfaceScalarField > meshPhi ( const GeometricField< Type, fvPatchField, volMesh > &  vf)
virtual

Implements ddtScheme< Type >.

Definition at line 1603 of file CrankNicolsonDdtScheme.C.

References Foam::dimVolume, mesh, Foam::name(), IOobject::NO_READ, IOobject::NO_WRITE, phi, and timeName.

Here is the call graph for this function:

◆ ddt0_()

CrankNicolsonDdtScheme< Type >::template DDt0Field< GeoField > & ddt0_ ( const word name,
const dimensionSet dims 
)

Definition at line 109 of file CrankNicolsonDdtScheme.C.

References Foam::dimTime, mesh, Foam::name(), runTime, and timeName.

Here is the call graph for this function:

◆ fvcDdtUfCorr() [3/4]

tmp< surfaceScalarField > fvcDdtUfCorr ( const GeometricField< scalar, fvPatchField, volMesh > &  U,
const GeometricField< scalar, fvsPatchField, surfaceMesh > &  Uf 
)

◆ fvcDdtPhiCorr() [3/4]

tmp< surfaceScalarField > fvcDdtPhiCorr ( const volScalarField U,
const surfaceScalarField phi 
)

◆ fvcDdtUfCorr() [4/4]

tmp< surfaceScalarField > fvcDdtUfCorr ( const volScalarField rho,
const volScalarField U,
const surfaceScalarField Uf 
)

◆ fvcDdtPhiCorr() [4/4]

tmp< surfaceScalarField > fvcDdtPhiCorr ( const volScalarField rho,
const volScalarField U,
const surfaceScalarField phi 
)

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