FreeFOAM The Cross-Platform CFD Toolkit
steadyStateDdtScheme.C
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) 1991-2010 OpenCFD Ltd.
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 <finiteVolume/fvcDiv.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,
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,
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,
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,
154  )
155  )
156  );
157 }
158 
159 
160 template<class Type>
163 (
165 )
166 {
167  tmp<fvMatrix<Type> > tfvm
168  (
169  new fvMatrix<Type>
170  (
171  vf,
173  )
174  );
175 
176  return tfvm;
177 }
178 
179 
180 template<class Type>
183 (
184  const dimensionedScalar& rho,
186 )
187 {
188  tmp<fvMatrix<Type> > tfvm
189  (
190  new fvMatrix<Type>
191  (
192  vf,
194  )
195  );
196 
197  return tfvm;
198 }
199 
200 
201 template<class Type>
204 (
205  const volScalarField& rho,
207 )
208 {
209  tmp<fvMatrix<Type> > tfvm
210  (
211  new fvMatrix<Type>
212  (
213  vf,
215  )
216  );
217 
218  return tfvm;
219 }
220 
221 
222 template<class Type>
225 (
226  const volScalarField& rA,
228  const fluxFieldType& phi
229 )
230 {
231  return tmp<fluxFieldType>
232  (
233  new fluxFieldType
234  (
235  IOobject
236  (
237  "ddtPhiCorr("
238  + rA.name() + ',' + U.name() + ',' + phi.name() + ')',
239  mesh().time().timeName(),
240  mesh()
241  ),
242  mesh(),
244  (
245  "0",
246  rA.dimensions()*phi.dimensions()/dimTime,
248  )
249  )
250  );
251 }
252 
253 
254 template<class Type>
257 (
258  const volScalarField& rA,
259  const volScalarField& rho,
261  const fluxFieldType& phi
262 )
263 {
264  return tmp<fluxFieldType>
265  (
266  new fluxFieldType
267  (
268  IOobject
269  (
270  "ddtPhiCorr("
271  + rA.name() + ',' + rho.name()
272  + ',' + U.name() + ',' + phi.name() + ')',
273  mesh().time().timeName(),
274  mesh()
275  ),
276  mesh(),
278  (
279  "0",
280  rA.dimensions()*rho.dimensions()*phi.dimensions()/dimTime,
282  )
283  )
284  );
285 }
286 
287 
288 template<class Type>
290 (
292 )
293 {
295  (
297  (
298  IOobject
299  (
300  "meshPhi",
301  mesh().time().timeName(),
302  mesh()
303  ),
304  mesh(),
306  )
307  );
308 }
309 
310 
311 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
312 
313 } // End namespace fv
314 
315 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
316 
317 } // End namespace Foam
318 
319 // ************************ vim: set sw=4 sts=4 et: ************************ //