ThermoCloud.C
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-2018 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 "ThermoCloud.H"
27 #include "ThermoParcel.H"
28 
29 #include "HeatTransferModel.H"
30 
31 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
32 
33 template<class CloudType>
35 {
36  heatTransferModel_.reset
37  (
39  (
40  this->subModelProperties(),
41  *this
42  ).ptr()
43  );
44 
45  TIntegrator_.reset
46  (
48  (
49  "T",
50  this->solution().integrationSchemes()
51  ).ptr()
52  );
53 
54  if (this->solution().coupled())
55  {
56  this->subModelProperties().lookup("radiation") >> radiation_;
57  }
58 
59  if (radiation_)
60  {
61  radAreaP_.reset
62  (
64  (
65  IOobject
66  (
67  this->name() + ":radAreaP",
68  this->db().time().timeName(),
69  this->db(),
70  IOobject::READ_IF_PRESENT,
71  IOobject::AUTO_WRITE
72  ),
73  this->mesh(),
75  )
76  );
77 
78  radT4_.reset
79  (
81  (
82  IOobject
83  (
84  this->name() + ":radT4",
85  this->db().time().timeName(),
86  this->db(),
87  IOobject::READ_IF_PRESENT,
88  IOobject::AUTO_WRITE
89  ),
90  this->mesh(),
92  )
93  );
94 
95  radAreaPT4_.reset
96  (
98  (
99  IOobject
100  (
101  this->name() + ":radAreaPT4",
102  this->db().time().timeName(),
103  this->db(),
104  IOobject::READ_IF_PRESENT,
105  IOobject::AUTO_WRITE
106  ),
107  this->mesh(),
109  )
110  );
111  }
112 }
113 
114 
115 template<class CloudType>
117 {
118  CloudType::cloudReset(c);
119 
120  heatTransferModel_.reset(c.heatTransferModel_.ptr());
121  TIntegrator_.reset(c.TIntegrator_.ptr());
122 
123  radiation_ = c.radiation_;
124 }
125 
126 
127 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
128 
129 template<class CloudType>
131 (
132  const word& cloudName,
133  const volScalarField& rho,
134  const volVectorField& U,
135  const dimensionedVector& g,
136  const SLGThermo& thermo,
137  bool readFields
138 )
139 :
140  CloudType
141  (
142  cloudName,
143  rho,
144  U,
145  thermo.thermo().mu(),
146  g,
147  false
148  ),
149  thermoCloud(),
150  cloudCopyPtr_(nullptr),
151  constProps_(this->particleProperties()),
152  thermo_(thermo),
153  T_(thermo.thermo().T()),
154  p_(thermo.thermo().p()),
155  heatTransferModel_(nullptr),
156  TIntegrator_(nullptr),
157  radiation_(false),
158  radAreaP_(nullptr),
159  radT4_(nullptr),
160  radAreaPT4_(nullptr),
161  hsTrans_
162  (
164  (
165  IOobject
166  (
167  this->name() + ":hsTrans",
168  this->db().time().timeName(),
169  this->db(),
170  IOobject::READ_IF_PRESENT,
171  IOobject::AUTO_WRITE
172  ),
173  this->mesh(),
175  )
176  ),
177  hsCoeff_
178  (
180  (
181  IOobject
182  (
183  this->name() + ":hsCoeff",
184  this->db().time().timeName(),
185  this->db(),
186  IOobject::READ_IF_PRESENT,
187  IOobject::AUTO_WRITE
188  ),
189  this->mesh(),
191  )
192  )
193 {
194  if (this->solution().active())
195  {
196  setModels();
197 
198  if (readFields)
199  {
200  parcelType::readFields(*this);
201  this->deleteLostParticles();
202  }
203  }
204 
205  if (this->solution().resetSourcesOnStartup())
206  {
207  resetSourceTerms();
208  }
209 }
210 
211 
212 template<class CloudType>
214 (
216  const word& name
217 )
218 :
219  CloudType(c, name),
220  thermoCloud(),
221  cloudCopyPtr_(nullptr),
222  constProps_(c.constProps_),
223  thermo_(c.thermo_),
224  T_(c.T()),
225  p_(c.p()),
226  heatTransferModel_(c.heatTransferModel_->clone()),
227  TIntegrator_(c.TIntegrator_->clone()),
228  radiation_(c.radiation_),
229  radAreaP_(nullptr),
230  radT4_(nullptr),
231  radAreaPT4_(nullptr),
232  hsTrans_
233  (
235  (
236  IOobject
237  (
238  this->name() + ":hsTrans",
239  this->db().time().timeName(),
240  this->db(),
241  IOobject::NO_READ,
242  IOobject::NO_WRITE,
243  false
244  ),
245  c.hsTrans()
246  )
247  ),
248  hsCoeff_
249  (
251  (
252  IOobject
253  (
254  this->name() + ":hsCoeff",
255  this->db().time().timeName(),
256  this->db(),
257  IOobject::NO_READ,
258  IOobject::NO_WRITE,
259  false
260  ),
261  c.hsCoeff()
262  )
263  )
264 {
265  if (radiation_)
266  {
267  radAreaP_.reset
268  (
270  (
271  IOobject
272  (
273  this->name() + ":radAreaP",
274  this->db().time().timeName(),
275  this->db(),
276  IOobject::NO_READ,
277  IOobject::NO_WRITE,
278  false
279  ),
280  c.radAreaP()
281  )
282  );
283 
284  radT4_.reset
285  (
287  (
288  IOobject
289  (
290  this->name() + ":radT4",
291  this->db().time().timeName(),
292  this->db(),
293  IOobject::NO_READ,
294  IOobject::NO_WRITE,
295  false
296  ),
297  c.radT4()
298  )
299  );
300 
301  radAreaPT4_.reset
302  (
304  (
305  IOobject
306  (
307  this->name() + ":radAreaPT4",
308  this->db().time().timeName(),
309  this->db(),
310  IOobject::NO_READ,
311  IOobject::NO_WRITE,
312  false
313  ),
314  c.radAreaPT4()
315  )
316  );
317  }
318 }
319 
320 
321 template<class CloudType>
323 (
324  const fvMesh& mesh,
325  const word& name,
326  const ThermoCloud<CloudType>& c
327 )
328 :
329  CloudType(mesh, name, c),
330  thermoCloud(),
331  cloudCopyPtr_(nullptr),
332  constProps_(),
333  thermo_(c.thermo()),
334  T_(c.T()),
335  p_(c.p()),
336  heatTransferModel_(nullptr),
337  TIntegrator_(nullptr),
338  radiation_(false),
339  radAreaP_(nullptr),
340  radT4_(nullptr),
341  radAreaPT4_(nullptr),
342  hsTrans_(nullptr),
343  hsCoeff_(nullptr)
344 {}
345 
346 
347 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
348 
349 template<class CloudType>
351 {}
352 
353 
354 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
355 
356 template<class CloudType>
358 (
359  parcelType& parcel,
360  const scalar lagrangianDt
361 )
362 {
363  CloudType::setParcelThermoProperties(parcel, lagrangianDt);
364 
365  parcel.T() = constProps_.T0();
366  parcel.Cp() = constProps_.Cp0();
367 }
368 
369 
370 template<class CloudType>
372 (
373  parcelType& parcel,
374  const scalar lagrangianDt,
375  const bool fullyDescribed
376 )
377 {
378  CloudType::checkParcelProperties(parcel, lagrangianDt, fullyDescribed);
379 }
380 
381 
382 template<class CloudType>
384 {
385  cloudCopyPtr_.reset
386  (
387  static_cast<ThermoCloud<CloudType>*>
388  (
389  clone(this->name() + "Copy").ptr()
390  )
391  );
392 }
393 
394 
395 template<class CloudType>
397 {
398  cloudReset(cloudCopyPtr_());
399  cloudCopyPtr_.clear();
400 }
401 
402 
403 template<class CloudType>
405 {
406  CloudType::resetSourceTerms();
407  hsTrans_->field() = 0.0;
408  hsCoeff_->field() = 0.0;
409 
410  if (radiation_)
411  {
412  radAreaP_->field() = 0.0;
413  radT4_->field() = 0.0;
414  radAreaPT4_->field() = 0.0;
415  }
416 }
417 
418 
419 template<class CloudType>
421 (
422  const ThermoCloud<CloudType>& cloudOldTime
423 )
424 {
425  CloudType::relaxSources(cloudOldTime);
426 
427  this->relax(hsTrans_(), cloudOldTime.hsTrans(), "h");
428  this->relax(hsCoeff_(), cloudOldTime.hsCoeff(), "h");
429 
430  if (radiation_)
431  {
432  this->relax(radAreaP_(), cloudOldTime.radAreaP(), "radiation");
433  this->relax(radT4_(), cloudOldTime.radT4(), "radiation");
434  this->relax(radAreaPT4_(), cloudOldTime.radAreaPT4(), "radiation");
435  }
436 }
437 
438 
439 template<class CloudType>
441 {
442  CloudType::scaleSources();
443 
444  this->scale(hsTrans_(), "h");
445  this->scale(hsCoeff_(), "h");
446 
447  if (radiation_)
448  {
449  this->scale(radAreaP_(), "radiation");
450  this->scale(radT4_(), "radiation");
451  this->scale(radAreaPT4_(), "radiation");
452  }
453 }
454 
455 
456 template<class CloudType>
458 {
459  CloudType::preEvolve();
460 
461  this->pAmbient() = thermo_.thermo().p().average().value();
462 }
463 
464 
465 template<class CloudType>
467 {
468  if (this->solution().canEvolve())
469  {
470  typename parcelType::trackingData td(*this);
471 
472  this->solve(*this, td);
473  }
474 }
475 
476 
477 template<class CloudType>
479 {
481 
482  this->updateMesh();
483 }
484 
485 
486 template<class CloudType>
488 {
489  CloudType::info();
490 
491  Info<< " Temperature min/max = " << Tmin() << ", " << Tmax()
492  << endl;
493 }
494 
495 
496 // ************************************************************************* //
autoPtr< integrationScheme > TIntegrator_
Temperature integration.
Definition: ThermoCloud.H:117
Switch radiation_
Include radiation.
Definition: ThermoCloud.H:123
virtual tmp< volScalarField > mu() const =0
Dynamic viscosity of mixture [kg/m/s].
DSMCCloud< dsmcParcel > CloudType
const SLGThermo & thermo_
SLG thermodynamics package.
Definition: ThermoCloud.H:98
void setParcelThermoProperties(parcelType &parcel, const scalar lagrangianDt)
Set parcel thermo properties.
Definition: ThermoCloud.C:358
Templated heat transfer model class.
Definition: ThermoCloud.H:53
volScalarField::Internal & hsTrans()
Sensible enthalpy transfer [J/kg].
Definition: ThermoCloudI.H:192
UEqn relax()
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const HashSet< word > &selectedFields, LIFOStack< regIOobject *> &storedObjects)
Read the selected GeometricFields of the specified type.
Definition: ReadFields.C:244
dimensionedSymmTensor sqr(const dimensionedVector &dv)
void restoreState()
Reset the current cloud to the previously stored state.
Definition: ThermoCloud.C:396
ThermoCloud(const word &cloudName, const volScalarField &rho, const volVectorField &U, const dimensionedVector &g, const SLGThermo &thermo, bool readFields=true)
Construct given carrier gas fields.
Definition: ThermoCloud.C:131
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
const volScalarField & p() const
Return const access to the carrier prressure field.
Definition: ThermoCloudI.H:71
virtual const volScalarField & T() const
Temperature [K].
Definition: basicThermo.C:489
volScalarField::Internal & radAreaPT4()
Radiation sum of parcel projected area*temperature^4 [m2K4].
Definition: ThermoCloudI.H:162
void cloudReset(ThermoCloud< CloudType > &c)
Reset state of cloud.
Definition: ThermoCloud.C:116
void preEvolve()
Pre-evolve.
Definition: ThermoCloud.C:457
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
Templated base class for thermodynamic cloud.
Definition: ThermoCloud.H:60
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
parcelType::constantProperties constProps_
Thermo parcel constant properties.
Definition: ThermoCloud.H:92
virtual volScalarField & p()
Pressure [Pa].
Definition: basicThermo.C:477
volScalarField::Internal & radAreaP()
Radiation sum of parcel projected areas [m^2].
Definition: ThermoCloudI.H:102
volScalarField::Internal & radT4()
Radiation sum of parcel temperature^4 [K4].
Definition: ThermoCloudI.H:132
T clone(const T &t)
Definition: List.H:54
dynamicFvMesh & mesh
virtual void autoMap(const mapPolyMesh &)
Remap the cells of particles corresponding to the.
Definition: ThermoCloud.C:478
A class for handling words, derived from string.
Definition: word.H:59
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
Definition: dimensionSets.H:52
rhoEqn solve()
word timeName
Definition: getTimeIndex.H:3
Thermo package for (S)olids (L)iquids and (G)ases Takes reference to thermo package, and provides:
Definition: SLGThermo.H:62
Base cloud calls templated on particle type.
Definition: Cloud.H:52
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:215
Virtual abstract base class for templated ThermoCloud.
Definition: thermoCloud.H:48
const SLGThermo & thermo() const
Return const access to thermo package.
Definition: ThermoCloudI.H:57
const volScalarField & T() const
Return const access to the carrier temperature field.
Definition: ThermoCloudI.H:64
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const dimensionSet dimEnergy
void checkParcelProperties(parcelType &parcel, const scalar lagrangianDt, const bool fullyDescribed)
Check parcel properties.
Definition: ThermoCloud.C:372
void resetSourceTerms()
Reset the cloud source terms.
Definition: ThermoCloud.C:404
void evolve()
Evolve the cloud.
Definition: ThermoCloud.C:466
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
void info()
Print cloud information.
Definition: ThermoCloud.C:487
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
dimensionedScalar pow4(const dimensionedScalar &ds)
void scaleSources()
Apply scaling to (transient) cloud sources.
Definition: ThermoCloud.C:440
messageStream Info
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:48
autoPtr< HeatTransferModel< ThermoCloud< CloudType > > > heatTransferModel_
Heat transfer model.
Definition: ThermoCloud.H:111
const dimensionedVector & g
volScalarField::Internal & hsCoeff()
Return coefficient for carrier phase hs equation.
Definition: ThermoCloudI.H:208
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
void storeState()
Store the current cloud state.
Definition: ThermoCloud.C:383
void setModels()
Set cloud sub-models.
Definition: ThermoCloud.C:34
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:57
const fluidThermo & thermo() const
Return reference to the thermo database.
Definition: SLGThermo.C:102
void relaxSources(const ThermoCloud< CloudType > &cloudOldTime)
Apply relaxation to (steady state) cloud sources.
Definition: ThermoCloud.C:421
virtual ~ThermoCloud()
Destructor.
Definition: ThermoCloud.C:350