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-2021 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().timeName(),
116  this->db(),
117  IOobject::READ_IF_PRESENT,
118  IOobject::AUTO_WRITE
119  ),
120  this->mesh(),
122  )
123  );
124  }
125 
126  if (this->solution().resetSourcesOnStartup())
127  {
128  resetSourceTerms();
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().timeName(),
158  this->db(),
159  IOobject::NO_READ,
160  IOobject::NO_WRITE,
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,
175  const ReactingCloud<CloudType>& c
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  const scalar lagrangianDt
200 )
201 {
202  CloudType::setParcelThermoProperties(parcel, lagrangianDt);
203 
204  parcel.Y() = this->composition().YMixture0();
205 }
206 
207 
208 template<class CloudType>
210 (
211  parcelType& parcel,
212  const scalar lagrangianDt,
213  const bool fullyDescribed
214 )
215 {
216  CloudType::checkParcelProperties(parcel, lagrangianDt, fullyDescribed);
217 
218  if (fullyDescribed)
219  {
220  checkSuppliedComposition
221  (
222  parcel.Y(),
223  this->composition().YMixture0(),
224  "YMixture"
225  );
226  }
227 
228  // derived information - store initial mass
229  parcel.mass0() = parcel.mass();
230 }
231 
232 
233 template<class CloudType>
235 {
236  cloudCopyPtr_.reset
237  (
238  static_cast<ReactingCloud<CloudType>*>
239  (
240  clone(this->name() + "Copy").ptr()
241  )
242  );
243 }
244 
245 
246 template<class CloudType>
248 {
249  cloudReset(cloudCopyPtr_());
250  cloudCopyPtr_.clear();
251 }
252 
253 
254 template<class CloudType>
256 {
257  CloudType::resetSourceTerms();
258  forAll(rhoTrans_, i)
259  {
260  rhoTrans_[i].field() = 0.0;
261  }
262 }
263 
264 
265 template<class CloudType>
267 (
268  const ReactingCloud<CloudType>& cloudOldTime
269 )
270 {
271  CloudType::relaxSources(cloudOldTime);
272 
273  typedef volScalarField::Internal dsfType;
274 
275  forAll(rhoTrans_, fieldi)
276  {
277  dsfType& rhoT = rhoTrans_[fieldi];
278  const dsfType& rhoT0 = cloudOldTime.rhoTrans()[fieldi];
279  this->relax(rhoT, rhoT0, "rho");
280  }
281 }
282 
283 
284 template<class CloudType>
286 {
287  CloudType::scaleSources();
288 
289  typedef volScalarField::Internal dsfType;
290 
291  forAll(rhoTrans_, fieldi)
292  {
293  dsfType& rhoT = rhoTrans_[fieldi];
294  this->scale(rhoT, "rho");
295  }
296 }
297 
298 
299 template<class CloudType>
301 {
302  if (this->solution().canEvolve())
303  {
304  typename parcelType::trackingData td(*this);
305 
306  this->solve(*this, td);
307  }
308 }
309 
310 
311 template<class CloudType>
313 {
315 
316  this->updateMesh();
317 }
318 
319 
320 template<class CloudType>
322 {
323  CloudType::info();
324 
325  this->phaseChange().info(Info);
326 }
327 
328 
329 // ************************************************************************* //
virtual ~ReactingCloud()
Destructor.
void resetSourceTerms()
Reset the cloud source terms.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
DSMCCloud< dsmcParcel > CloudType
Templated phase change model class.
Definition: ReactingCloud.H:57
void storeState()
Store the current cloud state.
parcelType::constantProperties constProps_
Parcel constant properties.
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
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
basicSpecieMixture & composition
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
void setParcelThermoProperties(parcelType &parcel, const scalar lagrangianDt)
Set parcel thermo properties.
void restoreState()
Reset the current cloud to the previously stored state.
twoPhaseChangeModel & phaseChange
void setModels()
Set cloud sub-models.
Definition: ReactingCloud.C:33
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
void evolve()
Evolve the cloud.
T clone(const T &t)
Definition: List.H:54
dynamicFvMesh & mesh
autoPtr< BasicCompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleMomentumTransportModel::transportModel &transport)
A class for handling words, derived from string.
Definition: word.H:59
Base-class for fluid thermodynamic properties.
Definition: fluidThermo.H:53
autoPtr< PhaseChangeModel< ReactingCloud< CloudType > > > phaseChangeModel_
Reacting phase change model.
word timeName
Definition: getTimeIndex.H:3
errorManip< error > abort(error &err)
Definition: errorManip.H:131
EaEqn relax()
void checkParcelProperties(parcelType &parcel, const scalar lagrangianDt, const bool fullyDescribed)
Check parcel properties.
Base cloud calls templated on particle type.
Definition: Cloud.H:52
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:221
volScalarField::Internal & rhoTrans(const label i)
Mass.
static const char nl
Definition: Ostream.H:260
const dimensionSet dimMass
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
virtual void autoMap(const mapPolyMesh &)
Remap the cells of particles corresponding to the.
void scaleSources()
Apply scaling to (transient) cloud sources.
Templated base class for reacting cloud.
Definition: ReactingCloud.H:72
rhoEqn solve()
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
void info()
Print cloud information.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
void cloudReset(ReactingCloud< CloudType > &c)
Reset state of cloud.
Definition: ReactingCloud.C:67
PtrList< volScalarField::Internal > rhoTrans_
Mass transfer fields - one per carrier phase specie.
messageStream Info
void checkSuppliedComposition(const scalarField &YSupplied, const scalarField &Y, const word &YName)
Check that size of a composition field is valid.
Definition: ReactingCloud.C:48
Selector class for relaxation factors, solver type and solution.
Definition: solution.H: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
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
void relaxSources(const ReactingCloud< CloudType > &cloudOldTime)
Apply relaxation to (steady state) cloud sources.