All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
CoEulerDdtScheme.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::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  //- Return the reciprocal of the Courant-number limited time-step
82  tmp<volScalarField> CorDeltaT() const;
83 
84  //- Return the reciprocal of the face-Courant-number limited time-step
85  tmp<surfaceScalarField> CofrDeltaT() const;
86 
87 
88 public:
89 
90  //- Runtime type information
91  TypeName("CoEuler");
92 
93 
94  // Constructors
95 
96  //- Construct from mesh and Istream
97  CoEulerDdtScheme(const fvMesh& mesh, Istream& is)
98  :
99  ddtScheme<Type>(mesh, is),
100  phiName_(is),
101  rhoName_(is),
102  maxCo_(readScalar(is))
103  {}
104 
105  //- Disallow default bitwise copy construction
106  CoEulerDdtScheme(const CoEulerDdtScheme&) = delete;
107 
108 
109  // Member Functions
110 
111  //- Return mesh reference
112  const fvMesh& mesh() const
113  {
114  return fv::ddtScheme<Type>::mesh();
115  }
116 
117  virtual tmp<VolField<Type>> fvcDdt
118  (
119  const dimensioned<Type>&
120  );
121 
122  virtual tmp<VolField<Type>> fvcDdt
123  (
124  const VolField<Type>&
125  );
126 
127  virtual tmp<VolField<Type>> fvcDdt
128  (
129  const dimensionedScalar&,
130  const VolField<Type>&
131  );
132 
133  virtual tmp<VolField<Type>> fvcDdt
134  (
135  const volScalarField&,
136  const VolField<Type>&
137  );
138 
139  virtual tmp<VolField<Type>> fvcDdt
140  (
141  const volScalarField& alpha,
142  const volScalarField& rho,
143  const VolField<Type>& psi
144  );
145 
146  virtual tmp<fvMatrix<Type>> fvmDdt
147  (
148  const VolField<Type>&
149  );
150 
151  virtual tmp<fvMatrix<Type>> fvmDdt
152  (
153  const dimensionedScalar&,
154  const VolField<Type>&
155  );
156 
157  virtual tmp<fvMatrix<Type>> fvmDdt
158  (
159  const volScalarField&,
160  const VolField<Type>&
161  );
162 
163  virtual tmp<fvMatrix<Type>> fvmDdt
164  (
165  const volScalarField& alpha,
166  const volScalarField& rho,
167  const VolField<Type>& psi
168  );
169 
171 
173  (
174  const VolField<Type>& U,
175  const SurfaceField<Type>& Uf
176  );
177 
179  (
180  const VolField<Type>& U,
181  const fluxFieldType& phi
182  );
183 
185  (
186  const volScalarField& rho,
187  const VolField<Type>& U,
188  const SurfaceField<Type>& rhoUf
189  );
190 
192  (
193  const volScalarField& rho,
194  const VolField<Type>& U,
195  const fluxFieldType& phi
196  );
197 
199  (
200  const volScalarField& alpha,
201  const volScalarField& rho,
202  const VolField<Type>& U,
203  const SurfaceField<Type>& Uf
204  )
205  {
207  return fluxFieldType::null();
208  }
209 
211  (
212  const volScalarField& alpha,
213  const volScalarField& rho,
214  const VolField<Type>& U,
215  const fluxFieldType& phi
216  )
217  {
219  return fluxFieldType::null();
220  }
221 
223  (
224  const VolField<Type>&
225  );
226 
227  virtual tmp<scalarField> meshPhi
228  (
229  const VolField<Type>&,
230  const label patchi
231  );
232 
233 
234  // Member Operators
235 
236  //- Disallow default bitwise assignment
237  void operator=(const CoEulerDdtScheme&) = delete;
238 };
239 
240 
241 template<>
243 (
244  const VolField<scalar>& U,
245  const SurfaceField<scalar>& Uf
246 );
247 
248 template<>
250 (
251  const volScalarField& U,
252  const surfaceScalarField& phi
253 );
254 
255 template<>
257 (
258  const volScalarField& rho,
259  const volScalarField& U,
260  const surfaceScalarField& Uf
261 );
262 
263 template<>
265 (
266  const volScalarField& rho,
267  const volScalarField& U,
268  const surfaceScalarField& phi
269 );
270 
271 template<>
273 (
274  const volScalarField& alpha,
275  const volScalarField& rho,
276  const volScalarField& U,
277  const surfaceScalarField& Uf
278 );
279 
280 template<>
282 (
283  const volScalarField& alpha,
284  const volScalarField& rho,
285  const volScalarField& U,
286  const surfaceScalarField& phi
287 );
288 
289 
290 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
291 
292 } // End namespace fv
293 
294 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
295 
296 } // End namespace Foam
297 
298 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
299 
300 #ifdef NoRepository
301  #include "CoEulerDdtScheme.C"
302 #endif
303 
304 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
305 
306 #endif
307 
308 // ************************************************************************* //
Generic GeometricField class.
static const GeometricField< Type, PatchField, GeoMesh > & null()
Return a null geometric field.
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
Courant number limited first-order Euler implicit/explicit ddt.
virtual tmp< fluxFieldType > fvcDdtUfCorr(const VolField< Type > &U, const SurfaceField< Type > &Uf)
virtual tmp< fvMatrix< Type > > fvmDdt(const VolField< Type > &)
TypeName("CoEuler")
Runtime type information.
const fvMesh & mesh() const
Return mesh reference.
virtual tmp< surfaceScalarField > meshPhi(const VolField< Type > &)
ddtScheme< Type >::fluxFieldType fluxFieldType
void operator=(const CoEulerDdtScheme &)=delete
Disallow default bitwise assignment.
CoEulerDdtScheme(const fvMesh &mesh, Istream &is)
Construct from mesh and Istream.
virtual tmp< VolField< Type > > fvcDdt(const dimensioned< Type > &)
virtual tmp< fluxFieldType > fvcDdtPhiCorr(const VolField< Type > &U, const fluxFieldType &phi)
Abstract base class for ddt schemes.
Definition: ddtScheme.H:68
const fvMesh & mesh() const
Return mesh reference.
Definition: ddtScheme.H:130
A class for managing temporary objects.
Definition: tmp.H:55
A class for handling words, derived from string.
Definition: word.H:62
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:381
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
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if successful.
Definition: doubleScalar.H:75
labelList fv(nPoints)