ThermoCloud.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) 2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
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
27Class
28 Foam::ThermoCloud
29
30Group
31 grpLagrangianIntermediateClouds
32
33Description
34 Templated base class for thermodynamic cloud
35
36 - Adds to kinematic cloud
37 - Heat transfer
38
39SourceFiles
40 ThermoCloudI.H
41 ThermoCloud.C
42
43\*---------------------------------------------------------------------------*/
44
45#ifndef ThermoCloud_H
46#define ThermoCloud_H
47
48#include "KinematicCloud.H"
49#include "thermoCloud.H"
50#include "SLGThermo.H"
51
52// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
53
54namespace Foam
55{
56
57// Forward declaration of classes
58
59template<class CloudType>
60class HeatTransferModel;
61
62/*---------------------------------------------------------------------------*\
63 Class ThermoCloud Declaration
64\*---------------------------------------------------------------------------*/
65
66template<class CloudType>
67class ThermoCloud
68:
69 public CloudType,
70 public thermoCloud
71{
72public:
73
74 // Public typedefs
75
76 //- Type of cloud this cloud was instantiated for
77 typedef CloudType cloudType;
78
79 //- Type of parcel the cloud was instantiated for
80 typedef typename CloudType::particleType parcelType;
81
82 //- Convenience typedef for this cloud type
84
85
86private:
87
88 // Private data
89
90 //- Cloud copy pointer
91 autoPtr<ThermoCloud<CloudType>> cloudCopyPtr_;
92
93
94 // Private member functions
95
96 //- No copy construct
97 ThermoCloud(const ThermoCloud&) = delete;
98
99 //- No copy assignment
100 void operator=(const ThermoCloud&) = delete;
101
102
103protected:
104
105 // Protected data
106
107 //- Thermo parcel constant properties
109
110
111 // References to the carrier gas fields
112
113 //- SLG thermodynamics package
114 const SLGThermo& thermo_;
115
116 //- Temperature [K]
117 const volScalarField& T_;
118
119 //- Pressure [Pa]
120 const volScalarField& p_;
121
122
123 // References to the cloud sub-models
124
125 //- Heat transfer model
128
129
130 // Reference to the particle integration schemes
131
132 //- Temperature integration
134
135
136 // Modelling options
137
138 //- Include radiation
140
141 //- Radiation sum of parcel projected areas
143
144 //- Radiation sum of parcel temperature^4
146
147 //- Radiation sum of parcel projected areas * temperature^4
149
150
151 // Sources
152
153 //- Sensible enthalpy transfer [J/kg]
155
156 //- Coefficient for carrier phase hs equation [W/K]
158
159
160 // Protected Member Functions
161
162 // Initialisation
163
164 //- Set cloud sub-models
165 void setModels();
166
167
168 // Cloud evolution functions
169
170 //- Reset state of cloud
172
173
174public:
175
176 // Constructors
177
178 //- Construct given carrier gas fields
180 (
181 const word& cloudName,
182 const volScalarField& rho,
183 const volVectorField& U,
184 const dimensionedVector& g,
185 const SLGThermo& thermo,
186 bool readFields = true
187 );
188
189 //- Copy constructor with new name
191
192 //- Copy constructor with new name - creates bare cloud
194 (
195 const fvMesh& mesh,
196 const word& name,
198 );
199
200 //- Construct and return clone based on (this) with new name
202 {
204 (
205 new ThermoCloud(*this, name)
206 );
207 }
208
209 //- Construct and return bare clone based on (this) with new name
210 virtual autoPtr<Cloud<parcelType>> cloneBare(const word& name) const
211 {
213 (
214 new ThermoCloud(this->mesh(), name, *this)
215 );
216 }
217
218
219 //- Destructor
220 virtual ~ThermoCloud() = default;
221
222
223 // Member Functions
224
225 // Access
226
227 //- Return a reference to the cloud copy
228 inline const ThermoCloud& cloudCopy() const;
229
230 //- Return the constant properties
231 inline const typename parcelType::constantProperties&
232 constProps() const;
233
234 //- Return access to the constant properties
236
237 //- Return const access to thermo package
238 inline const SLGThermo& thermo() const;
239
240 //- Return const access to the carrier temperature field
241 inline const volScalarField& T() const;
242
243 //- Return const access to the carrier pressure field
244 inline const volScalarField& p() const;
245
246
247 // Sub-models
248
249 //- Return reference to heat transfer model
251 heatTransfer() const;
252
253
254 // Integration schemes
255
256 //-Return reference to velocity integration
257 inline const integrationScheme& TIntegrator() const;
258
259
260 // Modelling options
261
262 //- Radiation flag
263 inline bool radiation() const;
264
265 //- Radiation sum of parcel projected areas [m2]
267
268 //- Radiation sum of parcel projected areas [m2]
269 inline const volScalarField::Internal&
270 radAreaP() const;
271
272 //- Radiation sum of parcel temperature^4 [K4]
274
275 //- Radiation sum of parcel temperature^4 [K4]
276 inline const volScalarField::Internal& radT4() const;
277
278 //- Radiation sum of parcel projected area*temperature^4 [m2K4]
280
281 //- Radiation sum of parcel temperature^4 [m2K4]
282 inline const volScalarField::Internal&
283 radAreaPT4() const;
284
285
286 // Sources
287
288 // Enthalpy
289
290 //- Sensible enthalpy transfer [J/kg]
292
293 //- Sensible enthalpy transfer [J/kg]
294 inline const volScalarField::Internal&
295 hsTrans() const;
296
297 //- Return coefficient for carrier phase hs equation
299
300 //- Return const coefficient for carrier phase hs equation
301 inline const volScalarField::Internal&
302 hsCoeff() const;
303
304 //- Return sensible enthalpy source term [J/kg/m3/s]
305 inline tmp<fvScalarMatrix> Sh(volScalarField& hs) const;
306
307
308 // Radiation - overrides thermoCloud virtual abstract members
309
310 //- Return tmp equivalent particulate emission
311 inline tmp<volScalarField> Ep() const;
312
313 //- Return tmp equivalent particulate absorption
314 inline tmp<volScalarField> ap() const;
315
316 //- Return tmp equivalent particulate scattering factor
317 inline tmp<volScalarField> sigmap() const;
318
319
320 // Check
321
322 //- Maximum temperature
323 inline scalar Tmax() const;
324
325 //- Minimum temperature
326 inline scalar Tmin() const;
327
328
329 // Cloud evolution functions
330
331 //- Set parcel thermo properties
333 (
334 parcelType& parcel,
335 const scalar lagrangianDt
336 );
337
338 //- Check parcel properties
340 (
341 parcelType& parcel,
342 const scalar lagrangianDt,
343 const bool fullyDescribed
344 );
345
346 //- Store the current cloud state
347 void storeState();
348
349 //- Reset the current cloud to the previously stored state
350 void restoreState();
351
352 //- Reset the cloud source terms
353 void resetSourceTerms();
354
355 //- Apply relaxation to (steady state) cloud sources
356 void relaxSources(const ThermoCloud<CloudType>& cloudOldTime);
357
358 //- Apply scaling to (transient) cloud sources
359 void scaleSources();
360
361 //- Pre-evolve
362 void preEvolve(const typename parcelType::trackingData& td);
363
364 //- Evolve the cloud
365 void evolve();
366
367
368 // Mapping
369
370 //- Remap the cells of particles corresponding to the
371 // mesh topology change with a default tracking data object
372 virtual void autoMap(const mapPolyMesh&);
373
374
375 // I-O
376
377 //- Print cloud information
378 void info();
379};
380
381
382// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
383
384} // End namespace Foam
385
386// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
387
388#include "ThermoCloudI.H"
389
390// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
391
392#ifdef NoRepository
393 #include "ThermoCloud.C"
394#endif
395
396// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
397
398#endif
399
400// ************************************************************************* //
const uniformDimensionedVectorField & g
ParticleType particleType
Definition: Cloud.H:114
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:75
const word & cloudName() const
Return the cloud type.
Definition: DSMCCloudI.H:37
const fvMesh & mesh() const
Return reference to the mesh.
Definition: DSMCCloudI.H:44
Class to hold DSMC particle constant properties.
Definition: DSMCParcel.H:82
Templated class to calculate the fluid-particle heat transfer coefficients based on a specified Nusse...
const word & name() const noexcept
Return the object name.
Definition: IOobjectI.H:65
autoPtr< IOobject > clone() const
Clone.
Definition: IOobject.H:459
Thermo package for (S)olids (L)iquids and (G)ases Takes reference to thermo package,...
Definition: SLGThermo.H:67
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:78
Templated base class for thermodynamic cloud.
Definition: ThermoCloud.H:70
const parcelType::constantProperties & constProps() const
Return the constant properties.
Definition: ThermoCloudI.H:45
volScalarField::Internal & radAreaP()
Radiation sum of parcel projected areas [m2].
Definition: ThermoCloudI.H:105
void setModels()
Set cloud sub-models.
Definition: ThermoCloud.C:37
volScalarField::Internal & hsCoeff()
Return coefficient for carrier phase hs equation.
Definition: ThermoCloudI.H:211
volScalarField::Internal & radAreaPT4()
Radiation sum of parcel projected area*temperature^4 [m2K4].
Definition: ThermoCloudI.H:165
virtual ~ThermoCloud()=default
Destructor.
void storeState()
Store the current cloud state.
Definition: ThermoCloud.C:376
const volScalarField & T() const
Return const access to the carrier temperature field.
Definition: ThermoCloudI.H:67
void setParcelThermoProperties(parcelType &parcel, const scalar lagrangianDt)
Set parcel thermo properties.
Definition: ThermoCloud.C:351
bool radiation() const
Radiation flag.
Definition: ThermoCloudI.H:97
void checkParcelProperties(parcelType &parcel, const scalar lagrangianDt, const bool fullyDescribed)
Check parcel properties.
Definition: ThermoCloud.C:365
scalar Tmin() const
Minimum temperature.
Definition: ThermoCloudI.H:397
const volScalarField & p_
Pressure [Pa].
Definition: ThermoCloud.H:119
volScalarField::Internal & radT4()
Radiation sum of parcel temperature^4 [K4].
Definition: ThermoCloudI.H:135
void relaxSources(const ThermoCloud< CloudType > &cloudOldTime)
Apply relaxation to (steady state) cloud sources.
Definition: ThermoCloud.C:414
autoPtr< volScalarField::Internal > hsCoeff_
Coefficient for carrier phase hs equation [W/K].
Definition: ThermoCloud.H:156
autoPtr< volScalarField::Internal > radT4_
Radiation sum of parcel temperature^4.
Definition: ThermoCloud.H:144
void scaleSources()
Apply scaling to (transient) cloud sources.
Definition: ThermoCloud.C:433
autoPtr< integrationScheme > TIntegrator_
Temperature integration.
Definition: ThermoCloud.H:132
CloudType::particleType parcelType
Type of parcel the cloud was instantiated for.
Definition: ThermoCloud.H:79
Switch radiation_
Include radiation.
Definition: ThermoCloud.H:138
const SLGThermo & thermo() const
Return const access to thermo package.
Definition: ThermoCloudI.H:60
const integrationScheme & TIntegrator() const
Return reference to velocity integration.
Definition: ThermoCloudI.H:90
CloudType cloudType
Type of cloud this cloud was instantiated for.
Definition: ThermoCloud.H:76
const ThermoCloud & cloudCopy() const
Return a reference to the cloud copy.
Definition: ThermoCloudI.H:37
virtual void autoMap(const mapPolyMesh &)
Remap the cells of particles corresponding to the.
Definition: ThermoCloud.C:474
scalar Tmax() const
Maximum temperature.
Definition: ThermoCloudI.H:376
tmp< volScalarField > sigmap() const
Return tmp equivalent particulate scattering factor.
Definition: ThermoCloudI.H:339
virtual autoPtr< Cloud< parcelType > > clone(const word &name)
Construct and return clone based on (this) with new name.
Definition: ThermoCloud.H:200
virtual autoPtr< Cloud< parcelType > > cloneBare(const word &name) const
Construct and return bare clone based on (this) with new name.
Definition: ThermoCloud.H:209
void cloudReset(ThermoCloud< CloudType > &c)
Reset state of cloud.
Definition: ThermoCloud.C:116
const HeatTransferModel< ThermoCloud< CloudType > > & heatTransfer() const
Return reference to heat transfer model.
Definition: ThermoCloudI.H:82
tmp< fvScalarMatrix > Sh(volScalarField &hs) const
Return sensible enthalpy source term [J/kg/m3/s].
Definition: ThermoCloudI.H:227
void evolve()
Evolve the cloud.
Definition: ThermoCloud.C:462
autoPtr< HeatTransferModel< ThermoCloud< CloudType > > > heatTransferModel_
Heat transfer model.
Definition: ThermoCloud.H:126
parcelType::constantProperties constProps_
Thermo parcel constant properties.
Definition: ThermoCloud.H:107
autoPtr< volScalarField::Internal > radAreaPT4_
Radiation sum of parcel projected areas * temperature^4.
Definition: ThermoCloud.H:147
ThermoCloud< CloudType > thermoCloudType
Convenience typedef for this cloud type.
Definition: ThermoCloud.H:82
volScalarField::Internal & hsTrans()
Sensible enthalpy transfer [J/kg].
Definition: ThermoCloudI.H:195
void preEvolve(const typename parcelType::trackingData &td)
Pre-evolve.
Definition: ThermoCloud.C:451
void info()
Print cloud information.
Definition: ThermoCloud.C:483
void restoreState()
Reset the current cloud to the previously stored state.
Definition: ThermoCloud.C:389
void resetSourceTerms()
Reset the cloud source terms.
Definition: ThermoCloud.C:397
const SLGThermo & thermo_
SLG thermodynamics package.
Definition: ThermoCloud.H:113
const volScalarField & T_
Temperature [K].
Definition: ThermoCloud.H:116
const volScalarField & p() const
Return const access to the carrier pressure field.
Definition: ThermoCloudI.H:74
tmp< volScalarField > Ep() const
Return tmp equivalent particulate emission.
Definition: ThermoCloudI.H:266
autoPtr< volScalarField::Internal > hsTrans_
Sensible enthalpy transfer [J/kg].
Definition: ThermoCloud.H:153
autoPtr< volScalarField::Internal > radAreaP_
Radiation sum of parcel projected areas.
Definition: ThermoCloud.H:141
tmp< volScalarField > ap() const
Return tmp equivalent particulate absorption.
Definition: ThermoCloudI.H:302
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: autoPtr.H:66
Class used to pass data into container.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
Base for a set of schemes which integrate simple ODEs which arise from semi-implcit rate expressions.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:162
Virtual abstract base class for templated ThermoCloud.
Definition: thermoCloud.H:51
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
A class for managing temporary objects.
Definition: tmp.H:65
A class for handling words, derived from Foam::string.
Definition: word.H:68
U
Definition: pEqn.H:72
Namespace for OpenFOAM.
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const wordHashSet &selectedFields, LIFOStack< regIOobject * > &storedObjects)
Read the selected GeometricFields of the templated type.