ThermoCloudI.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) 2011-2013 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 
27 
28 using namespace Foam::constant;
29 
30 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
31 
32 template<class CloudType>
33 inline const Foam::ThermoCloud<CloudType>&
35 {
36  return cloudCopyPtr_();
37 }
38 
39 
40 template<class CloudType>
41 inline const typename CloudType::particleType::constantProperties&
43 {
44  return constProps_;
45 }
46 
47 
48 template<class CloudType>
49 inline typename CloudType::particleType::constantProperties&
51 {
52  return constProps_;
53 }
54 
55 
56 template<class CloudType>
58 {
59  return thermo_;
60 }
61 
62 
63 template<class CloudType>
65 {
66  return T_;
67 }
68 
69 
70 template<class CloudType>
72 {
73  return p_;
74 }
75 
76 
77 template<class CloudType>
80 {
81  return heatTransferModel_;
82 }
83 
84 
85 template<class CloudType>
88 {
89  return TIntegrator_;
90 }
91 
92 
93 template<class CloudType>
95 {
96  return radiation_;
97 }
98 
99 
100 template<class CloudType>
103 {
104  if (!radiation_)
105  {
107  (
108  "inline Foam::DimensionedField<Foam::scalar, Foam::volMesh> "
109  "Foam::ThermoCloud<CloudType>::radAreaP()"
110  ) << "Radiation field requested, but radiation model not active"
111  << abort(FatalError);
112  }
113 
114  return radAreaP_();
115 }
116 
117 
118 template<class CloudType>
121 {
122  if (!radiation_)
123  {
125  (
126  "inline Foam::DimensionedField<Foam::scalar, Foam::volMesh> "
127  "Foam::ThermoCloud<CloudType>::radAreaP()"
128  ) << "Radiation field requested, but radiation model not active"
129  << abort(FatalError);
130  }
131 
132  return radAreaP_();
133 }
134 
135 
136 template<class CloudType>
139 {
140  if (!radiation_)
141  {
143  (
144  "inline Foam::DimensionedField<Foam::scalar, Foam::volMesh> "
145  "Foam::ThermoCloud<CloudType>::radT4()"
146  ) << "Radiation field requested, but radiation model not active"
147  << abort(FatalError);
148  }
149 
150  return radT4_();
151 }
152 
153 
154 template<class CloudType>
157 {
158  if (!radiation_)
159  {
161  (
162  "inline Foam::DimensionedField<Foam::scalar, Foam::volMesh> "
163  "Foam::ThermoCloud<CloudType>::radT4()"
164  ) << "Radiation field requested, but radiation model not active"
165  << abort(FatalError);
166  }
167 
168  return radT4_();
169 }
170 
171 
172 template<class CloudType>
175 {
176  if (!radiation_)
177  {
179  (
180  "inline Foam::DimensionedField<Foam::scalar, Foam::volMesh> "
181  "Foam::ThermoCloud<CloudType>::radAreaPT4()"
182  ) << "Radiation field requested, but radiation model not active"
183  << abort(FatalError);
184  }
185 
186  return radAreaPT4_();
187 }
188 
189 
190 template<class CloudType>
193 {
194  if (!radiation_)
195  {
197  (
198  "inline Foam::DimensionedField<Foam::scalar, Foam::volMesh> "
199  "Foam::ThermoCloud<CloudType>::radAreaPT4()"
200  ) << "Radiation field requested, but radiation model not active"
201  << abort(FatalError);
202  }
203 
204  return radAreaPT4_();
205 }
206 
207 
208 template<class CloudType>
211 {
212  return hsTrans_();
213 }
214 
215 
216 template<class CloudType>
219 {
220  return hsTrans_();
221 }
222 
223 
224 template<class CloudType>
227 {
228  return hsCoeff_();
229 }
230 
231 
232 template<class CloudType>
235 {
236  return hsCoeff_();
237 }
238 
239 
240 template<class CloudType>
243 {
244  if (debug)
245  {
246  Info<< "hsTrans min/max = " << min(hsTrans()).value() << ", "
247  << max(hsTrans()).value() << nl
248  << "hsCoeff min/max = " << min(hsCoeff()).value() << ", "
249  << max(hsCoeff()).value() << endl;
250  }
251 
252  if (this->solution().coupled())
253  {
254  if (this->solution().semiImplicit("h"))
255  {
256  const volScalarField Cp(thermo_.thermo().Cp());
258  Vdt(this->mesh().V()*this->db().time().deltaT());
259 
260  return
261  hsTrans()/Vdt
262  - fvm::SuSp(hsCoeff()/(Cp*Vdt), hs)
263  + hsCoeff()/(Cp*Vdt)*hs;
264  }
265  else
266  {
268  fvScalarMatrix& fvm = tfvm();
269 
270  fvm.source() = -hsTrans()/(this->db().time().deltaT());
271 
272  return tfvm;
273  }
274  }
275 
277 }
278 
279 
280 template<class CloudType>
282 {
284  (
285  new volScalarField
286  (
287  IOobject
288  (
289  this->name() + ":radiation:Ep",
290  this->db().time().timeName(),
291  this->db(),
294  false
295  ),
296  this->mesh(),
298  )
299  );
300 
301  if (radiation_)
302  {
303  scalarField& Ep = tEp().internalField();
304  const scalar dt = this->db().time().deltaTValue();
305  const scalarField& V = this->mesh().V();
306  const scalar epsilon = constProps_.epsilon0();
307  const scalarField& sumAreaPT4 = radAreaPT4_->field();
308 
309  Ep = sumAreaPT4*epsilon*physicoChemical::sigma.value()/V/dt;
310  }
311 
312  return tEp;
313 }
314 
315 
316 template<class CloudType>
318 {
320  (
321  new volScalarField
322  (
323  IOobject
324  (
325  this->name() + ":radiation:ap",
326  this->db().time().timeName(),
327  this->db(),
330  false
331  ),
332  this->mesh(),
333  dimensionedScalar("zero", dimless/dimLength, 0.0)
334  )
335  );
336 
337  if (radiation_)
338  {
339  scalarField& ap = tap().internalField();
340  const scalar dt = this->db().time().deltaTValue();
341  const scalarField& V = this->mesh().V();
342  const scalar epsilon = constProps_.epsilon0();
343  const scalarField& sumAreaP = radAreaP_->field();
344 
345  ap = sumAreaP*epsilon/V/dt;
346  }
347 
348  return tap;
349 }
350 
351 
352 template<class CloudType>
355 {
356  tmp<volScalarField> tsigmap
357  (
358  new volScalarField
359  (
360  IOobject
361  (
362  this->name() + ":radiation:sigmap",
363  this->db().time().timeName(),
364  this->db(),
367  false
368  ),
369  this->mesh(),
370  dimensionedScalar("zero", dimless/dimLength, 0.0)
371  )
372  );
373 
374  if (radiation_)
375  {
376  scalarField& sigmap = tsigmap().internalField();
377  const scalar dt = this->db().time().deltaTValue();
378  const scalarField& V = this->mesh().V();
379  const scalar epsilon = constProps_.epsilon0();
380  const scalar f = constProps_.f0();
381  const scalarField& sumAreaP = radAreaP_->field();
382 
383  sigmap *= sumAreaP*(1.0 - f)*(1.0 - epsilon)/V/dt;
384  }
385 
386  return tsigmap;
387 }
388 
389 
390 template<class CloudType>
391 inline Foam::scalar Foam::ThermoCloud<CloudType>::Tmax() const
392 {
393  scalar T = -GREAT;
394  scalar n = 0;
395  forAllConstIter(typename ThermoCloud<CloudType>, *this, iter)
396  {
397  const parcelType& p = iter();
398  T = max(T, p.T());
399  n++;
400  }
401 
402  reduce(T, maxOp<scalar>());
403  reduce(n, sumOp<label>());
404 
405  if (n > 0)
406  {
407  return T;
408  }
409  else
410  {
411  return 0.0;
412  }
413 }
414 
415 
416 template<class CloudType>
417 inline Foam::scalar Foam::ThermoCloud<CloudType>::Tmin() const
418 {
419  scalar T = GREAT;
420  scalar n = 0;
421  forAllConstIter(typename ThermoCloud<CloudType>, *this, iter)
422  {
423  const parcelType& p = iter();
424  T = min(T, p.T());
425  n++;
426  }
427 
428  reduce(T, minOp<scalar>());
429  reduce(n, sumOp<label>());
430 
431  if (n > 0)
432  {
433  return T;
434  }
435  else
436  {
437  return 0.0;
438  }
439 }
440 
441 
442 // ************************************************************************* //
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
Top level model for Integration schemes.
dimensionedScalar pow3(const dimensionedScalar &ds)
const volScalarField & p() const
Return const access to the carrier prressure field.
Definition: ThermoCloudI.H:71
tmp< volScalarField > ap() const
Return tmp equivalent particulate absorption.
Definition: ThermoCloudI.H:317
scalar Tmin() const
Minimum temperature.
Definition: ThermoCloudI.H:417
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
Collection of constants.
labelList f(nPoints)
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
const parcelType::constantProperties & constProps() const
Return the constant properties.
Definition: ThermoCloudI.H:42
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
scalar epsilon
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
scalar Tmax() const
Maximum temperature.
Definition: ThermoCloudI.H:391
tmp< volScalarField > Ep() const
Return tmp equivalent particulate emission.
Definition: ThermoCloudI.H:281
const dimensionSet dimEnergy
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:68
messageStream Info
dynamicFvMesh & mesh
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:48
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:217
DimensionedField< scalar, volMesh > & radAreaPT4()
Radiation sum of parcel projected area*temperature^4 [m2K4].
Definition: ThermoCloudI.H:174
Thermo package for (S)olids (L)iquids and (G)ases Takes reference to thermo package, and provides:
Definition: SLGThermo.H:62
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:39
label n
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
tmp< fvMatrix< Type > > SuSp(const DimensionedField< scalar, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
volScalarField & p
Definition: createFields.H:51
DimensionedField< scalar, volMesh > & radT4()
Radiation sum of parcel temperature^4 [K4].
Definition: ThermoCloudI.H:138
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
DimensionedField< scalar, volMesh > & radAreaP()
Radiation sum of parcel projected areas [m2].
Definition: ThermoCloudI.H:102
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
const volScalarField & T() const
Return const access to the carrier temperature field.
Definition: ThermoCloudI.H:64
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const ThermoCloud & cloudCopy() const
Return a reference to the cloud copy.
Definition: ThermoCloudI.H:34
const HeatTransferModel< ThermoCloud< CloudType > > & heatTransfer() const
Return reference to heat transfer model.
Definition: ThermoCloudI.H:79
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
const scalarIntegrationScheme & TIntegrator() const
Return reference to velocity integration.
Definition: ThermoCloudI.H:87
const SLGThermo & thermo() const
Return const access to thermo package.
Definition: ThermoCloudI.H:57
error FatalError
tmp< volScalarField > sigmap() const
Return tmp equivalent particulate scattering factor.
Definition: ThermoCloudI.H:354
Templated base class for thermodynamic cloud.
Definition: ThermoCloud.H:60
DimensionedField< scalar, volMesh > & hsTrans()
Sensible enthalpy transfer [J/kg].
Definition: ThermoCloudI.H:210
tmp< fvScalarMatrix > Sh(volScalarField &hs) const
Return sensible enthalpy source term [J/kg/m3/s].
Definition: ThermoCloudI.H:242
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
Field< Type > & source()
Definition: fvMatrix.H:291
bool radiation() const
Radiation flag.
Definition: ThermoCloudI.H:94
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
Templated heat transfer model class.
Definition: ThermoCloud.H:53
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A class for managing temporary objects.
Definition: PtrList.H:118
const Type & value() const
Return const reference to value.
word timeName
Definition: getTimeIndex.H:3
DimensionedField< scalar, volMesh > & hsCoeff()
Return coefficient for carrier phase hs equation.
Definition: ThermoCloudI.H:226