Second-oder Crank-Nicolson implicit ddt using the current and previous time-step fields as well as the previous time-step ddt. More...
Public Types | |
typedef ddtScheme< Type >::fluxFieldType | fluxFieldType |
Public Types inherited from ddtScheme< Type > | |
typedef SurfaceField< typename flux< Type >::type > | fluxFieldType |
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... | |
CrankNicolsonDdtScheme (const CrankNicolsonDdtScheme &)=delete | |
Disallow default bitwise copy construction. More... | |
const fvMesh & | mesh () const |
Return mesh reference. More... | |
scalar | ocCoeff () const |
Return the current off-centering coefficient. More... | |
virtual tmp< VolField< Type > > | fvcDdt (const dimensioned< Type > &) |
virtual tmp< VolField< Type > > | fvcDdt (const VolField< Type > &) |
virtual tmp< VolField< Type > > | fvcDdt (const dimensionedScalar &, const VolField< Type > &) |
virtual tmp< VolField< Type > > | fvcDdt (const volScalarField &, const VolField< Type > &) |
virtual tmp< VolField< Type > > | fvcDdt (const volScalarField &alpha, const volScalarField &rho, const VolField< Type > &psi) |
virtual tmp< fvMatrix< Type > > | fvmDdt (const VolField< Type > &) |
virtual tmp< fvMatrix< Type > > | fvmDdt (const dimensionedScalar &, const VolField< Type > &) |
virtual tmp< fvMatrix< Type > > | fvmDdt (const volScalarField &, const VolField< Type > &) |
virtual tmp< fvMatrix< Type > > | fvmDdt (const volScalarField &alpha, const volScalarField &rho, const VolField< Type > &psi) |
virtual tmp< fluxFieldType > | fvcDdtUfCorr (const VolField< Type > &U, const SurfaceField< Type > &Uf) |
virtual tmp< fluxFieldType > | fvcDdtPhiCorr (const VolField< Type > &U, const fluxFieldType &phi) |
virtual tmp< fluxFieldType > | fvcDdtUfCorr (const volScalarField &rho, const VolField< Type > &U, const SurfaceField< Type > &rhoUf) |
virtual tmp< fluxFieldType > | fvcDdtPhiCorr (const volScalarField &rho, const VolField< Type > &U, const fluxFieldType &phi) |
virtual tmp< fluxFieldType > | fvcDdtUfCorr (const volScalarField &alpha, const volScalarField &rho, const VolField< Type > &U, const SurfaceField< Type > &Uf) |
virtual tmp< fluxFieldType > | fvcDdtPhiCorr (const volScalarField &alpha, const volScalarField &rho, const VolField< Type > &U, const fluxFieldType &phi) |
virtual tmp< surfaceScalarField > | meshPhi (const VolField< Type > &) |
virtual tmp< scalarField > | meshPhi (const VolField< Type > &, const label patchi) |
void | operator= (const CrankNicolsonDdtScheme &)=delete |
Disallow default bitwise assignment. More... | |
template<class GeoField > | |
CrankNicolsonDdtScheme< Type >::template DDt0Field< GeoField > & | ddt0_ (const word &unTypedName, const dimensionSet &dims) |
tmp< surfaceScalarField > | fvcDdtUfCorr (const VolField< scalar > &U, const SurfaceField< scalar > &Uf) |
tmp< surfaceScalarField > | fvcDdtPhiCorr (const volScalarField &U, const surfaceScalarField &phi) |
tmp< surfaceScalarField > | fvcDdtUfCorr (const volScalarField &rho, const volScalarField &U, const surfaceScalarField &Uf) |
tmp< surfaceScalarField > | fvcDdtPhiCorr (const volScalarField &rho, const volScalarField &U, const surfaceScalarField &phi) |
tmp< surfaceScalarField > | fvcDdtUfCorr (const volScalarField &alpha, const volScalarField &rho, const volScalarField &U, const surfaceScalarField &Uf) |
tmp< surfaceScalarField > | fvcDdtPhiCorr (const volScalarField &alpha, const volScalarField &rho, const volScalarField &U, const surfaceScalarField &phi) |
Public Member Functions inherited from ddtScheme< Type > | |
virtual const word & | type () 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 &) | |
Construct from mesh and Istream. More... | |
ddtScheme (const ddtScheme &)=delete | |
Disallow default bitwise copy construction. More... | |
virtual | ~ddtScheme () |
Destructor. More... | |
const fvMesh & | mesh () const |
Return mesh reference. More... | |
virtual tmp< SurfaceField< Type > > | fvcDdt (const SurfaceField< Type > &) |
virtual tmp< surfaceScalarField > | fvcDdtPhiCoeff (const VolField< Type > &U, const fluxFieldType &phi, const fluxFieldType &phiCorr) |
virtual tmp< surfaceScalarField > | fvcDdtPhiCoeff (const VolField< Type > &U, const fluxFieldType &phi, const fluxFieldType &phiCorr, const volScalarField &rho) |
virtual tmp< surfaceScalarField > | fvcDdtPhiCoeff (const VolField< Type > &U, const fluxFieldType &phi) |
virtual tmp< surfaceScalarField > | fvcDdtPhiCoeff (const VolField< Type > &U, const fluxFieldType &phi, const volScalarField &rho) |
void | operator= (const ddtScheme &)=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... | |
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... | |
Protected Member Functions inherited from refCount | |
refCount () | |
Construct null initialising count to 0. More... | |
refCount (const refCount &)=delete | |
Disallow copy. More... | |
Protected Attributes inherited from ddtScheme< Type > | |
const fvMesh & | mesh_ |
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 stabilise 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);
Definition at line 110 of file CrankNicolsonDdtScheme.H.
typedef ddtScheme<Type>::fluxFieldType fluxFieldType |
Definition at line 283 of file CrankNicolsonDdtScheme.H.
CrankNicolsonDdtScheme | ( | const fvMesh & | mesh | ) |
Construct from mesh.
Definition at line 283 of file CrankNicolsonDdtScheme.C.
References CrankNicolsonDdtScheme< Type >::mesh(), polyMesh::moving(), and fvMesh::V00().
CrankNicolsonDdtScheme | ( | const fvMesh & | mesh, |
Istream & | is | ||
) |
Construct from mesh and Istream.
Definition at line 298 of file CrankNicolsonDdtScheme.C.
References dict, Foam::exit(), Foam::FatalIOError, FatalIOErrorInFunction, token::isNumber(), CrankNicolsonDdtScheme< Type >::mesh(), polyMesh::moving(), Function1< Type >::New(), token::number(), CrankNicolsonDdtScheme< Type >::ocCoeff(), Istream::putBack(), fvMesh::time(), Foam::unitFraction, Time::userUnits(), and fvMesh::V00().
|
delete |
Disallow default bitwise copy construction.
TypeName | ( | "CrankNicolson" | ) |
Runtime type information.
|
inline |
Return mesh reference.
Definition at line 219 of file CrankNicolsonDdtScheme.H.
References ddtScheme< Type >::mesh().
Referenced by CrankNicolsonDdtScheme< Type >::CrankNicolsonDdtScheme(), and CrankNicolsonDdtScheme< Type >::ocCoeff().
|
inline |
Return the current off-centering coefficient.
Definition at line 225 of file CrankNicolsonDdtScheme.H.
References CrankNicolsonDdtScheme< Type >::mesh().
Referenced by CrankNicolsonDdtScheme< Type >::CrankNicolsonDdtScheme().
|
virtual |
Implements ddtScheme< Type >.
Definition at line 350 of file CrankNicolsonDdtScheme.C.
References dimensioned< Type >::dimensions(), Foam::dimTime, Foam::evaluate(), dimensioned< Type >::name(), tmp< T >::ref(), and Foam::Zero.
Implements ddtScheme< Type >.
Definition at line 407 of file CrankNicolsonDdtScheme.C.
References GeometricField< Type, PatchField, GeoMesh >::boundaryField(), DimensionedField< Type, GeoMesh >::dimensions(), Foam::evaluate(), Foam::fv::ff(), IOobject::name(), GeometricField< Type, PatchField, GeoMesh >::New(), OldTimeField< FieldType >::oldTime(), and dimensioned< Type >::value().
|
virtual |
Implements ddtScheme< Type >.
Definition at line 483 of file CrankNicolsonDdtScheme.C.
References GeometricField< Type, PatchField, GeoMesh >::boundaryField(), DimensionedField< Type, GeoMesh >::dimensions(), Foam::evaluate(), Foam::fv::ff(), IOobject::name(), GeometricField< Type, PatchField, GeoMesh >::New(), OldTimeField< FieldType >::oldTime(), rho, and dimensioned< Type >::value().
|
virtual |
Implements ddtScheme< Type >.
Definition at line 560 of file CrankNicolsonDdtScheme.C.
References GeometricField< Type, PatchField, GeoMesh >::boundaryField(), DimensionedField< Type, GeoMesh >::dimensions(), Foam::evaluate(), Foam::fv::ff(), IOobject::name(), GeometricField< Type, PatchField, GeoMesh >::New(), OldTimeField< FieldType >::oldTime(), rho, and dimensioned< Type >::value().
|
virtual |
Implements ddtScheme< Type >.
Definition at line 646 of file CrankNicolsonDdtScheme.C.
References alpha(), GeometricField< Type, PatchField, GeoMesh >::boundaryField(), DimensionedField< Type, GeoMesh >::dimensions(), Foam::evaluate(), Foam::fv::ff(), IOobject::name(), GeometricField< Type, PatchField, GeoMesh >::New(), OldTimeField< FieldType >::oldTime(), rho, and dimensioned< Type >::value().
Implements ddtScheme< Type >.
Definition at line 758 of file CrankNicolsonDdtScheme.C.
References lduMatrix::diag(), DimensionedField< Type, GeoMesh >::dimensions(), Foam::dimTime, Foam::dimVolume, Foam::evaluate(), Foam::fv::ff(), IOobject::name(), OldTimeField< FieldType >::oldTime(), and fvMatrix< Type >::source().
|
virtual |
Implements ddtScheme< Type >.
Definition at line 840 of file CrankNicolsonDdtScheme.C.
References lduMatrix::diag(), DimensionedField< Type, GeoMesh >::dimensions(), Foam::dimTime, Foam::dimVolume, Foam::evaluate(), Foam::fv::ff(), IOobject::name(), OldTimeField< FieldType >::oldTime(), rho, and fvMatrix< Type >::source().
|
virtual |
Implements ddtScheme< Type >.
Definition at line 922 of file CrankNicolsonDdtScheme.C.
References lduMatrix::diag(), DimensionedField< Type, GeoMesh >::dimensions(), Foam::dimTime, Foam::dimVolume, Foam::evaluate(), Foam::fv::ff(), IOobject::name(), OldTimeField< FieldType >::oldTime(), rho, and fvMatrix< Type >::source().
|
virtual |
Implements ddtScheme< Type >.
Definition at line 1012 of file CrankNicolsonDdtScheme.C.
References alpha(), lduMatrix::diag(), DimensionedField< Type, GeoMesh >::dimensions(), Foam::dimTime, Foam::dimVolume, Foam::evaluate(), Foam::fv::ff(), IOobject::name(), OldTimeField< FieldType >::oldTime(), rho, and fvMatrix< Type >::source().
|
virtual |
Implements ddtScheme< Type >.
Definition at line 1127 of file CrankNicolsonDdtScheme.C.
References DimensionedField< Type, GeoMesh >::dimensions(), Foam::evaluate(), Foam::fvc::interpolate(), IOobject::name(), Foam::compressible::New(), OldTimeField< FieldType >::oldTime(), and U.
|
virtual |
Implements ddtScheme< Type >.
Definition at line 1180 of file CrankNicolsonDdtScheme.C.
References DimensionedField< Type, GeoMesh >::dimensions(), Foam::fvc::dotInterpolate(), Foam::evaluate(), IOobject::name(), Foam::compressible::New(), OldTimeField< FieldType >::oldTime(), and U.
|
virtual |
Implements ddtScheme< Type >.
Definition at line 1234 of file CrankNicolsonDdtScheme.C.
References Foam::abort(), DimensionedField< Type, GeoMesh >::dimensions(), Foam::dimVelocity, Foam::evaluate(), Foam::FatalError, FatalErrorInFunction, Foam::fvc::interpolate(), IOobject::name(), Foam::compressible::New(), OldTimeField< FieldType >::oldTime(), rho, and U.
|
virtual |
Implements ddtScheme< Type >.
Definition at line 1373 of file CrankNicolsonDdtScheme.C.
References Foam::abort(), DimensionedField< Type, GeoMesh >::dimensions(), Foam::dimVelocity, Foam::dimVolumetricFlux, Foam::fvc::dotInterpolate(), Foam::evaluate(), Foam::FatalError, FatalErrorInFunction, IOobject::name(), Foam::compressible::New(), OldTimeField< FieldType >::oldTime(), rho, and U.
|
inlinevirtual |
Implements ddtScheme< Type >.
Definition at line 311 of file CrankNicolsonDdtScheme.H.
References NotImplemented, and GeometricField< Type, PatchField, GeoMesh >::null().
|
inlinevirtual |
Implements ddtScheme< Type >.
Definition at line 323 of file CrankNicolsonDdtScheme.H.
References NotImplemented, and GeometricField< Type, PatchField, GeoMesh >::null().
|
virtual |
Implements ddtScheme< Type >.
Definition at line 1499 of file CrankNicolsonDdtScheme.C.
References Foam::dimVolume, Foam::evaluate(), Foam::name(), and GeometricField< Type, PatchField, GeoMesh >::New().
|
virtual |
Implements ddtScheme< Type >.
Definition at line 1522 of file CrankNicolsonDdtScheme.C.
References Foam::dimVolume, Foam::evaluate(), and patchi.
|
delete |
Disallow default bitwise assignment.
CrankNicolsonDdtScheme<Type>::template DDt0Field<GeoField>& ddt0_ | ( | const word & | unTypedName, |
const dimensionSet & | dims | ||
) |
Definition at line 106 of file CrankNicolsonDdtScheme.C.
References Foam::dimTime, Foam::name(), Time::startTime(), Time::startTimeIndex(), objectRegistry::time(), TimeState::timeIndex(), Time::timeName(), Foam::typedName(), dimensioned< Type >::value(), and Foam::Zero.
tmp< surfaceScalarField > fvcDdtUfCorr | ( | const VolField< scalar > & | U, |
const SurfaceField< scalar > & | Uf | ||
) |
tmp< surfaceScalarField > fvcDdtPhiCorr | ( | const volScalarField & | U, |
const surfaceScalarField & | phi | ||
) |
tmp< surfaceScalarField > fvcDdtUfCorr | ( | const volScalarField & | rho, |
const volScalarField & | U, | ||
const surfaceScalarField & | Uf | ||
) |
tmp< surfaceScalarField > fvcDdtPhiCorr | ( | const volScalarField & | rho, |
const volScalarField & | U, | ||
const surfaceScalarField & | phi | ||
) |
tmp< surfaceScalarField > fvcDdtUfCorr | ( | const volScalarField & | alpha, |
const volScalarField & | rho, | ||
const volScalarField & | U, | ||
const surfaceScalarField & | Uf | ||
) |
tmp< surfaceScalarField > fvcDdtPhiCorr | ( | const volScalarField & | alpha, |
const volScalarField & | rho, | ||
const volScalarField & | U, | ||
const surfaceScalarField & | phi | ||
) |