ThermoCloudI.H
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-2022 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 "fvmSup.H"
28 
29 using namespace Foam::constant;
30 
31 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
32 
33 template<class CloudType>
34 inline const Foam::ThermoCloud<CloudType>&
36 {
37  return cloudCopyPtr_();
38 }
39 
40 
41 template<class CloudType>
42 inline const typename CloudType::particleType::constantProperties&
44 {
45  return constProps_;
46 }
47 
48 
49 template<class CloudType>
50 inline typename CloudType::particleType::constantProperties&
52 {
53  return constProps_;
54 }
55 
56 
57 template<class CloudType>
58 inline const Foam::fluidThermo&
60 {
61  return carrierThermo_;
62 }
63 
64 
65 template<class CloudType>
67 {
68  return thermo_;
69 }
70 
71 
72 template<class CloudType>
74 {
75  return T_;
76 }
77 
78 
79 template<class CloudType>
81 {
82  return p_;
83 }
84 
85 
86 template<class CloudType>
89 {
90  return heatTransferModel_;
91 }
92 
93 
94 template<class CloudType>
97 {
98  return compositionModel_;
99 }
100 
101 
102 template<class CloudType>
103 inline const Foam::integrationScheme&
105 {
106  return TIntegrator_;
107 }
108 
109 
110 template<class CloudType>
112 {
113  return radiation_;
114 }
115 
116 
117 template<class CloudType>
120 {
121  if (!radiation_)
122  {
124  << "Radiation field requested, but radiation model not active"
125  << abort(FatalError);
126  }
127 
128  return radAreaP_();
129 }
130 
131 
132 template<class CloudType>
135 {
136  if (!radiation_)
137  {
139  << "Radiation field requested, but radiation model not active"
140  << abort(FatalError);
141  }
142 
143  return radAreaP_();
144 }
145 
146 
147 template<class CloudType>
150 {
151  if (!radiation_)
152  {
154  << "Radiation field requested, but radiation model not active"
155  << abort(FatalError);
156  }
157 
158  return radT4_();
159 }
160 
161 
162 template<class CloudType>
165 {
166  if (!radiation_)
167  {
169  << "Radiation field requested, but radiation model not active"
170  << abort(FatalError);
171  }
172 
173  return radT4_();
174 }
175 
176 
177 template<class CloudType>
180 {
181  if (!radiation_)
182  {
184  << "Radiation field requested, but radiation model not active"
185  << abort(FatalError);
186  }
187 
188  return radAreaPT4_();
189 }
190 
191 
192 template<class CloudType>
195 {
196  if (!radiation_)
197  {
199  << "Radiation field requested, but radiation model not active"
200  << abort(FatalError);
201  }
202 
203  return radAreaPT4_();
204 }
205 
206 
207 template<class CloudType>
210 {
211  return hsTrans_();
212 }
213 
214 
215 template<class CloudType>
218 {
219  return hsTrans_();
220 }
221 
222 
223 template<class CloudType>
226 {
227  return hsCoeff_();
228 }
229 
230 
231 template<class CloudType>
234 {
235  return hsCoeff_();
236 }
237 
238 
239 template<class CloudType>
242 {
243  if (debug)
244  {
245  Info<< "hsTrans min/max = " << min(hsTrans()).value() << ", "
246  << max(hsTrans()).value() << nl
247  << "hsCoeff min/max = " << min(hsCoeff()).value() << ", "
248  << max(hsCoeff()).value() << endl;
249  }
250 
251  if (this->solution().coupled())
252  {
253  if (this->solution().semiImplicit("h"))
254  {
255  const volScalarField Cp(carrierThermo_.Cp());
257  Vdt(this->mesh().V()*this->db().time().deltaT());
258 
259  if (hs.dimensions() == dimTemperature)
260  {
261  return
262  hsTrans()/Vdt
263  - fvm::SuSp(hsCoeff()/Vdt, hs)
264  + (hsCoeff()/Vdt)*hs;
265  }
266  else
267  {
268  return
269  hsTrans()/Vdt
270  - fvm::SuSp(hsCoeff()/(Cp*Vdt), hs)
271  + hsCoeff()/(Cp*Vdt)*hs;
272  }
273  }
274  else
275  {
277  fvScalarMatrix& fvm = tfvm.ref();
278 
279  fvm.source() = -hsTrans()/(this->db().time().deltaT());
280 
281  return tfvm;
282  }
283  }
284 
286 }
287 
288 
289 template<class CloudType>
291 {
293  (
295  (
296  this->name() + ":radiation:Ep",
297  this->mesh(),
299  )
300  );
301 
302  if (radiation_)
303  {
304  scalarField& Ep = tEp.ref().primitiveFieldRef();
305  const scalar dt = this->db().time().deltaTValue();
306  const scalarField& V = this->mesh().V();
307  const scalar epsilon = constProps_.epsilon0();
308  const scalarField& sumAreaPT4 = radAreaPT4_->field();
309 
310  Ep = sumAreaPT4*epsilon*physicoChemical::sigma.value()/V/dt;
311  }
312 
313  return tEp;
314 }
315 
316 
317 template<class CloudType>
319 {
321  (
323  (
324  this->name() + ":radiation:ap",
325  this->mesh(),
327  )
328  );
329 
330  if (radiation_)
331  {
332  scalarField& ap = tap.ref().primitiveFieldRef();
333  const scalar dt = this->db().time().deltaTValue();
334  const scalarField& V = this->mesh().V();
335  const scalar epsilon = constProps_.epsilon0();
336  const scalarField& sumAreaP = radAreaP_->field();
337 
338  ap = sumAreaP*epsilon/V/dt;
339  }
340 
341  return tap;
342 }
343 
344 
345 template<class CloudType>
348 {
349  tmp<volScalarField> tsigmap
350  (
352  (
353  this->name() + ":radiation:sigmap",
354  this->mesh(),
356  )
357  );
358 
359  if (radiation_)
360  {
361  scalarField& sigmap = tsigmap.ref().primitiveFieldRef();
362  const scalar dt = this->db().time().deltaTValue();
363  const scalarField& V = this->mesh().V();
364  const scalar epsilon = constProps_.epsilon0();
365  const scalar f = constProps_.f0();
366  const scalarField& sumAreaP = radAreaP_->field();
367 
368  sigmap = sumAreaP*(1.0 - f)*(1.0 - epsilon)/V/dt;
369  }
370 
371  return tsigmap;
372 }
373 
374 
375 template<class CloudType>
376 inline Foam::scalar Foam::ThermoCloud<CloudType>::Tmax() const
377 {
378  scalar T = -great;
379  label n = 0;
380  forAllConstIter(typename ThermoCloud<CloudType>, *this, iter)
381  {
382  const parcelType& p = iter();
383  T = max(T, p.T());
384  n++;
385  }
386 
387  reduce(T, maxOp<scalar>());
388  reduce(n, sumOp<label>());
389 
390  if (n > 0)
391  {
392  return T;
393  }
394  else
395  {
396  return 0.0;
397  }
398 }
399 
400 
401 template<class CloudType>
402 inline Foam::scalar Foam::ThermoCloud<CloudType>::Tmin() const
403 {
404  scalar T = great;
405  label n = 0;
406  forAllConstIter(typename ThermoCloud<CloudType>, *this, iter)
407  {
408  const parcelType& p = iter();
409  T = min(T, p.T());
410  n++;
411  }
412 
413  reduce(T, minOp<scalar>());
414  reduce(n, sumOp<label>());
415 
416  if (n > 0)
417  {
418  return T;
419  }
420  else
421  {
422  return 0.0;
423  }
424 }
425 
426 
427 // ************************************************************************* //
Collection of constants.
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
Templated heat transfer model class.
Definition: ThermoCloud.H:58
const ThermoCloud & cloudCopy() const
Return a reference to the cloud copy.
Definition: ThermoCloudI.H:35
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
tmp< volScalarField > sigmap() const
Return tmp equivalent particulate scattering factor.
Definition: ThermoCloudI.H:347
Thermo package for (S)olids (L)iquids and (G)ases Takes reference to thermo package, and provides:
Definition: parcelThermo.H:58
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, const Internal &, const PtrList< fvPatchField< scalar >> &)
Return a temporary field constructed from name,.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
const volScalarField & p() const
Return const access to the carrier pressure field.
Definition: ThermoCloudI.H:80
volScalarField::Internal & radAreaPT4()
Radiation sum of parcel projected area*temperature^4 [m2K4].
Definition: ThermoCloudI.H:179
tmp< fvScalarMatrix > Sh(const volScalarField &hs) const
Return sensible enthalpy source term [J/kg/m^3/s].
Definition: ThermoCloudI.H:241
tmp< volScalarField::Internal > hsTrans() const
Return sensible enthalpy transfer [J/kg].
Definition: ThermoCloudI.H:209
const dimensionSet dimless
fvMesh & mesh
const dimensionSet dimLength
Templated base class for thermodynamic cloud.
Definition: ThermoCloud.H:76
const HeatTransferModel< ThermoCloud< CloudType > > & heatTransfer() const
Return reference to heat transfer model.
Definition: ThermoCloudI.H:88
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m^2/K^4].
volScalarField::Internal & hsCoeffRef()
Access coefficient for carrier phase hs equation.
Definition: ThermoCloudI.H:233
volScalarField::Internal & radAreaP()
Radiation sum of parcel projected areas [m^2].
Definition: ThermoCloudI.H:119
const dimensionSet dimTime
volScalarField::Internal & radT4()
Radiation sum of parcel temperature^4 [K4].
Definition: ThermoCloudI.H:149
const dimensionSet & dimensions() const
Return dimensions.
tmp< volScalarField > ap() const
Return tmp equivalent particulate absorption.
Definition: ThermoCloudI.H:318
const parcelType::constantProperties & constProps() const
Return the constant properties.
Definition: ThermoCloudI.H:43
Base-class for fluid thermodynamic properties.
Definition: fluidThermo.H:53
const Type & value() const
Return const reference to value.
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
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
bool radiation() const
Radiation flag.
Definition: ThermoCloudI.H:111
errorManip< error > abort(error &err)
Definition: errorManip.H:131
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:221
static const char nl
Definition: Ostream.H:260
const integrationScheme & TIntegrator() const
Return reference to velocity integration.
Definition: ThermoCloudI.H:104
Field< Type > & source()
Definition: fvMatrix.H:300
labelList f(nPoints)
const dimensionSet dimEnergy
const volScalarField & T() const
Return const access to the carrier temperature field.
Definition: ThermoCloudI.H:73
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
const dimensionSet dimMass
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
Base for a set of schemes which integrate simple ODEs which arise from semi-implcit rate expressions...
dimensionedScalar pow3(const dimensionedScalar &ds)
scalar Tmax() const
Maximum temperature.
Definition: ThermoCloudI.H:376
tmp< volScalarField > Ep() const
Return tmp equivalent particulate emission.
Definition: ThermoCloudI.H:290
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
volScalarField::Internal & hsTransRef()
Access sensible enthalpy transfer [J/kg].
Definition: ThermoCloudI.H:217
scalar epsilon
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
scalar Cp(const scalar p, const scalar T) const
Definition: EtoHthermo.H:2
const fluidThermo & carrierThermo() const
Return const access to carrier thermo package.
Definition: ThermoCloudI.H:59
messageStream Info
label n
tmp< volScalarField::Internal > hsCoeff() const
Return coefficient for carrier phase hs equation.
Definition: ThermoCloudI.H:225
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:48
const parcelThermo & thermo() const
Return const access to thermo package.
Definition: ThermoCloudI.H:66
const CompositionModel< ThermoCloud< CloudType > > & composition() const
Return const access to reacting composition model.
Definition: ThermoCloudI.H:96
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:53
Templated reacting parcel composition model class Consists of carrier species (via thermo package)...
Definition: ThermoCloud.H:61
scalar Tmin() const
Minimum temperature.
Definition: ThermoCloudI.H:402
Calculate the matrix for implicit and explicit sources.
const dimensionSet dimTemperature