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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2017 OpenFOAM Foundation
9  Copyright (C) 2019-2022 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
30 
31 using namespace Foam::constant;
32 
33 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
34 
35 template<class CloudType>
36 inline const Foam::ThermoCloud<CloudType>&
38 {
39  return *cloudCopyPtr_;
40 }
41 
42 
43 template<class CloudType>
44 inline const typename CloudType::particleType::constantProperties&
46 {
47  return constProps_;
48 }
49 
50 
51 template<class CloudType>
52 inline typename CloudType::particleType::constantProperties&
54 {
55  return constProps_;
56 }
57 
58 
59 template<class CloudType>
61 {
62  return thermo_;
63 }
64 
65 
66 template<class CloudType>
68 {
69  return T_;
70 }
71 
72 
73 template<class CloudType>
75 {
76  return p_;
77 }
78 
79 
80 template<class CloudType>
83 {
84  return *heatTransferModel_;
85 }
86 
87 
88 template<class CloudType>
89 inline const Foam::integrationScheme&
91 {
92  return *TIntegrator_;
93 }
94 
95 
96 template<class CloudType>
98 {
99  return radiation_;
100 }
101 
102 
103 template<class CloudType>
106 {
107  if (!radiation_)
108  {
110  << "Radiation field requested, but radiation model not active"
111  << abort(FatalError);
112  }
114  return *radAreaP_;
115 }
116 
117 
118 template<class CloudType>
121 {
122  if (!radiation_)
123  {
125  << "Radiation field requested, but radiation model not active"
126  << abort(FatalError);
127  }
129  return *radAreaP_;
130 }
131 
132 
133 template<class CloudType>
136 {
137  if (!radiation_)
138  {
140  << "Radiation field requested, but radiation model not active"
141  << abort(FatalError);
142  }
144  return *radT4_;
145 }
146 
147 
148 template<class CloudType>
151 {
152  if (!radiation_)
153  {
155  << "Radiation field requested, but radiation model not active"
156  << abort(FatalError);
157  }
159  return *radT4_;
160 }
161 
162 
163 template<class CloudType>
166 {
167  if (!radiation_)
168  {
170  << "Radiation field requested, but radiation model not active"
171  << abort(FatalError);
172  }
174  return *radAreaPT4_;
175 }
176 
177 
178 template<class CloudType>
181 {
182  if (!radiation_)
183  {
185  << "Radiation field requested, but radiation model not active"
186  << abort(FatalError);
187  }
189  return *radAreaPT4_;
190 }
191 
192 
193 template<class CloudType>
196 {
197  return *hsTrans_;
198 }
199 
200 
201 template<class CloudType>
204 {
205  return *hsTrans_;
206 }
207 
208 
209 template<class CloudType>
212 {
213  return *hsCoeff_;
214 }
215 
216 
217 template<class CloudType>
220 {
221  return *hsCoeff_;
222 }
223 
224 
225 template<class CloudType>
228 {
229  DebugInfo
230  << "hsTrans min/max = " << min(hsTrans()).value() << ", "
231  << max(hsTrans()).value() << nl
232  << "hsCoeff min/max = " << min(hsCoeff()).value() << ", "
233  << max(hsCoeff()).value() << endl;
234 
235  if (this->solution().coupled())
236  {
237  if (this->solution().semiImplicit("h"))
238  {
239  const volScalarField Cp(thermo_.thermo().Cp());
241  Vdt(this->mesh().V()*this->db().time().deltaT());
242 
243  return
244  hsTrans()/Vdt
245  - fvm::SuSp(hsCoeff()/(Cp*Vdt), hs)
246  + hsCoeff()/(Cp*Vdt)*hs;
247  }
248  else
249  {
250  tmp<fvScalarMatrix> tfvm(new fvScalarMatrix(hs, dimEnergy/dimTime));
251  fvScalarMatrix& fvm = tfvm.ref();
252 
253  fvm.source() = -hsTrans()/(this->db().time().deltaT());
254 
255  return tfvm;
256  }
257  }
258 
260 }
261 
262 
263 template<class CloudType>
265 {
266  tmp<volScalarField> tEp
267  (
268  new volScalarField
269  (
270  IOobject
271  (
272  this->name() + ":radiation:Ep",
273  this->db().time().timeName(),
274  this->db(),
277  false
278  ),
279  this->mesh(),
281  )
282  );
283 
284  if (radiation_)
285  {
286  scalarField& Ep = tEp.ref().primitiveFieldRef();
287  const scalar dt = this->db().time().deltaTValue();
288  const scalarField& V = this->mesh().V();
289  const scalar epsilon = constProps_.epsilon0();
290  const scalarField& sumAreaPT4 = radAreaPT4_->field();
291 
292  Ep = sumAreaPT4*epsilon*physicoChemical::sigma.value()/V/dt;
293  }
294 
295  return tEp;
296 }
297 
298 
299 template<class CloudType>
301 {
302  tmp<volScalarField> tap
303  (
304  new volScalarField
305  (
306  IOobject
307  (
308  this->name() + ":radiation:ap",
309  this->db().time().timeName(),
310  this->db(),
313  false
314  ),
315  this->mesh(),
317  )
318  );
319 
320  if (radiation_)
321  {
322  scalarField& ap = tap.ref().primitiveFieldRef();
323  const scalar dt = this->db().time().deltaTValue();
324  const scalarField& V = this->mesh().V();
325  const scalar epsilon = constProps_.epsilon0();
326  const scalarField& sumAreaP = radAreaP_->field();
327 
328  ap = sumAreaP*epsilon/V/dt;
329  }
331  return tap;
332 }
333 
334 
335 template<class CloudType>
338 {
339  tmp<volScalarField> tsigmap
340  (
341  new volScalarField
342  (
343  IOobject
344  (
345  this->name() + ":radiation:sigmap",
346  this->db().time().timeName(),
347  this->db(),
350  false
351  ),
352  this->mesh(),
354  )
355  );
356 
357  if (radiation_)
358  {
359  scalarField& sigmap = tsigmap.ref().primitiveFieldRef();
360  const scalar dt = this->db().time().deltaTValue();
361  const scalarField& V = this->mesh().V();
362  const scalar epsilon = constProps_.epsilon0();
363  const scalar f = constProps_.f0();
364  const scalarField& sumAreaP = radAreaP_->field();
365 
366  sigmap = sumAreaP*(1.0 - f)*(1.0 - epsilon)/V/dt;
367  }
368 
369  return tsigmap;
370 }
371 
372 
373 template<class CloudType>
374 inline Foam::scalar Foam::ThermoCloud<CloudType>::Tmax() const
375 {
376  scalar val = -GREAT;
377  bool nonEmpty = false;
378 
379  for (const parcelType& p : *this)
380  {
381  val = max(val, p.T());
382  nonEmpty = true;
383  }
384 
385  if (returnReduceOr(nonEmpty))
386  {
387  return returnReduce(val, maxOp<scalar>());
388  }
389 
390  return 0;
391 }
392 
393 
394 template<class CloudType>
395 inline Foam::scalar Foam::ThermoCloud<CloudType>::Tmin() const
396 {
397  scalar val = GREAT;
398  bool nonEmpty = false;
399 
400  for (const parcelType& p : *this)
401  {
402  val = min(val, p.T());
403  nonEmpty = true;
404  }
405 
406  if (returnReduceOr(nonEmpty))
407  {
408  return returnReduce(val, minOp<scalar>());
409  }
410 
411  return 0;
412 }
413 
414 
415 // ************************************************************************* //
Different types of constants.
const Type & value() const noexcept
Return const reference to value.
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:37
Templated class to calculate the fluid-particle heat transfer coefficients based on a specified Nusse...
Definition: ThermoCloud.H:55
volScalarField::Internal & hsTrans()
Sensible enthalpy transfer [J/kg].
Definition: ThermoCloudI.H:188
const ThermoCloud & cloudCopy() const
Return a reference to the cloud copy.
Definition: ThermoCloudI.H:30
error FatalError
Error stream (stdout output on all processes), with additional &#39;FOAM FATAL ERROR&#39; header text and sta...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:578
tmp< volScalarField > sigmap() const
Return tmp equivalent particulate scattering factor.
Definition: ThermoCloudI.H:330
T & ref() const
Return non-const reference to the contents of a non-null managed pointer.
Definition: tmpI.H:210
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:40
constexpr char nl
The newline &#39;\n&#39; character (0x0a)
Definition: Ostream.H:49
bool coupled(solutionDict.getOrDefault("coupledEnergyField", false))
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:487
const volScalarField & p() const
Return const access to the carrier pressure field.
Definition: ThermoCloudI.H:67
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
volScalarField::Internal & radAreaPT4()
Radiation sum of parcel projected area*temperature^4 [m2K4].
Definition: ThermoCloudI.H:158
Ignore writing from objectRegistry::writeObject()
const dimensionSet dimless
Dimensionless.
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Perform reduction on a copy, using specified binary operation.
Templated base class for thermodynamic cloud.
Definition: ThermoCloud.H:62
const HeatTransferModel< ThermoCloud< CloudType > > & heatTransfer() const
Return reference to heat transfer model.
Definition: ThermoCloudI.H:75
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:84
volScalarField::Internal & radAreaP()
Radiation sum of parcel projected areas [m2].
Definition: ThermoCloudI.H:98
word timeName
Definition: getTimeIndex.H:3
volScalarField::Internal & radT4()
Radiation sum of parcel temperature^4 [K4].
Definition: ThermoCloudI.H:128
tmp< volScalarField > ap() const
Return tmp equivalent particulate absorption.
Definition: ThermoCloudI.H:293
dynamicFvMesh & mesh
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:52
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
const parcelType::constantProperties & constProps() const
Return the constant properties.
Definition: ThermoCloudI.H:38
static tmp< T > New(Args &&... args)
Construct tmp with forwarding arguments.
Definition: tmp.H:203
zeroField SuSp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
Thermo package for (S)olids (L)iquids and (G)ases Takes reference to thermo package, and provides:
Definition: SLGThermo.H:60
bool radiation() const
Radiation flag.
Definition: ThermoCloudI.H:90
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:26
errorManip< error > abort(error &err)
Definition: errorManip.H:139
#define DebugInfo
Report an information message using Foam::Info.
const integrationScheme & TIntegrator() const
Return reference to velocity integration.
Definition: ThermoCloudI.H:83
const volScalarField & Cp
Definition: EEqn.H:7
labelList f(nPoints)
const SLGThermo & thermo() const
Return const access to thermo package.
Definition: ThermoCloudI.H:53
const volScalarField & T() const
Return const access to the carrier temperature field.
Definition: ThermoCloudI.H:60
const dimensionSet dimEnergy
Base for a set of schemes which integrate simple ODEs which arise from semi-implcit rate expressions...
tmp< fvScalarMatrix > Sh(volScalarField &hs) const
Return sensible enthalpy source term [J/kg/m3/s].
Definition: ThermoCloudI.H:220
dimensionedScalar pow3(const dimensionedScalar &ds)
scalar Tmax() const
Maximum temperature.
Definition: ThermoCloudI.H:367
tmp< volScalarField > Ep() const
Return tmp equivalent particulate emission.
Definition: ThermoCloudI.H:257
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
scalar epsilon
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Definition: areaFieldsFwd.H:42
Nothing to be read.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:59
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
volScalarField & p
A class for managing temporary objects.
Definition: HashPtrTable.H:50
bool returnReduceOr(const bool value, const label comm=UPstream::worldComm)
Perform logical (or) MPI Allreduce on a copy. Uses UPstream::reduceOr.
volScalarField::Internal & hsCoeff()
Return coefficient for carrier phase hs equation.
Definition: ThermoCloudI.H:204
Defines the attributes of an object for which implicit objectRegistry management is supported...
Definition: IOobject.H:166
scalar Tmin() const
Minimum temperature.
Definition: ThermoCloudI.H:388
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:157