00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036 #ifndef ddtScheme_H
00037 #define ddtScheme_H
00038
00039 #include "tmp.H"
00040 #include "dimensionedType.H"
00041 #include "volFieldsFwd.H"
00042 #include "surfaceFieldsFwd.H"
00043 #include "typeInfo.H"
00044 #include "runTimeSelectionTables.H"
00045
00046
00047
00048 namespace Foam
00049 {
00050
00051 template<class Type>
00052 class fvMatrix;
00053
00054 class fvMesh;
00055
00056
00057
00058 namespace fv
00059 {
00060
00061
00062
00063
00064
00065 template<class Type>
00066 class ddtScheme
00067 :
00068 public refCount
00069 {
00070
00071 protected:
00072
00073
00074
00075 const fvMesh& mesh_;
00076
00077
00078
00079
00080
00081 ddtScheme(const ddtScheme&);
00082
00083
00084 void operator=(const ddtScheme&);
00085
00086
00087 public:
00088
00089
00090 virtual const word& type() const = 0;
00091
00092
00093
00094
00095 declareRunTimeSelectionTable
00096 (
00097 tmp,
00098 ddtScheme,
00099 Istream,
00100 (const fvMesh& mesh, Istream& schemeData),
00101 (mesh, schemeData)
00102 );
00103
00104
00105
00106
00107
00108 ddtScheme(const fvMesh& mesh)
00109 :
00110 mesh_(mesh)
00111 {}
00112
00113
00114 ddtScheme(const fvMesh& mesh, Istream&)
00115 :
00116 mesh_(mesh)
00117 {}
00118
00119
00120
00121
00122
00123 static tmp<ddtScheme<Type> > New
00124 (
00125 const fvMesh& mesh,
00126 Istream& schemeData
00127 );
00128
00129
00130
00131
00132 virtual ~ddtScheme();
00133
00134
00135
00136
00137
00138 const fvMesh& mesh() const
00139 {
00140 return mesh_;
00141 }
00142
00143 virtual tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
00144 (
00145 const dimensioned<Type>&
00146 ) = 0;
00147
00148 virtual tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
00149 (
00150 const GeometricField<Type, fvPatchField, volMesh>&
00151 ) = 0;
00152
00153 virtual tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
00154 (
00155 const dimensionedScalar&,
00156 const GeometricField<Type, fvPatchField, volMesh>&
00157 ) = 0;
00158
00159 virtual tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
00160 (
00161 const volScalarField&,
00162 const GeometricField<Type, fvPatchField, volMesh>&
00163 ) = 0;
00164
00165 virtual tmp<fvMatrix<Type> > fvmDdt
00166 (
00167 GeometricField<Type, fvPatchField, volMesh>&
00168 ) = 0;
00169
00170 virtual tmp<fvMatrix<Type> > fvmDdt
00171 (
00172 const dimensionedScalar&,
00173 GeometricField<Type, fvPatchField, volMesh>&
00174 ) = 0;
00175
00176 virtual tmp<fvMatrix<Type> > fvmDdt
00177 (
00178 const volScalarField&,
00179 GeometricField<Type, fvPatchField, volMesh>&
00180 ) = 0;
00181
00182
00183 typedef GeometricField
00184 <
00185 typename flux<Type>::type,
00186 fvPatchField,
00187 surfaceMesh
00188 > fluxFieldType;
00189
00190 tmp<surfaceScalarField> fvcDdtPhiCoeff
00191 (
00192 const GeometricField<Type, fvPatchField, volMesh>& U,
00193 const fluxFieldType& phi,
00194 const fluxFieldType& phiCorr
00195 );
00196
00197 tmp<surfaceScalarField> fvcDdtPhiCoeff
00198 (
00199 const GeometricField<Type, fvPatchField, volMesh>& U,
00200 const fluxFieldType& phi
00201 );
00202
00203 virtual tmp<fluxFieldType> fvcDdtPhiCorr
00204 (
00205 const volScalarField& rA,
00206 const GeometricField<Type, fvPatchField, volMesh>& U,
00207 const fluxFieldType& phi
00208 ) = 0;
00209
00210 tmp<surfaceScalarField> fvcDdtPhiCoeff
00211 (
00212 const volScalarField& rho,
00213 const GeometricField<Type, fvPatchField, volMesh>& rhoU,
00214 const fluxFieldType& phi
00215 );
00216
00217 virtual tmp<fluxFieldType> fvcDdtPhiCorr
00218 (
00219 const volScalarField& rA,
00220 const volScalarField& rho,
00221 const GeometricField<Type, fvPatchField, volMesh>& U,
00222 const fluxFieldType& phi
00223 ) = 0;
00224
00225
00226 virtual tmp<surfaceScalarField> meshPhi
00227 (
00228 const GeometricField<Type, fvPatchField, volMesh>&
00229 ) = 0;
00230 };
00231
00232
00233
00234
00235 }
00236
00237
00238
00239 }
00240
00241
00242
00243
00244
00245 #define makeFvDdtTypeScheme(SS, Type) \
00246 \
00247 defineNamedTemplateTypeNameAndDebug(SS<Type>, 0); \
00248 \
00249 ddtScheme<Type>::addIstreamConstructorToTable<SS<Type> > \
00250 add##SS##Type##IstreamConstructorToTable_;
00251
00252
00253 #define makeFvDdtScheme(SS) \
00254 \
00255 makeFvDdtTypeScheme(SS, scalar) \
00256 makeFvDdtTypeScheme(SS, vector) \
00257 makeFvDdtTypeScheme(SS, tensor) \
00258 makeFvDdtTypeScheme(SS, sphericalTensor) \
00259 \
00260 template<> \
00261 tmp<surfaceScalarField> SS<scalar>::fvcDdtPhiCorr \
00262 ( \
00263 const volScalarField& rA, \
00264 const volScalarField& U, \
00265 const surfaceScalarField& phi \
00266 ) \
00267 { \
00268 notImplemented(#SS"<scalar>::fvcDdtPhiCorr"); \
00269 return surfaceScalarField::null(); \
00270 } \
00271 \
00272 template<> \
00273 tmp<surfaceScalarField> SS<scalar>::fvcDdtPhiCorr \
00274 ( \
00275 const volScalarField& rA, \
00276 const volScalarField& rho, \
00277 const volScalarField& U, \
00278 const surfaceScalarField& phi \
00279 ) \
00280 { \
00281 notImplemented(#SS"<scalar>::fvcDdtPhiCorr"); \
00282 return surfaceScalarField::null(); \
00283 }
00284
00285
00286
00287
00288 #ifdef NoRepository
00289 # include "ddtScheme.C"
00290 #endif
00291
00292
00293
00294 #endif
00295
00296