FreeFOAM The Cross-Platform CFD Toolkit
KinematicCloudI_.H
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 <finiteVolume/fvmSup.H>
27 
28 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
29 
30 template<class ParcelType>
32 {
33  return parcelTypeId_;
34 }
35 
36 
37 template<class ParcelType>
39 {
40  return mesh_;
41 }
42 
43 
44 template<class ParcelType>
45 inline const Foam::IOdictionary&
47 {
48  return particleProperties_;
49 }
50 
51 
52 template<class ParcelType>
53 inline const typename ParcelType::constantProperties&
55 {
56  return constProps_;
57 }
58 
59 
60 template<class ParcelType>
62 {
63  return active_;
64 }
65 
66 
67 template<class ParcelType>
69 {
70  return coupled_;
71 }
72 
73 
74 template <class ParcelType>
75 inline const Foam::Switch
77 {
78  return cellValueSourceCorrection_;
79 }
80 
81 
82 template<class ParcelType>
83 inline const Foam::volScalarField&
85 {
86  return rho_;
87 }
88 
89 
90 template<class ParcelType>
92 {
93  return U_;
94 }
95 
96 
97 template<class ParcelType>
99 {
100  return mu_;
101 }
102 
103 
104 template<class ParcelType>
105 inline const Foam::dimensionedVector&
107 {
108  return g_;
109 }
110 
111 
112 template<class ParcelType>
113 inline const Foam::particleForces&
115 {
116  return forces_;
117 }
118 
119 
120 template<class ParcelType>
121 inline const Foam::dictionary&
123 {
124  return interpolationSchemes_;
125 }
126 
127 
128 template<class ParcelType>
131 {
132  return dispersionModel_;
133 }
134 
135 
136 template<class ParcelType>
139 {
140  return dispersionModel_();
141 }
142 
143 
144 template<class ParcelType>
147 {
148  return dragModel_;
149 }
150 
151 
152 template<class ParcelType>
155 {
156  return injectionModel_;
157 }
158 
159 
160 template<class ParcelType>
163 {
164  return patchInteractionModel_;
165 }
166 
167 
168 template<class ParcelType>
171 {
172  return injectionModel_();
173 }
174 
175 
176 template<class ParcelType>
179 {
180  return postProcessingModel_();
181 }
182 
183 
184 template<class ParcelType>
185 inline const Foam::vectorIntegrationScheme&
187 {
188  return UIntegrator_;
189 }
190 
191 
192 template<class ParcelType>
194 {
195  scalar sysMass = 0.0;
196  forAllConstIter(typename KinematicCloud<ParcelType>, *this, iter)
197  {
198  const ParcelType& p = iter();
199  sysMass += p.mass()*p.nParticle();
200  }
201 
202  return sysMass;
203 }
204 
205 
206 template<class ParcelType>
208 {
209  return rndGen_;
210 }
211 
212 
213 template<class ParcelType>
216 {
217  return UTrans_;
218 }
219 
220 
221 template<class ParcelType>
224 {
226  (
228  (
229  IOobject
230  (
231  this->name() + "SU",
232  this->db().time().timeName(),
233  this->mesh(),
234  IOobject::NO_READ,
235  IOobject::AUTO_WRITE
236  ),
237  this->mesh(),
239  (
240  "zero",
242  vector::zero
243  )
244  )
245  );
246 
247  vectorField& SU = tSU().field();
248  SU = UTrans_/(mesh_.V()*this->db().time().deltaT());
249 
250  return tSU;
251 }
252 
253 
254 template<class ParcelType>
257 {
258  tmp<volScalarField> ttheta
259  (
260  new volScalarField
261  (
262  IOobject
263  (
264  this->name() + "Theta",
265  this->db().time().timeName(),
266  this->db(),
267  IOobject::NO_READ,
268  IOobject::NO_WRITE,
269  false
270  ),
271  mesh_,
272  dimensionedScalar("zero", dimless, 0.0)
273  )
274  );
275 
276  scalarField& theta = ttheta().internalField();
277  forAllConstIter(typename KinematicCloud<ParcelType>, *this, iter)
278  {
279  const ParcelType& p = iter();
280  const label cellI = p.cell();
281 
282  theta[cellI] += p.nParticle()*p.volume();
283  }
284 
285  theta /= mesh().V();
286 
287  return ttheta;
288 }
289 
290 
291 template<class ParcelType>
294 {
295  tmp<volScalarField> talpha
296  (
297  new volScalarField
298  (
299  IOobject
300  (
301  this->name() + "Alpha",
302  this->db().time().timeName(),
303  this->db(),
304  IOobject::NO_READ,
305  IOobject::NO_WRITE,
306  false
307  ),
308  mesh_,
309  dimensionedScalar("zero", dimless, 0.0)
310  )
311  );
312 
313  scalarField& alpha = talpha().internalField();
314  forAllConstIter(typename KinematicCloud<ParcelType>, *this, iter)
315  {
316  const ParcelType& p = iter();
317  const label cellI = p.cell();
318 
319  alpha[cellI] += p.nParticle()*p.mass();
320  }
321 
322  alpha /= (mesh().V()*rho_);
323 
324  return talpha;
325 }
326 
327 
328 template<class ParcelType>
331 {
332  tmp<volScalarField> trhoEff
333  (
334  new volScalarField
335  (
336  IOobject
337  (
338  this->name() + "RhoEff",
339  this->db().time().timeName(),
340  this->db(),
341  IOobject::NO_READ,
342  IOobject::NO_WRITE,
343  false
344  ),
345  mesh_,
346  dimensionedScalar("zero", dimDensity, 0.0)
347  )
348  );
349 
350  scalarField& rhoEff = trhoEff().internalField();
351  forAllConstIter(typename KinematicCloud<ParcelType>, *this, iter)
352  {
353  const ParcelType& p = iter();
354  const label cellI = p.cell();
355 
356  rhoEff[cellI] += p.nParticle()*p.mass();
357  }
358 
359  rhoEff /= mesh().V();
360 
361  return trhoEff;
362 }
363 
364 
365 // ************************ vim: set sw=4 sts=4 et: ************************ //