steadyStateDdtScheme.C
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-2018 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 \*---------------------------------------------------------------------------*/
25 
26 #include "steadyStateDdtScheme.H"
27 #include "fvcDiv.H"
28 #include "fvMatrices.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 
35 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36 
37 namespace fv
38 {
39 
40 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
41 
42 template<class Type>
43 tmp<GeometricField<Type, fvPatchField, volMesh>>
45 (
46  const dimensioned<Type>& dt
47 )
48 {
50  (
52  (
53  IOobject
54  (
55  "ddt("+dt.name()+')',
56  mesh().time().timeName(),
57  mesh()
58  ),
59  mesh(),
61  (
62  "0",
63  dt.dimensions()/dimTime,
64  Zero
65  )
66  )
67  );
68 }
69 
70 
71 template<class Type>
74 (
76 )
77 {
79  (
81  (
82  IOobject
83  (
84  "ddt("+vf.name()+')',
85  mesh().time().timeName(),
86  mesh()
87  ),
88  mesh(),
90  (
91  "0",
92  vf.dimensions()/dimTime,
93  Zero
94  )
95  )
96  );
97 }
98 
99 
100 template<class Type>
103 (
104  const dimensionedScalar& rho,
106 )
107 {
109  (
111  (
112  IOobject
113  (
114  "ddt("+rho.name()+','+vf.name()+')',
115  mesh().time().timeName(),
116  mesh()
117  ),
118  mesh(),
120  (
121  "0",
122  rho.dimensions()*vf.dimensions()/dimTime,
123  Zero
124  )
125  )
126  );
127 }
128 
129 
130 template<class Type>
133 (
134  const volScalarField& rho,
136 )
137 {
139  (
141  (
142  IOobject
143  (
144  "ddt("+rho.name()+','+vf.name()+')',
145  mesh().time().timeName(),
146  mesh()
147  ),
148  mesh(),
150  (
151  "0",
152  rho.dimensions()*vf.dimensions()/dimTime,
153  Zero
154  )
155  )
156  );
157 }
158 
159 
160 template<class Type>
163 (
164  const volScalarField& alpha,
165  const volScalarField& rho,
167 )
168 {
170  (
172  (
173  IOobject
174  (
175  "ddt("+alpha.name()+','+rho.name()+','+vf.name()+')',
176  mesh().time().timeName(),
177  mesh()
178  ),
179  mesh(),
181  (
182  "0",
183  rho.dimensions()*vf.dimensions()/dimTime,
184  Zero
185  )
186  )
187  );
188 }
189 
190 
191 template<class Type>
194 (
196 )
197 {
198  tmp<fvMatrix<Type>> tfvm
199  (
200  new fvMatrix<Type>
201  (
202  vf,
204  )
205  );
206 
207  return tfvm;
208 }
209 
210 
211 template<class Type>
214 (
215  const dimensionedScalar& rho,
217 )
218 {
219  tmp<fvMatrix<Type>> tfvm
220  (
221  new fvMatrix<Type>
222  (
223  vf,
225  )
226  );
227 
228  return tfvm;
229 }
230 
231 
232 template<class Type>
235 (
236  const volScalarField& rho,
238 )
239 {
240  tmp<fvMatrix<Type>> tfvm
241  (
242  new fvMatrix<Type>
243  (
244  vf,
246  )
247  );
248 
249  return tfvm;
250 }
251 
252 
253 template<class Type>
256 (
257  const volScalarField& alpha,
258  const volScalarField& rho,
260 )
261 {
262  tmp<fvMatrix<Type>> tfvm
263  (
264  new fvMatrix<Type>
265  (
266  vf,
267  alpha.dimensions()*rho.dimensions()*vf.dimensions()*dimVol/dimTime
268  )
269  );
270 
271  return tfvm;
272 }
273 
274 
275 template<class Type>
278 (
281 )
282 {
283  return tmp<fluxFieldType>
284  (
285  new fluxFieldType
286  (
287  IOobject
288  (
289  "ddtCorr(" + U.name() + ',' + Uf.name() + ')',
290  mesh().time().timeName(),
291  mesh()
292  ),
293  mesh(),
295  (
296  "0",
298  Zero
299  )
300  )
301  );
302 }
303 
304 
305 template<class Type>
308 (
310  const fluxFieldType& phi
311 )
312 {
313  return tmp<fluxFieldType>
314  (
315  new fluxFieldType
316  (
317  IOobject
318  (
319  "ddtCorr(" + U.name() + ',' + phi.name() + ')',
320  mesh().time().timeName(),
321  mesh()
322  ),
323  mesh(),
325  (
326  "0",
327  phi.dimensions()/dimTime,
328  Zero
329  )
330  )
331  );
332 }
333 
334 
335 template<class Type>
338 (
339  const volScalarField& rho,
342 )
343 {
344  return tmp<fluxFieldType>
345  (
346  new fluxFieldType
347  (
348  IOobject
349  (
350  "ddtCorr("
351  + rho.name()
352  + ',' + U.name() + ',' + Uf.name() + ')',
353  mesh().time().timeName(),
354  mesh()
355  ),
356  mesh(),
358  (
359  "0",
361  Zero
362  )
363  )
364  );
365 }
366 
367 
368 template<class Type>
371 (
372  const volScalarField& rho,
374  const fluxFieldType& phi
375 )
376 {
377  return tmp<fluxFieldType>
378  (
379  new fluxFieldType
380  (
381  IOobject
382  (
383  "ddtCorr("
384  + rho.name()
385  + ',' + U.name() + ',' + phi.name() + ')',
386  mesh().time().timeName(),
387  mesh()
388  ),
389  mesh(),
391  (
392  "0",
393  phi.dimensions()/dimTime,
394  Zero
395  )
396  )
397  );
398 }
399 
400 
401 template<class Type>
403 (
405 )
406 {
408  (
409  "meshPhi",
410  mesh(),
412  );
413 }
414 
415 
416 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
417 
418 } // End namespace fv
419 
420 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
421 
422 } // End namespace Foam
423 
424 // ************************************************************************* //
const word & name() const
Return name.
Definition: IOobject.H:295
virtual tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt(const dimensioned< Type > &)
const dimensionSet dimVol(dimVolume)
Definition: dimensionSets.H:59
Generic GeometricField class.
Generic dimensioned Type class.
static tmp< GeometricField< scalar, fvsPatchField, surfaceMesh > > New(const word &name, const Mesh &, const dimensionSet &, const word &patchFieldType=fvsPatchField< scalar >::calculatedType())
Return a temporary field constructed from name, mesh, dimensionSet.
virtual tmp< fluxFieldType > fvcDdtPhiCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const fluxFieldType &phi)
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
const dimensionSet & dimensions() const
Return dimensions.
dynamicFvMesh & mesh
labelList fv(nPoints)
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
static const zero Zero
Definition: zero.H:97
Calculate the divergence of the given field.
const word & name() const
Return const reference to name.
virtual tmp< fluxFieldType > fvcDdtUfCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)
virtual tmp< surfaceScalarField > meshPhi(const GeometricField< Type, fvPatchField, volMesh > &)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
A special matrix type and solver, designed for finite volume solutions of scalar equations.
const dimensionSet & dimensions() const
Return const reference to dimensions.
virtual tmp< fvMatrix< Type > > fvmDdt(const GeometricField< Type, fvPatchField, volMesh > &)
A class for managing temporary objects.
Definition: PtrList.H:53
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:57
Namespace for OpenFOAM.