ReactingCloud.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2013 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 
28 #include "CompositionModel.H"
29 #include "PhaseChangeModel.H"
30 
31 // * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * * //
32 
33 template<class CloudType>
35 {
36  compositionModel_.reset
37  (
39  (
40  this->subModelProperties(),
41  *this
42  ).ptr()
43  );
44 
45  phaseChangeModel_.reset
46  (
48  (
49  this->subModelProperties(),
50  *this
51  ).ptr()
52  );
53 }
54 
55 
56 template<class CloudType>
58 (
59  const scalarField& YSupplied,
60  const scalarField& Y,
61  const word& YName
62 )
63 {
64  if (YSupplied.size() != Y.size())
65  {
67  (
68  "ReactingCloud<CloudType>::checkSuppliedComposition"
69  "("
70  "const scalarField&, "
71  "const scalarField&, "
72  "const word&"
73  ")"
74  ) << YName << " supplied, but size is not compatible with "
75  << "parcel composition: " << nl << " "
76  << YName << "(" << YSupplied.size() << ") vs required composition "
77  << YName << "(" << Y.size() << ")" << nl
78  << abort(FatalError);
79  }
80 }
81 
82 
83 template<class CloudType>
85 {
86  CloudType::cloudReset(c);
87 
88  compositionModel_.reset(c.compositionModel_.ptr());
89  phaseChangeModel_.reset(c.phaseChangeModel_.ptr());
90 }
91 
92 
93 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
94 
95 template<class CloudType>
97 (
98  const word& cloudName,
99  const volScalarField& rho,
100  const volVectorField& U,
101  const dimensionedVector& g,
102  const SLGThermo& thermo,
103  bool readFields
104 )
105 :
106  CloudType(cloudName, rho, U, g, thermo, false),
107  reactingCloud(),
108  cloudCopyPtr_(NULL),
109  constProps_(this->particleProperties()),
110  compositionModel_(NULL),
111  phaseChangeModel_(NULL),
112  rhoTrans_(thermo.carrier().species().size())
113 {
114  if (this->solution().active())
115  {
116  setModels();
117 
118  if (readFields)
119  {
120  parcelType::readFields(*this, this->composition());
121  }
122  }
123 
124  // Set storage for mass source fields and initialise to zero
125  forAll(rhoTrans_, i)
126  {
127  const word& specieName = thermo.carrier().species()[i];
128  rhoTrans_.set
129  (
130  i,
132  (
133  IOobject
134  (
135  this->name() + ":rhoTrans_" + specieName,
136  this->db().time().timeName(),
137  this->db(),
138  IOobject::READ_IF_PRESENT,
139  IOobject::AUTO_WRITE
140  ),
141  this->mesh(),
142  dimensionedScalar("zero", dimMass, 0.0)
143  )
144  );
145  }
146 
147  if (this->solution().resetSourcesOnStartup())
148  {
149  resetSourceTerms();
150  }
151 }
152 
153 
154 template<class CloudType>
156 (
158  const word& name
159 )
160 :
161  CloudType(c, name),
162  reactingCloud(),
163  cloudCopyPtr_(NULL),
164  constProps_(c.constProps_),
165  compositionModel_(c.compositionModel_->clone()),
166  phaseChangeModel_(c.phaseChangeModel_->clone()),
167  rhoTrans_(c.rhoTrans_.size())
168 {
169  forAll(c.rhoTrans_, i)
170  {
171  const word& specieName = this->thermo().carrier().species()[i];
172  rhoTrans_.set
173  (
174  i,
176  (
177  IOobject
178  (
179  this->name() + ":rhoTrans_" + specieName,
180  this->db().time().timeName(),
181  this->db(),
182  IOobject::NO_READ,
183  IOobject::NO_WRITE,
184  false
185  ),
186  c.rhoTrans_[i]
187  )
188  );
189  }
190 }
191 
192 
193 template<class CloudType>
195 (
196  const fvMesh& mesh,
197  const word& name,
198  const ReactingCloud<CloudType>& c
199 )
200 :
201  CloudType(mesh, name, c),
202  reactingCloud(),
203  cloudCopyPtr_(NULL),
204  constProps_(),
205  compositionModel_(c.compositionModel_->clone()),
206 // compositionModel_(NULL),
207  phaseChangeModel_(NULL),
208  rhoTrans_(0)
209 {}
210 
211 
212 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
213 
214 template<class CloudType>
216 {}
217 
218 
219 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
220 
221 template<class CloudType>
223 (
224  parcelType& parcel,
225  const scalar lagrangianDt
226 )
227 {
228  CloudType::setParcelThermoProperties(parcel, lagrangianDt);
229 
230  parcel.pc() = this->thermo().thermo().p()[parcel.cell()];
231  parcel.Y() = composition().YMixture0();
232 }
233 
234 
235 template<class CloudType>
237 (
238  parcelType& parcel,
239  const scalar lagrangianDt,
240  const bool fullyDescribed
241 )
242 {
243  CloudType::checkParcelProperties(parcel, lagrangianDt, fullyDescribed);
244 
245  if (fullyDescribed)
246  {
247  checkSuppliedComposition
248  (
249  parcel.Y(),
250  composition().YMixture0(),
251  "YMixture"
252  );
253  }
254 
255  // derived information - store initial mass
256  parcel.mass0() = parcel.mass();
257 }
258 
259 
260 template<class CloudType>
262 {
263  cloudCopyPtr_.reset
264  (
265  static_cast<ReactingCloud<CloudType>*>
266  (
267  clone(this->name() + "Copy").ptr()
268  )
269  );
270 }
271 
272 
273 template<class CloudType>
275 {
276  cloudReset(cloudCopyPtr_());
277  cloudCopyPtr_.clear();
278 }
279 
280 
281 template<class CloudType>
283 {
284  CloudType::resetSourceTerms();
285  forAll(rhoTrans_, i)
286  {
287  rhoTrans_[i].field() = 0.0;
288  }
289 }
290 
291 
292 template<class CloudType>
294 (
295  const ReactingCloud<CloudType>& cloudOldTime
296 )
297 {
298  CloudType::relaxSources(cloudOldTime);
299 
300  typedef DimensionedField<scalar, volMesh> dsfType;
301 
302  forAll(rhoTrans_, fieldI)
303  {
304  dsfType& rhoT = rhoTrans_[fieldI];
305  const dsfType& rhoT0 = cloudOldTime.rhoTrans()[fieldI];
306  this->relax(rhoT, rhoT0, "rho");
307  }
308 }
309 
310 
311 template<class CloudType>
313 {
314  CloudType::scaleSources();
315 
316  typedef DimensionedField<scalar, volMesh> dsfType;
317 
318  forAll(rhoTrans_, fieldI)
319  {
320  dsfType& rhoT = rhoTrans_[fieldI];
321  this->scale(rhoT, "rho");
322  }
323 }
324 
325 
326 template<class CloudType>
328 {
329  if (this->solution().canEvolve())
330  {
331  typename parcelType::template
332  TrackingData<ReactingCloud<CloudType> > td(*this);
333 
334  this->solve(td);
335  }
336 }
337 
338 
339 template<class CloudType>
341 {
342  typedef typename particle::TrackingData<ReactingCloud<CloudType> > tdType;
343 
344  tdType td(*this);
345 
346  Cloud<parcelType>::template autoMap<tdType>(td, mapper);
347 
348  this->updateMesh();
349 }
350 
351 
352 template<class CloudType>
354 {
355  CloudType::info();
356 
357  this->phaseChange().info(Info);
358 }
359 
360 
361 template<class CloudType>
363 {
364  if (this->size())
365  {
366  CloudType::particleType::writeFields(*this, this->composition());
367  }
368 }
369 
370 
371 // ************************************************************************* //
UEqn relax()
const basicSpecieMixture & carrier() const
Return reference to the gaseous components.
Definition: SLGThermo.C:108
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...
rhoEqn solve()
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.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
void cloudReset(ReactingCloud< CloudType > &c)
Reset state of cloud.
Definition: ReactingCloud.C:84
Virtual abstract base class for templated ReactingCloud.
Definition: reactingCloud.H:48
A class for handling words, derived from string.
Definition: word.H:59
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
void info()
Print cloud information.
DSMCCloud< dsmcParcel > CloudType
messageStream Info
PtrList< DimensionedField< scalar, volMesh > > rhoTrans_
Mass transfer fields - one per carrier phase specie.
dynamicFvMesh & mesh
void resetSourceTerms()
Reset the cloud source terms.
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:48
Templated phase change model class.
Definition: ReactingCloud.H:55
autoPtr< BasicCompressibleTurbulenceModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleTurbulenceModel::transportModel &transport, const word &propertiesName)
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:217
void evolve()
Evolve the cloud.
Thermo package for (S)olids (L)iquids and (G)ases Takes reference to thermo package, and provides:
Definition: SLGThermo.H:62
DimensionedField< scalar, volMesh > & rhoTrans(const label i)
Mass.
Templated base class for reacting cloud.
Definition: ReactingCloud.H:62
autoPtr< PhaseChangeModel< ReactingCloud< CloudType > > > phaseChangeModel_
Reacting phase change model.
static const char nl
Definition: Ostream.H:260
basicMultiComponentMixture & composition
Definition: createFields.H:35
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
void readFields(const Mesh &mesh, const IOobjectList &objects, PtrList< GeoField > &fields)
void restoreState()
Reset the current cloud to the previously stored state.
void checkParcelProperties(parcelType &parcel, const scalar lagrangianDt, const bool fullyDescribed)
Check parcel properties.
void setParcelThermoProperties(parcelType &parcel, const scalar lagrangianDt)
Set parcel thermo properties.
#define forAll(list, i)
Definition: UList.H:421
virtual void writeFields() const
Write the field data for the cloud.
parcelType::constantProperties constProps_
Parcel constant properties.
autoPtr< CompositionModel< ReactingCloud< CloudType > > > compositionModel_
Reacting composition model.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
errorManip< error > abort(error &err)
Definition: errorManip.H:131
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
void scaleSources()
Apply scaling to (transient) cloud sources.
error FatalError
Templated reacting parcel composition model class Consists of carrier species (via thermo package)...
Definition: ReactingCloud.H:52
virtual ~ReactingCloud()
Destructor.
Base cloud calls templated on particle type.
Definition: Cloud.H:52
void setModels()
Set cloud sub-models.
Definition: ReactingCloud.C:34
psiReactionThermo & thermo
Definition: createFields.H:32
const speciesTable & species() const
Return the table of species.
void checkSuppliedComposition(const scalarField &YSupplied, const scalarField &Y, const word &YName)
Check that size of a composition field is valid.
Definition: ReactingCloud.C:58
void storeState()
Store the current cloud state.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
void relaxSources(const ReactingCloud< CloudType > &cloudOldTime)
Apply relaxation to (steady state) cloud sources.
word timeName
Definition: getTimeIndex.H:3