OpenFOAM logo
Open Source CFD Toolkit

CrankNicholsonDdtScheme.H

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------*\
00002   =========                 |
00003   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
00004    \\    /   O peration     |
00005     \\  /    A nd           | Copyright (C) 1991-2005 OpenCFD Ltd.
00006      \\/     M anipulation  |
00007 -------------------------------------------------------------------------------
00008 License
00009     This file is part of OpenFOAM.
00010 
00011     OpenFOAM is free software; you can redistribute it and/or modify it
00012     under the terms of the GNU General Public License as published by the
00013     Free Software Foundation; either version 2 of the License, or (at your
00014     option) any later version.
00015 
00016     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
00017     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
00019     for more details.
00020 
00021     You should have received a copy of the GNU General Public License
00022     along with OpenFOAM; if not, write to the Free Software Foundation,
00023     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
00024 
00025 Class
00026     CrankNicholsonDdtScheme
00027 
00028 Description
00029     Second-oder CrankNicholson implicit ddt using the current and
00030     two previous time-step as well as the previous time-step ddt.
00031 
00032 SourceFiles
00033     CrankNicholsonDdtScheme.C
00034 
00035 \*---------------------------------------------------------------------------*/
00036 
00037 #ifndef CrankNicholsonDdtScheme_H
00038 #define CrankNicholsonDdtScheme_H
00039 
00040 #include "ddtScheme.H"
00041 
00042 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00043 
00044 namespace Foam
00045 {
00046 
00047 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00048 
00049 namespace fv
00050 {
00051 
00052 /*---------------------------------------------------------------------------*\
00053                        Class CrankNicholsonDdtScheme Declaration
00054 \*---------------------------------------------------------------------------*/
00055 
00056 template<class Type>
00057 class CrankNicholsonDdtScheme
00058 :
00059     public fv::ddtScheme<Type>
00060 {
00061     // Private Data
00062 
00063         //- Off-centering coefficient, 1 -> CN, less than one blends with EI
00064         scalar ocCoeff_;
00065 
00066 
00067     // Private Member Functions
00068 
00069         //- Disallow default bitwise copy construct
00070         CrankNicholsonDdtScheme(const CrankNicholsonDdtScheme&);
00071 
00072         //- Disallow default bitwise assignment
00073         void operator=(const CrankNicholsonDdtScheme&);
00074 
00075         template<class GeoField>
00076         GeoField& ddt0_(const word& name, const dimensionSet& dims);
00077 
00078         template<class GeoField>
00079         bool evaluate(const GeoField& ddt0) const;
00080 
00081         template<class GeoField>
00082         dimensionedScalar rDtCoef_(const GeoField& ddt0) const;
00083 
00084         template<class GeoField>
00085         dimensionedScalar rDtCoef0_(const GeoField& ddt0) const;
00086 
00087         template<class GeoField>
00088         tmp<GeoField> offCentre_(const GeoField& ddt0) const;
00089 
00090 
00091 public:
00092 
00093     //- Runtime type information
00094     TypeName("CrankNicholson");
00095 
00096 
00097     // Constructors
00098 
00099         //- Construct from mesh
00100         CrankNicholsonDdtScheme(const fvMesh& mesh)
00101         :
00102             ddtScheme<Type>(mesh),
00103             ocCoeff_(1.0)
00104         {}
00105 
00106         //- Construct from mesh and Istream
00107         CrankNicholsonDdtScheme(const fvMesh& mesh, Istream& is)
00108         :
00109             ddtScheme<Type>(mesh, is),
00110             ocCoeff_(readScalar(is))
00111         {
00112             if (ocCoeff_ < 0 || ocCoeff_ > 1)
00113             {
00114                 FatalIOErrorIn
00115                 (
00116                     "CrankNicholsonDdtScheme(const fvMesh& mesh, Istream& is)",
00117                     is
00118                 )   << "coefficient = " << ocCoeff_
00119                     << " should be >= 0 and <= 1"
00120                     << exit(FatalIOError);
00121             }
00122         }
00123 
00124 
00125     // Member Functions
00126 
00127         //- Return mesh reference
00128         const fvMesh& mesh() const
00129         {
00130             return fv::ddtScheme<Type>::mesh();
00131         }
00132 
00133         tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
00134         (
00135             const dimensioned<Type>&
00136         );
00137 
00138         tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
00139         (
00140             const GeometricField<Type, fvPatchField, volMesh>&
00141         );
00142 
00143         tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
00144         (
00145             const dimensionedScalar&,
00146             const GeometricField<Type, fvPatchField, volMesh>&
00147         );
00148 
00149         tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
00150         (
00151             const volScalarField&,
00152             const GeometricField<Type, fvPatchField, volMesh>&
00153         );
00154 
00155         tmp<fvMatrix<Type> > fvmDdt
00156         (
00157             GeometricField<Type, fvPatchField, volMesh>&
00158         );
00159 
00160         tmp<fvMatrix<Type> > fvmDdt
00161         (
00162             const dimensionedScalar&,
00163             GeometricField<Type, fvPatchField, volMesh>&
00164         );
00165 
00166         tmp<fvMatrix<Type> > fvmDdt
00167         (
00168             const volScalarField&,
00169             GeometricField<Type, fvPatchField, volMesh>&
00170         );
00171 
00172         typedef typename ddtScheme<Type>::fluxFieldType fluxFieldType;
00173 
00174         tmp<fluxFieldType> fvcDdtPhiCorr
00175         (
00176             const volScalarField& rA,
00177             const GeometricField<Type, fvPatchField, volMesh>& U,
00178             const fluxFieldType& phi
00179         );
00180 
00181         tmp<fluxFieldType> fvcDdtPhiCorr
00182         (
00183             const volScalarField& rA,
00184             const volScalarField& rho,
00185             const GeometricField<Type, fvPatchField, volMesh>& U,
00186             const fluxFieldType& phi
00187         );
00188 
00189 
00190         tmp<surfaceScalarField> meshPhi
00191         (
00192             const GeometricField<Type, fvPatchField, volMesh>&
00193         );
00194 };
00195 
00196 
00197 template<>
00198 tmp<surfaceScalarField> CrankNicholsonDdtScheme<scalar>::fvcDdtPhiCorr
00199 (
00200     const volScalarField& rA,
00201     const volScalarField& U,
00202     const surfaceScalarField& phi
00203 );
00204 
00205 
00206 template<>
00207 tmp<surfaceScalarField> CrankNicholsonDdtScheme<scalar>::fvcDdtPhiCorr
00208 (
00209     const volScalarField& rA,
00210     const volScalarField& rho,
00211     const volScalarField& U,
00212     const surfaceScalarField& phi
00213 );
00214 
00215 
00216 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00217 
00218 } // End namespace fv
00219 
00220 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00221 
00222 } // End namespace Foam
00223 
00224 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00225 
00226 #ifdef NoRepository
00227 #   include "CrankNicholsonDdtScheme.C"
00228 #endif
00229 
00230 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00231 
00232 #endif
00233 
00234 // ************************************************************************* //
For further information go to www.openfoam.org