ReactingCloud.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 "ReactingCloud.H"
27 #include "CompositionModel.H"
28 #include "PhaseChangeModel.H"
29 
30 // * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * * //
31 
32 template<class CloudType>
34 {
35  phaseChangeModel_.reset
36  (
38  (
39  this->subModelProperties(),
40  *this
41  ).ptr()
42  );
43 }
44 
45 
46 template<class CloudType>
48 (
49  const scalarField& YSupplied,
50  const scalarField& Y,
51  const word& YName
52 )
53 {
54  if (YSupplied.size() != Y.size())
55  {
57  << YName << " supplied, but size is not compatible with "
58  << "parcel composition: " << nl << " "
59  << YName << "(" << YSupplied.size() << ") vs required composition "
60  << YName << "(" << Y.size() << ")" << nl
61  << abort(FatalError);
62  }
63 }
64 
65 
66 template<class CloudType>
68 {
69  CloudType::cloudReset(c);
70 
71  phaseChangeModel_.reset(c.phaseChangeModel_.ptr());
72 }
73 
74 
75 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
76 
77 template<class CloudType>
79 (
80  const word& cloudName,
81  const volScalarField& rho,
82  const volVectorField& U,
83  const dimensionedVector& g,
84  const fluidThermo& carrierThermo,
85  const bool readFields
86 )
87 :
88  CloudType(cloudName, rho, U, g, carrierThermo, false),
89  cloudCopyPtr_(nullptr),
90  constProps_(this->particleProperties()),
91  phaseChangeModel_(nullptr)
92 {
93  setModels();
94 
95  rhoTrans_.setSize(this->composition().carrier().species().size());
96 
97  if (readFields)
98  {
99  parcelType::readFields(*this, this->composition());
100  this->deleteLostParticles();
101  }
102 
103  // Set storage for mass source fields and initialise to zero
104  forAll(rhoTrans_, i)
105  {
106  const word& specieName = this->composition().carrier().species()[i];
107  rhoTrans_.set
108  (
109  i,
111  (
112  IOobject
113  (
114  this->name() + ":rhoTrans_" + specieName,
115  this->db().time().name(),
116  this->db(),
119  ),
120  this->mesh(),
122  )
123  );
124  }
125 
126  if (this->solution().resetSourcesOnStartup())
127  {
129  }
130 }
131 
132 
133 template<class CloudType>
135 (
137  const word& name
138 )
139 :
140  CloudType(c, name),
141  cloudCopyPtr_(nullptr),
142  constProps_(c.constProps_),
143  phaseChangeModel_(c.phaseChangeModel_->clone()),
144  rhoTrans_(c.rhoTrans_.size())
145 {
146  forAll(c.rhoTrans_, i)
147  {
148  const word& specieName = this->composition().carrier().species()[i];
149  rhoTrans_.set
150  (
151  i,
153  (
154  IOobject
155  (
156  this->name() + ":rhoTrans_" + specieName,
157  this->db().time().name(),
158  this->db(),
161  false
162  ),
163  c.rhoTrans_[i]
164  )
165  );
166  }
167 }
168 
169 
170 template<class CloudType>
172 (
173  const fvMesh& mesh,
174  const word& name,
176 )
177 :
178  CloudType(mesh, name, c),
179  cloudCopyPtr_(nullptr),
180  constProps_(),
181  phaseChangeModel_(nullptr),
182  rhoTrans_(0)
183 {}
184 
185 
186 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
187 
188 template<class CloudType>
190 {}
191 
192 
193 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
194 
195 template<class CloudType>
197 (
198  parcelType& parcel
199 )
200 {
201  CloudType::setParcelThermoProperties(parcel);
202 
203  parcel.Y() = this->composition().YMixture0();
204 }
205 
206 
207 template<class CloudType>
209 (
210  parcelType& parcel,
211  const label injectori
212 )
213 {
214  CloudType::checkParcelProperties(parcel, injectori);
215 
216  if (injectori != -1 && this->injectors()[injectori].fullyDescribed())
217  {
218  checkSuppliedComposition
219  (
220  parcel.Y(),
221  this->composition().YMixture0(),
222  "YMixture"
223  );
224  }
225 }
226 
227 
228 template<class CloudType>
230 {
231  cloudCopyPtr_.reset
232  (
233  static_cast<ReactingCloud<CloudType>*>
234  (
235  clone(this->name() + "Copy").ptr()
236  )
237  );
238 }
239 
240 
241 template<class CloudType>
243 {
244  cloudReset(cloudCopyPtr_());
245  cloudCopyPtr_.clear();
246 }
247 
248 
249 template<class CloudType>
251 {
252  CloudType::resetSourceTerms();
253  forAll(rhoTrans_, i)
254  {
255  rhoTrans_[i].primitiveFieldRef() = 0.0;
256  }
257 }
258 
259 
260 template<class CloudType>
262 (
263  const ReactingCloud<CloudType>& cloudOldTime
264 )
265 {
266  CloudType::relaxSources(cloudOldTime);
267 
268  typedef volScalarField::Internal dsfType;
269 
270  forAll(rhoTrans_, fieldi)
271  {
272  dsfType& rhoT = rhoTrans_[fieldi];
273  const dsfType& rhoT0 = cloudOldTime.rhoTrans()[fieldi];
274  this->relax(rhoT, rhoT0, "rho");
275  }
276 }
277 
278 
279 template<class CloudType>
281 {
282  CloudType::scaleSources();
283 
284  typedef volScalarField::Internal dsfType;
285 
286  forAll(rhoTrans_, fieldi)
287  {
288  dsfType& rhoT = rhoTrans_[fieldi];
289  this->scale(rhoT, "rho");
290  }
291 }
292 
293 
294 template<class CloudType>
296 {
297  if (this->solution().canEvolve())
298  {
299  typename parcelType::trackingData td(*this);
300 
301  this->solve(*this, td);
302  }
303 }
304 
305 
306 template<class CloudType>
308 {
309  CloudType::info();
310 
311  this->phaseChange().info(Info);
312 }
313 
314 
315 // ************************************************************************* //
EaEqn relax()
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
void deleteLostParticles()
Remove lost particles from cloud and delete.
Definition: Cloud.C:233
label size() const
Return the number of particles in the cloud.
Definition: Cloud.H:191
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:80
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:225
const fvMesh & mesh() const
Return references to the mesh.
Definition: DSMCCloudI.H:41
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.
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
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
Templated phase change model class.
Templated base class for reacting cloud.
Definition: ReactingCloud.H:76
PtrList< volScalarField::Internal > rhoTrans_
Mass transfer fields - one per carrier phase specie.
void setModels()
Set cloud sub-models.
Definition: ReactingCloud.C:33
void storeState()
Store the current cloud state.
virtual ~ReactingCloud()
Destructor.
void cloudReset(ReactingCloud< CloudType > &c)
Reset state of cloud.
Definition: ReactingCloud.C:67
void scaleSources()
Apply scaling to (transient) cloud sources.
void checkSuppliedComposition(const scalarField &YSupplied, const scalarField &Y, const word &YName)
Check that size of a composition field is valid.
Definition: ReactingCloud.C:48
ReactingCloud(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: ReactingCloud.C:79
volScalarField::Internal & rhoTrans(const label i)
Mass.
void evolve()
Evolve the cloud.
void relaxSources(const ReactingCloud< CloudType > &cloudOldTime)
Apply relaxation to (steady state) cloud sources.
void info()
Print cloud information.
void restoreState()
Reset the current cloud to the previously stored state.
void checkParcelProperties(parcelType &parcel, const label injectori)
Check parcel properties.
void resetSourceTerms()
Reset the cloud source terms.
void setParcelThermoProperties(parcelType &parcel)
Set parcel thermo properties.
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
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
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
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
errorManip< error > abort(error &err)
Definition: errorManip.H:131
messageStream Info
T clone(const T &t)
Definition: List.H:55
const dimensionSet dimMass
error FatalError
static const char nl
Definition: Ostream.H:266
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.
PtrList< volScalarField > & Y
const word cloudName(propsDict.lookup("cloudName"))