Home
Downloads
Documentation
Installation
User Guide
man-pages
API Documentation
README
Release Notes
Changes
License
Support
SourceForge Project
Main Page
Related Pages
Namespaces
Classes
Files
File List
File Members
All
Classes
Namespaces
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
Friends
Macros
Pages
src
finiteVolume
finiteVolume
ddtSchemes
steadyStateDdtScheme
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
>
28
#include <
finiteVolume/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> >
44
steadyStateDdtScheme<Type>::fvcDdt
45
(
46
const
dimensioned<Type>
& dt
47
)
48
{
49
return
tmp<GeometricField<Type, fvPatchField, volMesh>
>
50
(
51
new
GeometricField<Type, fvPatchField, volMesh>
52
(
53
IOobject
54
(
55
"ddt("
+dt.
name
()+
')'
,
56
mesh
().time().timeName(),
57
mesh
()
58
),
59
mesh
(),
60
dimensioned<Type>
61
(
62
"0"
,
63
dt.
dimensions
()/
dimTime
,
64
pTraits<Type>::zero
65
)
66
)
67
);
68
}
69
70
71
template
<
class
Type>
72
tmp<GeometricField<Type, fvPatchField, volMesh>
>
73
steadyStateDdtScheme<Type>::fvcDdt
74
(
75
const
GeometricField<Type, fvPatchField, volMesh>
& vf
76
)
77
{
78
return
tmp<GeometricField<Type, fvPatchField, volMesh>
>
79
(
80
new
GeometricField<Type, fvPatchField, volMesh>
81
(
82
IOobject
83
(
84
"ddt("
+vf.
name
()+
')'
,
85
mesh
().time().timeName(),
86
mesh
()
87
),
88
mesh
(),
89
dimensioned<Type>
90
(
91
"0"
,
92
vf.
dimensions
()/
dimTime
,
93
pTraits<Type>::zero
94
)
95
)
96
);
97
}
98
99
100
template
<
class
Type>
101
tmp<GeometricField<Type, fvPatchField, volMesh>
>
102
steadyStateDdtScheme<Type>::fvcDdt
103
(
104
const
dimensionedScalar
&
rho
,
105
const
GeometricField<Type, fvPatchField, volMesh>
& vf
106
)
107
{
108
return
tmp<GeometricField<Type, fvPatchField, volMesh>
>
109
(
110
new
GeometricField<Type, fvPatchField, volMesh>
111
(
112
IOobject
113
(
114
"ddt("
+rho.
name
()+
','
+vf.
name
()+
')'
,
115
mesh
().time().timeName(),
116
mesh
()
117
),
118
mesh
(),
119
dimensioned<Type>
120
(
121
"0"
,
122
rho.
dimensions
()*vf.
dimensions
()/
dimTime
,
123
pTraits<Type>::zero
124
)
125
)
126
);
127
}
128
129
130
template
<
class
Type>
131
tmp<GeometricField<Type, fvPatchField, volMesh>
>
132
steadyStateDdtScheme<Type>::fvcDdt
133
(
134
const
volScalarField
& rho,
135
const
GeometricField<Type, fvPatchField, volMesh>
& vf
136
)
137
{
138
return
tmp<GeometricField<Type, fvPatchField, volMesh>
>
139
(
140
new
GeometricField<Type, fvPatchField, volMesh>
141
(
142
IOobject
143
(
144
"ddt("
+rho.
name
()+
','
+vf.
name
()+
')'
,
145
mesh
().time().timeName(),
146
mesh
()
147
),
148
mesh
(),
149
dimensioned<Type>
150
(
151
"0"
,
152
rho.
dimensions
()*vf.
dimensions
()/
dimTime
,
153
pTraits<Type>::zero
154
)
155
)
156
);
157
}
158
159
160
template
<
class
Type>
161
tmp<fvMatrix<Type>
>
162
steadyStateDdtScheme<Type>::fvmDdt
163
(
164
GeometricField<Type, fvPatchField, volMesh>
& vf
165
)
166
{
167
tmp<fvMatrix<Type>
> tfvm
168
(
169
new
fvMatrix<Type>
170
(
171
vf,
172
vf.
dimensions
()*
dimVol
/
dimTime
173
)
174
);
175
176
return
tfvm;
177
}
178
179
180
template
<
class
Type>
181
tmp<fvMatrix<Type>
>
182
steadyStateDdtScheme<Type>::fvmDdt
183
(
184
const
dimensionedScalar
& rho,
185
GeometricField<Type, fvPatchField, volMesh>
& vf
186
)
187
{
188
tmp<fvMatrix<Type>
> tfvm
189
(
190
new
fvMatrix<Type>
191
(
192
vf,
193
rho.
dimensions
()*vf.
dimensions
()*
dimVol
/
dimTime
194
)
195
);
196
197
return
tfvm;
198
}
199
200
201
template
<
class
Type>
202
tmp<fvMatrix<Type>
>
203
steadyStateDdtScheme<Type>::fvmDdt
204
(
205
const
volScalarField
& rho,
206
GeometricField<Type, fvPatchField, volMesh>
& vf
207
)
208
{
209
tmp<fvMatrix<Type>
> tfvm
210
(
211
new
fvMatrix<Type>
212
(
213
vf,
214
rho.
dimensions
()*vf.
dimensions
()*
dimVol
/
dimTime
215
)
216
);
217
218
return
tfvm;
219
}
220
221
222
template
<
class
Type>
223
tmp<typename steadyStateDdtScheme<Type>::fluxFieldType
>
224
steadyStateDdtScheme<Type>::fvcDdtPhiCorr
225
(
226
const
volScalarField
& rA,
227
const
GeometricField<Type, fvPatchField, volMesh>
&
U
,
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
(),
243
dimensioned<typename flux<Type>::type
>
244
(
245
"0"
,
246
rA.
dimensions
()*phi.
dimensions
()/
dimTime
,
247
pTraits<typename flux<Type>::type
>
::zero
248
)
249
)
250
);
251
}
252
253
254
template
<
class
Type>
255
tmp<typename steadyStateDdtScheme<Type>::fluxFieldType
>
256
steadyStateDdtScheme<Type>::fvcDdtPhiCorr
257
(
258
const
volScalarField
& rA,
259
const
volScalarField
& rho,
260
const
GeometricField<Type, fvPatchField, volMesh>
& U,
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
(),
277
dimensioned<typename flux<Type>::type
>
278
(
279
"0"
,
280
rA.
dimensions
()*rho.
dimensions
()*phi.
dimensions
()/
dimTime
,
281
pTraits<typename flux<Type>::type
>
::zero
282
)
283
)
284
);
285
}
286
287
288
template
<
class
Type>
289
tmp<surfaceScalarField>
steadyStateDdtScheme<Type>::meshPhi
290
(
291
const
GeometricField<Type, fvPatchField, volMesh>
& vf
292
)
293
{
294
return
tmp<surfaceScalarField>
295
(
296
new
surfaceScalarField
297
(
298
IOobject
299
(
300
"meshPhi"
,
301
mesh
().time().
timeName
(),
302
mesh
()
303
),
304
mesh
(),
305
dimensionedScalar
(
"0"
,
dimVolume
/
dimTime
, 0.0)
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: ************************ //