localEulerDdtScheme.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::localEulerDdtScheme
26 
27 Description
28  Local time-step first-order Euler implicit/explicit ddt.
29 
30  The reciprocal of the local time-step field is looked-up from the database.
31 
32  This scheme should only be used for steady-state computations using
33  transient codes where local time-stepping is preferably to under-relaxation
34  for transport consistency reasons.
35 
36 See also
37  Foam::fv::CoEulerDdtScheme
38 
39 SourceFiles
40  localEulerDdt.C
41  localEulerDdtScheme.C
42  localEulerDdtSchemes.C
43 
44 \*---------------------------------------------------------------------------*/
45 
46 #ifndef localEulerDdtScheme_H
47 #define localEulerDdtScheme_H
48 
49 #include "ddtScheme.H"
50 
51 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
52 
53 namespace Foam
54 {
55 
56 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
57 
58 namespace fv
59 {
60 
61 /*---------------------------------------------------------------------------*\
62  Class localEulerDdt Declaration
63 \*---------------------------------------------------------------------------*/
64 
65 class localEulerDdt
66 {
67 public:
68 
69  //- Name of the reciprocal local time-step field
70  static word rDeltaTName;
71 
72  //- Name of the reciprocal local sub-cycling time-step field
73  static word rSubDeltaTName;
74 
75 
76  // Constructors
77 
79  {}
80 
81 
82  // Member Functions
83 
84  //- Return true if LTS is enabled
85  static bool enabled(const fvMesh& mesh);
86 
87  //- Return the reciprocal of the local time-step
88  // looked-up from the objectRegistry
89  static const volScalarField& localRDeltaT(const fvMesh& mesh);
90 
91  //- Calculate and return the reciprocal of the local sub-cycling
92  // time-step
94  (
95  const fvMesh& mesh,
97  );
98 };
99 
100 
101 /*---------------------------------------------------------------------------*\
102  Class localEulerDdtScheme Declaration
103 \*---------------------------------------------------------------------------*/
104 
105 template<class Type>
107 :
108  public localEulerDdt,
109  public fv::ddtScheme<Type>
110 {
111  // Private Member Functions
112 
113  //- Disallow default bitwise copy construct
115 
116  //- Disallow default bitwise assignment
117  void operator=(const localEulerDdtScheme&);
118 
119  //- Return the reciprocal of the local time-step
120  const volScalarField& localRDeltaT() const;
121 
122 
123 public:
124 
125  //- Runtime type information
126  TypeName("localEuler");
127 
128 
129  // Constructors
130 
131  //- Construct from mesh
133  :
134  ddtScheme<Type>(mesh)
135  {}
136 
137  //- Construct from mesh and Istream
139  :
140  ddtScheme<Type>(mesh, is)
141  {}
142 
143 
144  // Member Functions
145 
146  //- Return mesh reference
147  const fvMesh& mesh() const
148  {
149  return fv::ddtScheme<Type>::mesh();
150  }
151 
153  (
154  const dimensioned<Type>&
155  );
156 
158  (
160  );
161 
163  (
164  const dimensionedScalar&,
166  );
167 
169  (
170  const volScalarField&,
172  );
173 
175  (
176  const volScalarField& alpha,
177  const volScalarField& rho,
179  );
180 
181  tmp<fvMatrix<Type>> fvmDdt
182  (
184  );
185 
186  tmp<fvMatrix<Type>> fvmDdt
187  (
188  const dimensionedScalar&,
190  );
191 
192  tmp<fvMatrix<Type>> fvmDdt
193  (
194  const volScalarField&,
196  );
197 
198  tmp<fvMatrix<Type>> fvmDdt
199  (
200  const volScalarField& alpha,
201  const volScalarField& rho,
203  );
206 
207  tmp<fluxFieldType> fvcDdtUfCorr
208  (
211  );
212 
213  tmp<fluxFieldType> fvcDdtPhiCorr
214  (
216  const fluxFieldType& phi
217  );
218 
219  tmp<fluxFieldType> fvcDdtUfCorr
220  (
221  const volScalarField& rho,
224  );
225 
226  tmp<fluxFieldType> fvcDdtPhiCorr
227  (
228  const volScalarField& rho,
230  const fluxFieldType& phi
231  );
232 
234  (
236  );
237 };
238 
239 
240 template<>
242 (
245 );
246 
247 template<>
249 (
250  const volScalarField& U,
251  const surfaceScalarField& phi
252 );
253 
254 template<>
256 (
257  const volScalarField& rho,
258  const volScalarField& U,
259  const surfaceScalarField& Uf
260 );
261 
262 template<>
264 (
265  const volScalarField& rho,
266  const volScalarField& U,
267  const surfaceScalarField& phi
268 );
269 
270 
271 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
272 
273 } // End namespace fv
274 
275 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
276 
277 } // End namespace Foam
278 
279 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
280 
281 #ifdef NoRepository
282  #include "localEulerDdtScheme.C"
283 #endif
284 
285 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
286 
287 #endif
288 
289 // ************************************************************************* //
Local time-step first-order Euler implicit/explicit ddt.
surfaceScalarField & phi
U
Definition: pEqn.H:83
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
const fvMesh & mesh() const
Return mesh reference.
Definition: ddtScheme.H:135
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
static tmp< volScalarField > localRSubDeltaT(const fvMesh &mesh, const label nAlphaSubCycles)
Calculate and return the reciprocal of the local sub-cycling.
Definition: localEulerDdt.C:57
tmp< fluxFieldType > fvcDdtUfCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)
Uf
Definition: pEqn.H:67
Generic dimensioned Type class.
Abstract base class for ddt schemes.
Definition: ddtScheme.H:64
dynamicFvMesh & mesh
static bool enabled(const fvMesh &mesh)
Return true if LTS is enabled.
Definition: localEulerDdt.C:36
A class for handling words, derived from string.
Definition: word.H:59
labelList fv(nPoints)
tmp< fluxFieldType > fvcDdtPhiCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const fluxFieldType &phi)
static const volScalarField & localRDeltaT(const fvMesh &mesh)
Return the reciprocal of the local time-step.
Definition: localEulerDdt.C:45
#define TypeName(TypeNameString)
Declare a ClassName() with extra virtual type info.
Definition: typeInfo.H:70
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
tmp< surfaceScalarField > meshPhi(const volVectorField &U)
Definition: fvcMeshPhi.C:34
const volScalarField & psi
A class for managing temporary objects.
Definition: PtrList.H:54
static word rDeltaTName
Name of the reciprocal local time-step field.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
static word rSubDeltaTName
Name of the reciprocal local sub-cycling time-step field.
Namespace for OpenFOAM.
label nAlphaSubCycles(readLabel(alphaControls.lookup("nAlphaSubCycles")))