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-2024 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 "integrationScheme.H"
28 #include "HeatTransferModel.H"
29 #include "CompositionModel.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  compositionModel_.reset
46  (
48  (
49  this->subModelProperties(),
50  *this
51  ).ptr()
52  );
53 
54  TIntegrator_.reset
55  (
57  (
58  "T",
59  this->solution().integrationSchemes()
60  ).ptr()
61  );
62 
63  if (this->solution().coupled())
64  {
65  this->subModelProperties().lookup("radiation") >> radiation_;
66  }
67 
68  if (radiation_)
69  {
70  radAreaP_.reset
71  (
73  (
74  IOobject
75  (
76  this->name() + ":radAreaP",
77  this->db().time().name(),
78  this->db(),
79  IOobject::READ_IF_PRESENT,
80  IOobject::AUTO_WRITE
81  ),
82  this->mesh(),
84  )
85  );
86 
87  radT4_.reset
88  (
90  (
91  IOobject
92  (
93  this->name() + ":radT4",
94  this->db().time().name(),
95  this->db(),
96  IOobject::READ_IF_PRESENT,
97  IOobject::AUTO_WRITE
98  ),
99  this->mesh(),
101  )
102  );
103 
104  radAreaPT4_.reset
105  (
107  (
108  IOobject
109  (
110  this->name() + ":radAreaPT4",
111  this->db().time().name(),
112  this->db(),
113  IOobject::READ_IF_PRESENT,
114  IOobject::AUTO_WRITE
115  ),
116  this->mesh(),
118  )
119  );
120  }
121 }
122 
123 
124 template<class CloudType>
126 {
127  CloudType::cloudReset(c);
128 
129  heatTransferModel_.reset(c.heatTransferModel_.ptr());
130  compositionModel_.reset(c.compositionModel_.ptr());
131 
132  TIntegrator_.reset(c.TIntegrator_.ptr());
133 
134  radiation_ = c.radiation_;
135 }
136 
137 
138 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
139 
140 template<class CloudType>
142 (
143  const word& cloudName,
144  const volScalarField& rho,
145  const volVectorField& U,
146  const dimensionedVector& g,
147  const fluidThermo& carrierThermo,
148  const bool readFields
149 )
150 :
151  CloudType(cloudName, rho, U, g, carrierThermo, false),
152  cloudCopyPtr_(nullptr),
153  constProps_(this->particleProperties()),
154  carrierThermo_(carrierThermo),
155  thermo_(carrierThermo_),
156  T_(carrierThermo.T()),
157  p_(carrierThermo.p()),
158  heatTransferModel_(nullptr),
159  compositionModel_(nullptr),
160  TIntegrator_(nullptr),
161  radiation_(false),
162  radAreaP_(nullptr),
163  radT4_(nullptr),
164  radAreaPT4_(nullptr),
165  hsTrans_
166  (
167  new volScalarField::Internal
168  (
169  IOobject
170  (
171  this->name() + ":hsTrans",
172  this->db().time().name(),
173  this->db(),
174  IOobject::READ_IF_PRESENT,
175  IOobject::AUTO_WRITE
176  ),
177  this->mesh(),
179  )
180  ),
181  hsCoeff_
182  (
183  new volScalarField::Internal
184  (
185  IOobject
186  (
187  this->name() + ":hsCoeff",
188  this->db().time().name(),
189  this->db(),
190  IOobject::READ_IF_PRESENT,
191  IOobject::AUTO_WRITE
192  ),
193  this->mesh(),
195  )
196  )
197 {
198  setModels();
199 
200  if (readFields)
201  {
202  parcelType::readFields(*this);
203  this->deleteLostParticles();
204  }
205 
206  if (this->solution().resetSourcesOnStartup())
207  {
209  }
210 }
211 
212 
213 template<class CloudType>
215 (
217  const word& name
218 )
219 :
220  CloudType(c, name),
221  cloudCopyPtr_(nullptr),
222  constProps_(c.constProps_),
223  carrierThermo_(c.carrierThermo_),
224  thermo_(c.thermo_),
225  T_(c.T()),
226  p_(c.p()),
227  heatTransferModel_(c.heatTransferModel_->clone()),
228  compositionModel_(c.compositionModel_->clone()),
229  TIntegrator_(c.TIntegrator_->clone()),
230  radiation_(c.radiation_),
231  radAreaP_(nullptr),
232  radT4_(nullptr),
233  radAreaPT4_(nullptr),
234  hsTrans_
235  (
236  new volScalarField::Internal
237  (
238  IOobject
239  (
240  this->name() + ":hsTrans",
241  this->db().time().name(),
242  this->db(),
243  IOobject::NO_READ,
244  IOobject::NO_WRITE,
245  false
246  ),
247  c.hsTrans()
248  )
249  ),
250  hsCoeff_
251  (
252  new volScalarField::Internal
253  (
254  IOobject
255  (
256  this->name() + ":hsCoeff",
257  this->db().time().name(),
258  this->db(),
259  IOobject::NO_READ,
260  IOobject::NO_WRITE,
261  false
262  ),
263  c.hsCoeff()
264  )
265  )
266 {
267  if (radiation_)
268  {
269  radAreaP_.reset
270  (
272  (
273  IOobject
274  (
275  this->name() + ":radAreaP",
276  this->db().time().name(),
277  this->db(),
280  false
281  ),
282  c.radAreaP()
283  )
284  );
285 
286  radT4_.reset
287  (
289  (
290  IOobject
291  (
292  this->name() + ":radT4",
293  this->db().time().name(),
294  this->db(),
297  false
298  ),
299  c.radT4()
300  )
301  );
302 
303  radAreaPT4_.reset
304  (
306  (
307  IOobject
308  (
309  this->name() + ":radAreaPT4",
310  this->db().time().name(),
311  this->db(),
314  false
315  ),
316  c.radAreaPT4()
317  )
318  );
319  }
320 }
321 
322 
323 template<class CloudType>
325 (
326  const fvMesh& mesh,
327  const word& name,
329 )
330 :
331  CloudType(mesh, name, c),
332  cloudCopyPtr_(nullptr),
333  constProps_(),
334  carrierThermo_(c.carrierThermo_),
335  thermo_(c.thermo()),
336  T_(c.T()),
337  p_(c.p()),
338  heatTransferModel_(nullptr),
339  compositionModel_(c.compositionModel_->clone()),
340  TIntegrator_(nullptr),
341  radiation_(false),
342  radAreaP_(nullptr),
343  radT4_(nullptr),
344  radAreaPT4_(nullptr),
345  hsTrans_(nullptr),
346  hsCoeff_(nullptr)
347 {}
348 
349 
350 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
351 
352 template<class CloudType>
354 {}
355 
356 
357 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
358 
359 template<class CloudType>
361 (
362  parcelType& parcel
363 )
364 {
365  CloudType::setParcelThermoProperties(parcel);
366 
367  parcel.T() = constProps_.T0();
368  parcel.Cp() = constProps_.Cp0();
369 }
370 
371 
372 template<class CloudType>
374 (
375  parcelType& parcel,
376  const label injectori
377 )
378 {
379  CloudType::checkParcelProperties(parcel, injectori);
380 }
381 
382 
383 template<class CloudType>
385 {
386  cloudCopyPtr_.reset
387  (
388  static_cast<ThermoCloud<CloudType>*>
389  (
390  clone(this->name() + "Copy").ptr()
391  )
392  );
393 }
394 
395 
396 template<class CloudType>
398 {
399  cloudReset(cloudCopyPtr_());
400  cloudCopyPtr_.clear();
401 }
402 
403 
404 template<class CloudType>
406 {
407  CloudType::resetSourceTerms();
408  hsTrans_->primitiveFieldRef() = 0.0;
409  hsCoeff_->primitiveFieldRef() = 0.0;
410 
411  if (radiation_)
412  {
413  radAreaP_->primitiveFieldRef() = 0.0;
414  radT4_->primitiveFieldRef() = 0.0;
415  radAreaPT4_->primitiveFieldRef() = 0.0;
416  }
417 }
418 
419 
420 template<class CloudType>
422 (
423  const ThermoCloud<CloudType>& cloudOldTime
424 )
425 {
426  CloudType::relaxSources(cloudOldTime);
427 
428  this->relax(hsTrans_(), cloudOldTime.hsTrans_(), "h");
429  this->relax(hsCoeff_(), cloudOldTime.hsCoeff_(), "h");
430 
431  if (radiation_)
432  {
433  this->relax(radAreaP_(), cloudOldTime.radAreaP(), "radiation");
434  this->relax(radT4_(), cloudOldTime.radT4(), "radiation");
435  this->relax(radAreaPT4_(), cloudOldTime.radAreaPT4(), "radiation");
436  }
437 }
438 
439 
440 template<class CloudType>
442 {
443  CloudType::scaleSources();
444 
445  this->scale(hsTrans_(), "h");
446  this->scale(hsCoeff_(), "h");
447 
448  if (radiation_)
449  {
450  this->scale(radAreaP_(), "radiation");
451  this->scale(radT4_(), "radiation");
452  this->scale(radAreaPT4_(), "radiation");
453  }
454 }
455 
456 
457 template<class CloudType>
459 {
460  CloudType::preEvolve();
461 
462  this->pAmbient() = carrierThermo_.p().average().value();
463 }
464 
465 
466 template<class CloudType>
468 {
469  if (this->solution().canEvolve())
470  {
471  typename parcelType::trackingData td(*this);
472 
473  this->solve(*this, td);
474  }
475 }
476 
477 
478 template<class CloudType>
480 {
481  CloudType::info();
482 
483  Info<< " Temperature min/max = " << Tmin() << ", " << Tmax()
484  << endl;
485 }
486 
487 
488 template<class CloudType>
490 {
491  if (compositionModel_.valid())
492  {
493  CloudType::particleType::writeFields(*this, this->composition());
494  }
495  else
496  {
497  CloudType::particleType::writeFields(*this);
498  }
499 }
500 
501 
502 // ************************************************************************* //
EaEqn relax()
void deleteLostParticles()
Remove lost particles from cloud and delete.
Definition: Cloud.C:233
Templated reacting parcel composition model class Consists of carrier species (via thermo package),...
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:80
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:225
void info() const
Print cloud information.
Definition: DSMCCloud.C:996
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Generic GeometricField class.
Templated heat transfer model class.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:312
const word & name() const
Return name.
Definition: IOobject.H:310
Templated base class for thermodynamic cloud.
Definition: ThermoCloud.H:83
volScalarField::Internal & radAreaP()
Radiation sum of parcel projected areas [m^2].
Definition: ThermoCloudI.H:119
void setModels()
Set cloud sub-models.
Definition: ThermoCloud.C:34
volScalarField::Internal & radAreaPT4()
Radiation sum of parcel projected area*temperature^4 [m2K4].
Definition: ThermoCloudI.H:179
void storeState()
Store the current cloud state.
Definition: ThermoCloud.C:384
volScalarField::Internal & radT4()
Radiation sum of parcel temperature^4 [K4].
Definition: ThermoCloudI.H:149
void relaxSources(const ThermoCloud< CloudType > &cloudOldTime)
Apply relaxation to (steady state) cloud sources.
Definition: ThermoCloud.C:422
autoPtr< volScalarField::Internal > hsCoeff_
Coefficient for carrier phase hs equation [W/K].
Definition: ThermoCloud.H:167
autoPtr< volScalarField::Internal > radT4_
Radiation sum of parcel temperature^4.
Definition: ThermoCloud.H:155
void scaleSources()
Apply scaling to (transient) cloud sources.
Definition: ThermoCloud.C:441
ThermoCloud(const word &cloudName, const volScalarField &rho, const volVectorField &U, const dimensionedVector &g, const fluidThermo &carrierThermo, const bool readFields=true)
Construct given carrier fields and thermo.
Definition: ThermoCloud.C:142
Switch radiation_
Include radiation.
Definition: ThermoCloud.H:149
virtual void writeFields() const
Write the field data for the cloud.
Definition: ThermoCloud.C:489
void cloudReset(ThermoCloud< CloudType > &c)
Reset state of cloud.
Definition: ThermoCloud.C:125
void evolve()
Evolve the cloud.
Definition: ThermoCloud.C:467
autoPtr< volScalarField::Internal > radAreaPT4_
Radiation sum of parcel projected areas * temperature^4.
Definition: ThermoCloud.H:158
void info()
Print cloud information.
Definition: ThermoCloud.C:479
void restoreState()
Reset the current cloud to the previously stored state.
Definition: ThermoCloud.C:397
void checkParcelProperties(parcelType &parcel, const label injectori)
Check parcel properties.
Definition: ThermoCloud.C:374
void preEvolve()
Pre-evolve.
Definition: ThermoCloud.C:458
void resetSourceTerms()
Reset the cloud source terms.
Definition: ThermoCloud.C:405
virtual ~ThermoCloud()
Destructor.
Definition: ThermoCloud.C:353
autoPtr< volScalarField::Internal > hsTrans_
Sensible enthalpy transfer [J/kg].
Definition: ThermoCloud.H:164
autoPtr< volScalarField::Internal > radAreaP_
Radiation sum of parcel projected areas.
Definition: ThermoCloud.H:152
void setParcelThermoProperties(parcelType &parcel)
Set parcel thermo properties.
Definition: ThermoCloud.C:361
Base-class for fluid thermodynamic properties.
Definition: fluidThermo.H:57
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
const Time & time() const
Return time.
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:51
A class for handling words, derived from string.
Definition: word.H:62
U
Definition: pEqn.H:72
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
const dimensionedScalar c
Speed of light in a vacuum.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
const dimensionSet dimEnergy
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
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
messageStream Info
const dimensionSet dimLength
const dimensionSet dimTemperature
T clone(const T &t)
Definition: List.H:55
dimensionedScalar pow4(const dimensionedScalar &ds)
const dimensionSet dimArea
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.
volScalarField & p
fluidMulticomponentThermo & thermo
Definition: createFields.H:31
const word cloudName(propsDict.lookup("cloudName"))