All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
SLTSDdtScheme.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::SLTSDdtScheme
26 
27 Description
28  Stabilised local time-step first-order Euler implicit/explicit ddt.
29  The time-step is adjusted locally so that an advective equations remains
30  diagonally dominant.
31 
32  This scheme should only be used for steady-state computations
33  using transient codes where local time-stepping is preferably to
34  under-relaxation for transport consistency reasons.
35 
36 See also
37  Foam::fv::CoEulerDdtScheme
38 
39 SourceFiles
40  SLTSDdtScheme.C
41 
42 \*---------------------------------------------------------------------------*/
43 
44 #ifndef SLTSDdtScheme_H
45 #define SLTSDdtScheme_H
46 
47 #include "ddtScheme.H"
48 
49 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
50 
51 namespace Foam
52 {
53 
54 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
55 
56 namespace fv
57 {
58 
59 /*---------------------------------------------------------------------------*\
60  Class SLTSDdtScheme Declaration
61 \*---------------------------------------------------------------------------*/
62 
63 template<class Type>
64 class SLTSDdtScheme
65 :
66  public fv::ddtScheme<Type>
67 {
68  // Private Data
69 
70  //- Name of the flux field used to calculate the local time-step
71  word phiName_;
72 
73  //- Name of the density field used to obtain the volumetric flux
74  // from the mass flux if required
75  word rhoName_;
76 
77  //- Under-relaxation factor
78  scalar alpha_;
79 
80 
81  // Private Member Functions
82 
83  //- Calculate a relaxed diagonal from the given flux field
84  void relaxedDiag(scalarField& rD, const surfaceScalarField& phi) const;
85 
86  //- Return the reciprocal of the stabilised local time-step
87  tmp<volScalarField> SLrDeltaT() const;
88 
89 
90 public:
91 
92  //- Runtime type information
93  TypeName("SLTS");
94 
95 
96  // Constructors
97 
98  //- Construct from mesh and Istream
99  SLTSDdtScheme(const fvMesh& mesh, Istream& is)
100  :
101  ddtScheme<Type>(mesh, is),
102  phiName_(is),
103  rhoName_(is),
104  alpha_(readScalar(is))
105  {}
106 
107  //- Disallow default bitwise copy construction
108  SLTSDdtScheme(const SLTSDdtScheme&) = delete;
109 
110 
111  // Member Functions
112 
113  //- Return mesh reference
114  const fvMesh& mesh() const
115  {
116  return fv::ddtScheme<Type>::mesh();
117  }
118 
119  virtual tmp<VolField<Type>> fvcDdt
120  (
121  const dimensioned<Type>&
122  );
123 
124  virtual tmp<VolField<Type>> fvcDdt
125  (
126  const VolField<Type>&
127  );
128 
129  virtual tmp<VolField<Type>> fvcDdt
130  (
131  const dimensionedScalar&,
132  const VolField<Type>&
133  );
134 
135  virtual tmp<VolField<Type>> fvcDdt
136  (
137  const volScalarField&,
138  const VolField<Type>&
139  );
140 
141  virtual tmp<VolField<Type>> fvcDdt
142  (
143  const volScalarField& alpha,
144  const volScalarField& rho,
145  const VolField<Type>& psi
146  );
147 
148  virtual tmp<fvMatrix<Type>> fvmDdt
149  (
150  const VolField<Type>&
151  );
152 
153  virtual tmp<fvMatrix<Type>> fvmDdt
154  (
155  const dimensionedScalar&,
156  const VolField<Type>&
157  );
158 
159  virtual tmp<fvMatrix<Type>> fvmDdt
160  (
161  const volScalarField&,
162  const VolField<Type>&
163  );
164 
165  virtual tmp<fvMatrix<Type>> fvmDdt
166  (
167  const volScalarField& alpha,
168  const volScalarField& rho,
169  const VolField<Type>& psi
170  );
171 
173 
175  (
176  const VolField<Type>& U,
177  const SurfaceField<Type>& Uf
178  );
179 
181  (
182  const VolField<Type>& U,
183  const fluxFieldType& phi
184  );
185 
187  (
188  const volScalarField& rho,
189  const VolField<Type>& U,
190  const SurfaceField<Type>& rhoUf
191  );
192 
194  (
195  const volScalarField& rho,
196  const VolField<Type>& U,
197  const fluxFieldType& phi
198  );
199 
201  (
202  const volScalarField& alpha,
203  const volScalarField& rho,
204  const VolField<Type>& U,
205  const SurfaceField<Type>& Uf
206  )
207  {
209  return fluxFieldType::null();
210  }
211 
213  (
214  const volScalarField& alpha,
215  const volScalarField& rho,
216  const VolField<Type>& U,
217  const fluxFieldType& phi
218  )
219  {
221  return fluxFieldType::null();
222  }
223 
225  (
226  const VolField<Type>&
227  );
228 
229  virtual tmp<scalarField> meshPhi
230  (
231  const VolField<Type>&,
232  const label patchi
233  );
234 
235 
236  // Member Operators
237 
238  //- Disallow default bitwise assignment
239  void operator=(const SLTSDdtScheme&) = delete;
240 };
241 
242 
243 template<>
245 (
246  const VolField<scalar>& U,
247  const SurfaceField<scalar>& Uf
248 );
249 
250 template<>
252 (
253  const volScalarField& U,
254  const surfaceScalarField& phi
255 );
256 
257 template<>
259 (
260  const volScalarField& rho,
261  const volScalarField& U,
262  const surfaceScalarField& Uf
263 );
264 
265 template<>
267 (
268  const volScalarField& rho,
269  const volScalarField& U,
270  const surfaceScalarField& phi
271 );
272 
273 template<>
275 (
276  const volScalarField& alpha,
277  const volScalarField& rho,
278  const volScalarField& U,
279  const surfaceScalarField& Uf
280 );
281 
282 template<>
284 (
285  const volScalarField& alpha,
286  const volScalarField& rho,
287  const volScalarField& U,
288  const surfaceScalarField& phi
289 );
290 
291 
292 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
293 
294 } // End namespace fv
295 
296 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
297 
298 } // End namespace Foam
299 
300 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
301 
302 #ifdef NoRepository
303  #include "SLTSDdtScheme.C"
304 #endif
305 
306 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
307 
308 #endif
309 
310 // ************************************************************************* //
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
Stabilised local time-step first-order Euler implicit/explicit ddt. The time-step is adjusted locally...
Definition: SLTSDdtScheme.H:66
virtual tmp< fluxFieldType > fvcDdtUfCorr(const VolField< Type > &U, const SurfaceField< Type > &Uf)
virtual tmp< fvMatrix< Type > > fvmDdt(const VolField< Type > &)
void operator=(const SLTSDdtScheme &)=delete
Disallow default bitwise assignment.
const fvMesh & mesh() const
Return mesh reference.
SLTSDdtScheme(const fvMesh &mesh, Istream &is)
Construct from mesh and Istream.
Definition: SLTSDdtScheme.H:98
TypeName("SLTS")
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 > &)
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)