localEulerDdtScheme.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2023 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 face time-step field
73  static word rDeltaTfName;
74 
75  //- Name of the reciprocal local sub-cycling time-step field
76  static word rSubDeltaTName;
77 
78 
79  // Constructors
80 
82  {}
83 
84 
85  // Member Functions
86 
87  //- Return true if LTS is enabled
88  static bool enabled(const fvMesh& mesh);
89 
90  //- Return the reciprocal of the local time-step
91  // looked-up from the objectRegistry
92  static const volScalarField& localRDeltaT(const fvMesh& mesh);
93 
94  //- Return the reciprocal of the local face time-step
95  // looked-up from the objectRegistry
96  static const surfaceScalarField& localRDeltaTf(const fvMesh& mesh);
97 
98  //- Calculate and return the reciprocal of the local sub-cycling
99  // time-step
101  (
102  const fvMesh& mesh,
103  const label nAlphaSubCycles
104  );
105 };
106 
107 
108 /*---------------------------------------------------------------------------*\
109  Class localEulerDdtScheme Declaration
110 \*---------------------------------------------------------------------------*/
111 
112 template<class Type>
114 :
115  public localEulerDdt,
116  public fv::ddtScheme<Type>
117 {
118  // Private Member Functions
119 
120  //- Return the reciprocal of the local time-step
121  const volScalarField& localRDeltaT() const;
122 
123  //- Return the reciprocal of the local face time-step
124  const surfaceScalarField& localRDeltaTf() const;
125 
126 
127 public:
128 
129  //- Runtime type information
130  TypeName("localEuler");
131 
132 
133  // Constructors
134 
135  //- Construct from mesh
137  :
138  ddtScheme<Type>(mesh)
139  {}
140 
141  //- Construct from mesh and Istream
143  :
144  ddtScheme<Type>(mesh, is)
145  {}
146 
147  //- Disallow default bitwise copy construction
148  localEulerDdtScheme(const localEulerDdtScheme&) = delete;
149 
150 
151  // Member Functions
152 
153  using ddtScheme<Type>::mesh;
154 
155  virtual tmp<VolField<Type>> fvcDdt
156  (
157  const dimensioned<Type>&
158  );
159 
160  virtual tmp<VolField<Type>> fvcDdt
161  (
162  const VolField<Type>&
163  );
164 
165  virtual tmp<VolField<Type>> fvcDdt
166  (
167  const dimensionedScalar&,
168  const VolField<Type>&
169  );
170 
171  virtual tmp<VolField<Type>> fvcDdt
172  (
173  const volScalarField&,
174  const VolField<Type>&
175  );
176 
177  virtual tmp<VolField<Type>> fvcDdt
178  (
179  const volScalarField& alpha,
180  const volScalarField& rho,
181  const VolField<Type>& psi
182  );
183 
185  (
186  const SurfaceField<Type>&
187  );
188 
189  virtual tmp<fvMatrix<Type>> fvmDdt
190  (
191  const VolField<Type>&
192  );
193 
194  virtual tmp<fvMatrix<Type>> fvmDdt
195  (
196  const dimensionedScalar&,
197  const VolField<Type>&
198  );
199 
200  virtual tmp<fvMatrix<Type>> fvmDdt
201  (
202  const volScalarField&,
203  const VolField<Type>&
204  );
205 
206  virtual tmp<fvMatrix<Type>> fvmDdt
207  (
208  const volScalarField& alpha,
209  const volScalarField& rho,
210  const VolField<Type>& psi
211  );
212 
214 
216  (
217  const VolField<Type>& U,
218  const SurfaceField<Type>& Uf
219  );
220 
222  (
223  const VolField<Type>& U,
224  const fluxFieldType& phi
225  );
226 
228  (
229  const volScalarField& rho,
230  const VolField<Type>& U,
231  const SurfaceField<Type>& rhoUf
232  );
233 
235  (
236  const volScalarField& rho,
237  const VolField<Type>& U,
238  const fluxFieldType& phi
239  );
240 
242  (
243  const volScalarField& alpha,
244  const volScalarField& rho,
245  const VolField<Type>& U,
246  const SurfaceField<Type>& Uf
247  );
248 
250  (
251  const volScalarField& alpha,
252  const volScalarField& rho,
253  const VolField<Type>& U,
254  const fluxFieldType& phi
255  );
256 
258  (
259  const VolField<Type>&
260  );
261 
262  virtual tmp<scalarField> meshPhi
263  (
264  const VolField<Type>&,
265  const label patchi
266  );
267 
268 
269  // Member Operators
270 
271  //- Disallow default bitwise assignment
272  void operator=(const localEulerDdtScheme&) = delete;
273 };
274 
275 
276 template<>
278 (
279  const VolField<scalar>& U,
280  const SurfaceField<scalar>& Uf
281 );
282 
283 template<>
285 (
286  const volScalarField& U,
287  const surfaceScalarField& phi
288 );
289 
290 template<>
292 (
293  const volScalarField& rho,
294  const volScalarField& U,
295  const surfaceScalarField& Uf
296 );
297 
298 template<>
300 (
301  const volScalarField& rho,
302  const volScalarField& U,
303  const surfaceScalarField& phi
304 );
305 
306 template<>
308 (
309  const volScalarField& alpha,
310  const volScalarField& rho,
311  const volScalarField& U,
312  const surfaceScalarField& Uf
313 );
314 
315 template<>
317 (
318  const volScalarField& alpha,
319  const volScalarField& rho,
320  const volScalarField& U,
321  const surfaceScalarField& phi
322 );
323 
324 template<>
326 (
327  const volScalarField& alpha,
328  const volScalarField& rho,
329  const volScalarField& U,
330  const surfaceScalarField& Uf
331 );
332 
333 template<>
335 (
336  const volScalarField& alpha,
337  const volScalarField& rho,
338  const volScalarField& U,
339  const surfaceScalarField& phi
340 );
341 
342 
343 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
344 
345 } // End namespace fv
346 
347 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
348 
349 } // End namespace Foam
350 
351 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
352 
353 #ifdef NoRepository
354  #include "localEulerDdtScheme.C"
355 #endif
356 
357 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
358 
359 #endif
360 
361 // ************************************************************************* //
const label nAlphaSubCycles(alphaControls.lookup< label >("nAlphaSubCycles"))
Generic GeometricField class.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:60
Generic dimensioned Type class.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
Abstract base class for ddt schemes.
Definition: ddtScheme.H:68
const fvMesh & mesh() const
Return mesh reference.
Definition: ddtScheme.H:130
Local time-step first-order Euler implicit/explicit ddt.
virtual tmp< fluxFieldType > fvcDdtUfCorr(const VolField< Type > &U, const SurfaceField< Type > &Uf)
localEulerDdtScheme(const fvMesh &mesh)
Construct from mesh.
void operator=(const localEulerDdtScheme &)=delete
Disallow default bitwise assignment.
virtual tmp< fvMatrix< Type > > fvmDdt(const VolField< Type > &)
TypeName("localEuler")
Runtime type information.
virtual tmp< surfaceScalarField > meshPhi(const VolField< Type > &)
ddtScheme< Type >::fluxFieldType fluxFieldType
virtual tmp< fluxFieldType > fvcDdtPhiCorr(const VolField< Type > &U, const fluxFieldType &phi)
virtual tmp< VolField< Type > > fvcDdt(const dimensioned< Type > &)
static word rSubDeltaTName
Name of the reciprocal local sub-cycling time-step field.
static const surfaceScalarField & localRDeltaTf(const fvMesh &mesh)
Return the reciprocal of the local face time-step.
Definition: localEulerDdt.C:58
static word rDeltaTfName
Name of the reciprocal local face time-step field.
static tmp< volScalarField > localRSubDeltaT(const fvMesh &mesh, const label nAlphaSubCycles)
Calculate and return the reciprocal of the local sub-cycling.
Definition: localEulerDdt.C:70
static word rDeltaTName
Name of the reciprocal local time-step field.
static const volScalarField & localRDeltaT(const fvMesh &mesh)
Return the reciprocal of the local time-step.
Definition: localEulerDdt.C:46
static bool enabled(const fvMesh &mesh)
Return true if LTS is enabled.
Definition: localEulerDdt.C:37
A class for managing temporary objects.
Definition: tmp.H:55
A class for handling words, derived from string.
Definition: word.H:62
label patchi
const volScalarField & psi
U
Definition: pEqn.H:72
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
Namespace for OpenFOAM.
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
labelList fv(nPoints)