CoEulerDdtScheme.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 Class
25  Foam::fv::CoEulerDdtScheme
26 
27 Description
28  Courant number limited first-order Euler implicit/explicit ddt.
29 
30  The time-step is adjusted locally so that the local Courant number
31  does not exceed the specified value.
32 
33  This scheme should only be used for steady-state computations
34  using transient codes where local time-stepping is preferable to
35  under-relaxation for transport consistency reasons.
36 
37 SourceFiles
38  CoEulerDdtScheme.C
39 
40 \*---------------------------------------------------------------------------*/
41 
42 #ifndef CoEulerDdtScheme_H
43 #define CoEulerDdtScheme_H
44 
45 #include "ddtScheme.H"
46 
47 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
48 
49 namespace Foam
50 {
51 
52 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
53 
54 namespace fv
55 {
56 
57 /*---------------------------------------------------------------------------*\
58  Class CoEulerDdtScheme Declaration
59 \*---------------------------------------------------------------------------*/
60 
61 template<class Type>
62 class CoEulerDdtScheme
63 :
64  public fv::ddtScheme<Type>
65 {
66  // Private Data
67 
68  //- Name of the flux field used to calculate the local time-step
69  word phiName_;
70 
71  //- Name of the density field used to obtain the volumetric flux
72  // from the mass flux if required
73  word rhoName_;
74 
75  //- Maximum local Courant number
76  scalar maxCo_;
77 
78 
79  // Private Member Functions
80 
81  //- Disallow default bitwise copy construct
83 
84  //- Disallow default bitwise assignment
85  void operator=(const CoEulerDdtScheme&);
86 
87  //- Return the reciprocal of the Courant-number limited time-step
88  tmp<volScalarField> CorDeltaT() const;
89 
90  //- Return the reciprocal of the face-Courant-number limited time-step
91  tmp<surfaceScalarField> CofrDeltaT() const;
92 
93 
94 public:
95 
96  //- Runtime type information
97  TypeName("CoEuler");
98 
99 
100  // Constructors
101 
102  //- Construct from mesh and Istream
103  CoEulerDdtScheme(const fvMesh& mesh, Istream& is)
104  :
105  ddtScheme<Type>(mesh, is),
106  phiName_(is),
107  rhoName_(is),
108  maxCo_(readScalar(is))
109  {}
110 
111 
112  // Member Functions
113 
114  //- Return mesh reference
115  const fvMesh& mesh() const
116  {
117  return fv::ddtScheme<Type>::mesh();
118  }
119 
121  (
122  const dimensioned<Type>&
123  );
124 
126  (
128  );
129 
131  (
132  const dimensionedScalar&,
134  );
135 
137  (
138  const volScalarField&,
140  );
141 
143  (
144  const volScalarField& alpha,
145  const volScalarField& rho,
147  );
148 
150  (
152  );
153 
155  (
156  const dimensionedScalar&,
158  );
159 
161  (
162  const volScalarField&,
164  );
165 
167  (
168  const volScalarField& alpha,
169  const volScalarField& rho,
171  );
174 
176  (
179  );
180 
182  (
184  const fluxFieldType& phi
185  );
186 
188  (
189  const volScalarField& rho,
192  );
193 
195  (
196  const volScalarField& rho,
198  const fluxFieldType& phi
199  );
200 
202  (
204  );
205 };
206 
207 
208 template<>
210 (
213 );
214 
215 template<>
217 (
218  const volScalarField& U,
219  const surfaceScalarField& phi
220 );
221 
222 template<>
224 (
225  const volScalarField& rho,
226  const volScalarField& U,
227  const surfaceScalarField& Uf
228 );
229 
230 template<>
232 (
233  const volScalarField& rho,
234  const volScalarField& U,
235  const surfaceScalarField& phi
236 );
237 
238 
239 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
240 
241 } // End namespace fv
242 
243 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
244 
245 } // End namespace Foam
246 
247 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
248 
249 #ifdef NoRepository
250  #include "CoEulerDdtScheme.C"
251 #endif
252 
253 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
254 
255 #endif
256 
257 // ************************************************************************* //
surfaceScalarField & phi
const fvMesh & mesh() const
Return mesh reference.
tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt(const dimensioned< Type > &)
U
Definition: pEqn.H:83
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
Uf
Definition: pEqn.H:67
Generic GeometricField class.
Generic dimensioned Type class.
TypeName("CoEuler")
Runtime type information.
Abstract base class for ddt schemes.
Definition: ddtScheme.H:64
const fvMesh & mesh() const
Return mesh reference.
Definition: ddtScheme.H:135
A class for handling words, derived from string.
Definition: word.H:59
labelList fv(nPoints)
tmp< surfaceScalarField > meshPhi(const GeometricField< Type, fvPatchField, volMesh > &)
tmp< fluxFieldType > fvcDdtUfCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)
Courant number limited first-order Euler implicit/explicit ddt.
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
tmp< fvMatrix< Type > > fvmDdt(const GeometricField< Type, fvPatchField, volMesh > &)
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
ddtScheme< Type >::fluxFieldType fluxFieldType
const volScalarField & psi
tmp< fluxFieldType > fvcDdtPhiCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const fluxFieldType &phi)
A class for managing temporary objects.
Definition: PtrList.H:53
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Namespace for OpenFOAM.