All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
chemistryModel.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-2017 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 "chemistryModel.H"
27 #include "reactingMixture.H"
28 #include "UniformField.H"
30 
31 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
32 
33 template<class CompType, class ThermoType>
35 (
36  const fvMesh& mesh,
37  const word& phaseName
38 )
39 :
40  CompType(mesh, phaseName),
41  ODESystem(),
42  Y_(this->thermo().composition().Y()),
43  reactions_
44  (
45  dynamic_cast<const reactingMixture<ThermoType>&>(this->thermo())
46  ),
47  specieThermo_
48  (
49  dynamic_cast<const reactingMixture<ThermoType>&>
50  (this->thermo()).speciesData()
51  ),
52 
53  nSpecie_(Y_.size()),
54  nReaction_(reactions_.size()),
55  Treact_(CompType::template lookupOrDefault<scalar>("Treact", 0.0)),
56  RR_(nSpecie_),
57  c_(nSpecie_),
58  dcdt_(nSpecie_)
59 {
60  // Create the fields for the chemistry sources
61  forAll(RR_, fieldi)
62  {
63  RR_.set
64  (
65  fieldi,
67  (
68  IOobject
69  (
70  "RR." + Y_[fieldi].name(),
71  mesh.time().timeName(),
72  mesh,
73  IOobject::NO_READ,
74  IOobject::NO_WRITE
75  ),
76  mesh,
78  )
79  );
80  }
81 
82  Info<< "chemistryModel: Number of species = " << nSpecie_
83  << " and reactions = " << nReaction_ << endl;
84 }
85 
86 
87 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
88 
89 template<class CompType, class ThermoType>
91 {}
92 
93 
94 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
95 
96 template<class CompType, class ThermoType>
98 (
99  const scalarField& c,
100  const scalar T,
101  const scalar p,
102  scalarField& dcdt
103 ) const
104 {
105  scalar pf, cf, pr, cr;
106  label lRef, rRef;
107 
108  dcdt = Zero;
109 
110  forAll(reactions_, i)
111  {
112  const Reaction<ThermoType>& R = reactions_[i];
113 
114  scalar omegai = omega
115  (
116  R, c, T, p, pf, cf, lRef, pr, cr, rRef
117  );
118 
119  forAll(R.lhs(), s)
120  {
121  const label si = R.lhs()[s].index;
122  const scalar sl = R.lhs()[s].stoichCoeff;
123  dcdt[si] -= sl*omegai;
124  }
125 
126  forAll(R.rhs(), s)
127  {
128  const label si = R.rhs()[s].index;
129  const scalar sr = R.rhs()[s].stoichCoeff;
130  dcdt[si] += sr*omegai;
131  }
132  }
133 }
134 
135 
136 template<class CompType, class ThermoType>
138 (
139  const label index,
140  const scalarField& c,
141  const scalar T,
142  const scalar p,
143  scalar& pf,
144  scalar& cf,
145  label& lRef,
146  scalar& pr,
147  scalar& cr,
148  label& rRef
149 ) const
150 {
151  const Reaction<ThermoType>& R = reactions_[index];
152  scalar w = omega(R, c, T, p, pf, cf, lRef, pr, cr, rRef);
153  return(w);
154 }
155 
156 
157 template<class CompType, class ThermoType>
159 (
160  const Reaction<ThermoType>& R,
161  const scalarField& c,
162  const scalar T,
163  const scalar p,
164  scalar& pf,
165  scalar& cf,
166  label& lRef,
167  scalar& pr,
168  scalar& cr,
169  label& rRef
170 ) const
171 {
172  const scalar kf = R.kf(p, T, c);
173  const scalar kr = R.kr(kf, p, T, c);
174 
175  pf = 1.0;
176  pr = 1.0;
177 
178  const label Nl = R.lhs().size();
179  const label Nr = R.rhs().size();
180 
181  label slRef = 0;
182  lRef = R.lhs()[slRef].index;
183 
184  pf = kf;
185  for (label s = 1; s < Nl; s++)
186  {
187  const label si = R.lhs()[s].index;
188 
189  if (c[si] < c[lRef])
190  {
191  const scalar exp = R.lhs()[slRef].exponent;
192  pf *= pow(max(0.0, c[lRef]), exp);
193  lRef = si;
194  slRef = s;
195  }
196  else
197  {
198  const scalar exp = R.lhs()[s].exponent;
199  pf *= pow(max(0.0, c[si]), exp);
200  }
201  }
202  cf = max(0.0, c[lRef]);
203 
204  {
205  const scalar exp = R.lhs()[slRef].exponent;
206  if (exp < 1.0)
207  {
208  if (cf > SMALL)
209  {
210  pf *= pow(cf, exp - 1.0);
211  }
212  else
213  {
214  pf = 0.0;
215  }
216  }
217  else
218  {
219  pf *= pow(cf, exp - 1.0);
220  }
221  }
222 
223  label srRef = 0;
224  rRef = R.rhs()[srRef].index;
225 
226  // Find the matrix element and element position for the rhs
227  pr = kr;
228  for (label s = 1; s < Nr; s++)
229  {
230  const label si = R.rhs()[s].index;
231  if (c[si] < c[rRef])
232  {
233  const scalar exp = R.rhs()[srRef].exponent;
234  pr *= pow(max(0.0, c[rRef]), exp);
235  rRef = si;
236  srRef = s;
237  }
238  else
239  {
240  const scalar exp = R.rhs()[s].exponent;
241  pr *= pow(max(0.0, c[si]), exp);
242  }
243  }
244  cr = max(0.0, c[rRef]);
245 
246  {
247  const scalar exp = R.rhs()[srRef].exponent;
248  if (exp < 1.0)
249  {
250  if (cr>SMALL)
251  {
252  pr *= pow(cr, exp - 1.0);
253  }
254  else
255  {
256  pr = 0.0;
257  }
258  }
259  else
260  {
261  pr *= pow(cr, exp - 1.0);
262  }
263  }
264 
265  return pf*cf - pr*cr;
266 }
267 
268 
269 template<class CompType, class ThermoType>
271 (
272  const scalar time,
273  const scalarField& c,
274  scalarField& dcdt
275 ) const
276 {
277  const scalar T = c[nSpecie_];
278  const scalar p = c[nSpecie_ + 1];
279 
280  for (label i = 0; i < nSpecie_; i++)
281  {
282  c_[i] = max(0.0, c[i]);
283  }
284 
285  omega(c_, T, p, dcdt);
286 
287  // Constant pressure
288  // dT/dt = ...
289  scalar rho = 0.0;
290  scalar cSum = 0.0;
291  for (label i = 0; i < nSpecie_; i++)
292  {
293  const scalar W = specieThermo_[i].W();
294  cSum += c_[i];
295  rho += W*c_[i];
296  }
297  scalar cp = 0.0;
298  for (label i=0; i<nSpecie_; i++)
299  {
300  cp += c_[i]*specieThermo_[i].cp(p, T);
301  }
302  cp /= rho;
303 
304  scalar dT = 0.0;
305  for (label i = 0; i < nSpecie_; i++)
306  {
307  const scalar hi = specieThermo_[i].ha(p, T);
308  dT += hi*dcdt[i];
309  }
310  dT /= rho*cp;
311 
312  dcdt[nSpecie_] = -dT;
313 
314  // dp/dt = ...
315  dcdt[nSpecie_ + 1] = 0.0;
316 }
317 
318 
319 template<class CompType, class ThermoType>
321 (
322  const scalar t,
323  const scalarField& c,
324  scalarField& dcdt,
325  scalarSquareMatrix& dfdc
326 ) const
327 {
328  const scalar T = c[nSpecie_];
329  const scalar p = c[nSpecie_ + 1];
330 
331  forAll(c_, i)
332  {
333  c_[i] = max(c[i], 0.0);
334  }
335 
336  dfdc = Zero;
337 
338  // Length of the first argument must be nSpecie_
339  omega(c_, T, p, dcdt);
340 
341  forAll(reactions_, ri)
342  {
343  const Reaction<ThermoType>& R = reactions_[ri];
344 
345  const scalar kf0 = R.kf(p, T, c_);
346  const scalar kr0 = R.kr(kf0, p, T, c_);
347 
348  forAll(R.lhs(), j)
349  {
350  const label sj = R.lhs()[j].index;
351  scalar kf = kf0;
352  forAll(R.lhs(), i)
353  {
354  const label si = R.lhs()[i].index;
355  const scalar el = R.lhs()[i].exponent;
356  if (i == j)
357  {
358  if (el < 1.0)
359  {
360  if (c_[si] > SMALL)
361  {
362  kf *= el*pow(c_[si] + VSMALL, el - 1.0);
363  }
364  else
365  {
366  kf = 0.0;
367  }
368  }
369  else
370  {
371  kf *= el*pow(c_[si], el - 1.0);
372  }
373  }
374  else
375  {
376  kf *= pow(c_[si], el);
377  }
378  }
379 
380  forAll(R.lhs(), i)
381  {
382  const label si = R.lhs()[i].index;
383  const scalar sl = R.lhs()[i].stoichCoeff;
384  dfdc(si, sj) -= sl*kf;
385  }
386  forAll(R.rhs(), i)
387  {
388  const label si = R.rhs()[i].index;
389  const scalar sr = R.rhs()[i].stoichCoeff;
390  dfdc(si, sj) += sr*kf;
391  }
392  }
393 
394  forAll(R.rhs(), j)
395  {
396  const label sj = R.rhs()[j].index;
397  scalar kr = kr0;
398  forAll(R.rhs(), i)
399  {
400  const label si = R.rhs()[i].index;
401  const scalar er = R.rhs()[i].exponent;
402  if (i == j)
403  {
404  if (er < 1.0)
405  {
406  if (c_[si] > SMALL)
407  {
408  kr *= er*pow(c_[si] + VSMALL, er - 1.0);
409  }
410  else
411  {
412  kr = 0.0;
413  }
414  }
415  else
416  {
417  kr *= er*pow(c_[si], er - 1.0);
418  }
419  }
420  else
421  {
422  kr *= pow(c_[si], er);
423  }
424  }
425 
426  forAll(R.lhs(), i)
427  {
428  const label si = R.lhs()[i].index;
429  const scalar sl = R.lhs()[i].stoichCoeff;
430  dfdc(si, sj) += sl*kr;
431  }
432  forAll(R.rhs(), i)
433  {
434  const label si = R.rhs()[i].index;
435  const scalar sr = R.rhs()[i].stoichCoeff;
436  dfdc(si, sj) -= sr*kr;
437  }
438  }
439  }
440 
441  // Calculate the dcdT elements numerically
442  const scalar delta = 1.0e-3;
443 
444  omega(c_, T + delta, p, dcdt_);
445  for (label i=0; i<nSpecie_; i++)
446  {
447  dfdc(i, nSpecie_) = dcdt_[i];
448  }
449 
450  omega(c_, T - delta, p, dcdt_);
451  for (label i=0; i<nSpecie_; i++)
452  {
453  dfdc(i, nSpecie_) = 0.5*(dfdc(i, nSpecie_) - dcdt_[i])/delta;
454  }
455 
456  dfdc(nSpecie_, nSpecie_) = 0;
457  dfdc(nSpecie_ + 1, nSpecie_) = 0;
458 }
459 
460 
461 template<class CompType, class ThermoType>
464 {
466  (
467  new volScalarField
468  (
469  IOobject
470  (
471  "tc",
472  this->time().timeName(),
473  this->mesh(),
474  IOobject::NO_READ,
475  IOobject::NO_WRITE,
476  false
477  ),
478  this->mesh(),
479  dimensionedScalar("zero", dimTime, SMALL),
480  extrapolatedCalculatedFvPatchScalarField::typeName
481  )
482  );
483 
484  scalarField& tc = ttc.ref();
485 
486  tmp<volScalarField> trho(this->thermo().rho());
487  const scalarField& rho = trho();
488 
489  const scalarField& T = this->thermo().T();
490  const scalarField& p = this->thermo().p();
491 
492  const label nReaction = reactions_.size();
493 
494  scalar pf, cf, pr, cr;
495  label lRef, rRef;
496 
497  if (this->chemistry_)
498  {
499  forAll(rho, celli)
500  {
501  const scalar rhoi = rho[celli];
502  const scalar Ti = T[celli];
503  const scalar pi = p[celli];
504 
505  scalar cSum = 0.0;
506 
507  for (label i=0; i<nSpecie_; i++)
508  {
509  c_[i] = rhoi*Y_[i][celli]/specieThermo_[i].W();
510  cSum += c_[i];
511  }
512 
513  forAll(reactions_, i)
514  {
515  const Reaction<ThermoType>& R = reactions_[i];
516 
517  omega(R, c_, Ti, pi, pf, cf, lRef, pr, cr, rRef);
518 
519  forAll(R.rhs(), s)
520  {
521  tc[celli] += R.rhs()[s].stoichCoeff*pf*cf;
522  }
523  }
524 
525  tc[celli] = nReaction*cSum/tc[celli];
526  }
527  }
528 
529  ttc.ref().correctBoundaryConditions();
530 
531  return ttc;
532 }
533 
534 
535 template<class CompType, class ThermoType>
538 {
539  tmp<volScalarField> tQdot
540  (
541  new volScalarField
542  (
543  IOobject
544  (
545  "Qdot",
546  this->mesh_.time().timeName(),
547  this->mesh_,
548  IOobject::NO_READ,
549  IOobject::NO_WRITE,
550  false
551  ),
552  this->mesh_,
554  )
555  );
556 
557  if (this->chemistry_)
558  {
559  scalarField& Qdot = tQdot.ref();
560 
561  forAll(Y_, i)
562  {
563  forAll(Qdot, celli)
564  {
565  const scalar hi = specieThermo_[i].Hc();
566  Qdot[celli] -= hi*RR_[i][celli];
567  }
568  }
569  }
570 
571  return tQdot;
572 }
573 
574 
575 template<class CompType, class ThermoType>
578 (
579  const label ri,
580  const label si
581 ) const
582 {
583  scalar pf, cf, pr, cr;
584  label lRef, rRef;
585 
587  (
589  (
590  IOobject
591  (
592  "RR",
593  this->mesh().time().timeName(),
594  this->mesh(),
595  IOobject::NO_READ,
596  IOobject::NO_WRITE
597  ),
598  this->mesh(),
600  )
601  );
602 
604 
605  tmp<volScalarField> trho(this->thermo().rho());
606  const scalarField& rho = trho();
607 
608  const scalarField& T = this->thermo().T();
609  const scalarField& p = this->thermo().p();
610 
611  forAll(rho, celli)
612  {
613  const scalar rhoi = rho[celli];
614  const scalar Ti = T[celli];
615  const scalar pi = p[celli];
616 
617  for (label i=0; i<nSpecie_; i++)
618  {
619  const scalar Yi = Y_[i][celli];
620  c_[i] = rhoi*Yi/specieThermo_[i].W();
621  }
622 
623  const scalar w = omegaI
624  (
625  ri,
626  c_,
627  Ti,
628  pi,
629  pf,
630  cf,
631  lRef,
632  pr,
633  cr,
634  rRef
635  );
636 
637  RR[celli] = w*specieThermo_[si].W();
638  }
639 
640  return tRR;
641 }
642 
643 
644 template<class CompType, class ThermoType>
646 {
647  if (!this->chemistry_)
648  {
649  return;
650  }
651 
652  tmp<volScalarField> trho(this->thermo().rho());
653  const scalarField& rho = trho();
654 
655  const scalarField& T = this->thermo().T();
656  const scalarField& p = this->thermo().p();
657 
658  forAll(rho, celli)
659  {
660  const scalar rhoi = rho[celli];
661  const scalar Ti = T[celli];
662  const scalar pi = p[celli];
663 
664  for (label i=0; i<nSpecie_; i++)
665  {
666  const scalar Yi = Y_[i][celli];
667  c_[i] = rhoi*Yi/specieThermo_[i].W();
668  }
669 
670  omega(c_, Ti, pi, dcdt_);
671 
672  for (label i=0; i<nSpecie_; i++)
673  {
674  RR_[i][celli] = dcdt_[i]*specieThermo_[i].W();
675  }
676  }
677 }
678 
679 
680 template<class CompType, class ThermoType>
681 template<class DeltaTType>
683 (
684  const DeltaTType& deltaT
685 )
686 {
688 
689  scalar deltaTMin = GREAT;
690 
691  if (!this->chemistry_)
692  {
693  return deltaTMin;
694  }
695 
696  tmp<volScalarField> trho(this->thermo().rho());
697  const scalarField& rho = trho();
698 
699  const scalarField& T = this->thermo().T();
700  const scalarField& p = this->thermo().p();
701 
702  scalarField c0(nSpecie_);
703 
704  forAll(rho, celli)
705  {
706  scalar Ti = T[celli];
707 
708  if (Ti > Treact_)
709  {
710  const scalar rhoi = rho[celli];
711  scalar pi = p[celli];
712 
713  for (label i=0; i<nSpecie_; i++)
714  {
715  c_[i] = rhoi*Y_[i][celli]/specieThermo_[i].W();
716  c0[i] = c_[i];
717  }
718 
719  // Initialise time progress
720  scalar timeLeft = deltaT[celli];
721 
722  // Calculate the chemical source terms
723  while (timeLeft > SMALL)
724  {
725  scalar dt = timeLeft;
726  this->solve(c_, Ti, pi, dt, this->deltaTChem_[celli]);
727  timeLeft -= dt;
728  }
729 
730  deltaTMin = min(this->deltaTChem_[celli], deltaTMin);
731 
732  for (label i=0; i<nSpecie_; i++)
733  {
734  RR_[i][celli] =
735  (c_[i] - c0[i])*specieThermo_[i].W()/deltaT[celli];
736  }
737  }
738  else
739  {
740  for (label i=0; i<nSpecie_; i++)
741  {
742  RR_[i][celli] = 0;
743  }
744  }
745  }
746 
747  return deltaTMin;
748 }
749 
750 
751 template<class CompType, class ThermoType>
753 (
754  const scalar deltaT
755 )
756 {
757  // Don't allow the time-step to change more than a factor of 2
758  return min
759  (
761  2*deltaT
762  );
763 }
764 
765 
766 template<class CompType, class ThermoType>
768 (
769  const scalarField& deltaT
770 )
771 {
772  return this->solve<scalarField>(deltaT);
773 }
774 
775 
776 // ************************************************************************* //
virtual tmp< volScalarField::Internal > calculateRR(const label reactionI, const label speciei) const
Return reaction rate of the speciei in reactionI.
scalar delta
virtual scalar kr(const scalar kfwd, const scalar p, const scalar T, const scalarField &c) const
Reverse rate constant from the given forward rate constant.
Definition: Reaction.C:404
virtual scalar kf(const scalar p, const scalar T, const scalarField &c) const
Forward rate constant.
Definition: Reaction.C:392
A class representing the concept of a uniform field which stores only the single value and providing ...
Definition: UniformField.H:47
virtual scalar omegaI(label iReaction, const scalarField &c, const scalar T, const scalar p, scalar &pf, scalar &cf, label &lRef, scalar &pr, scalar &cr, label &rRef) const
Return the reaction rate for iReaction and the reference.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
Foam::reactingMixture.
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 List< specieCoeffs > & lhs() const
Definition: ReactionI.H:51
Abstract base class for the systems of ordinary differential equations.
Definition: ODESystem.H:46
virtual tmp< volScalarField > tc() const
Return the chemical time scale.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
virtual void jacobian(const scalar t, const scalarField &c, scalarField &dcdt, scalarSquareMatrix &dfdc) const
Calculate the Jacobian of the system.
basicMultiComponentMixture & composition
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
virtual ~chemistryModel()
Destructor.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
tmp< volScalarField > trho
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:644
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
Simple extension of ReactionThermo to handle reaction kinetics in addition to the equilibrium thermod...
Definition: Reaction.H:53
scalar Qdot
Definition: solveChemistry.H:2
psiReactionThermo & thermo
Definition: createFields.H:31
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
dynamicFvMesh & mesh
dimensionedScalar exp(const dimensionedScalar &ds)
A class for handling words, derived from string.
Definition: word.H:59
rhoEqn solve()
word timeName
Definition: getTimeIndex.H:3
static const zero Zero
Definition: zero.H:91
const volScalarField & cp
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(psi *p+alphal *rhol0+((alphav *psiv+alphal *psil) - psi) *pSat, rhoMin);# 1 "/home/ubuntu/OpenFOAM-5.0/applications/solvers/multiphase/cavitatingFoam/alphavPsi.H" 1{ alphav=max(min((rho - rholSat)/(rhovSat - rholSat), scalar(1)), scalar(0));alphal=1.0 - alphav;Info<< "max-min alphav: "<< max(alphav).value()<< " "<< min(alphav).value()<< endl;psiModel-> correct()
Definition: pEqn.H:72
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const dimensionSet dimEnergy
virtual void calculate()
Calculates the reaction rates.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Extends base chemistry model by adding a thermo package, and ODE functions. Introduces chemistry equa...
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
PtrList< volScalarField > & Y
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...
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
virtual void derivatives(const scalar t, const scalarField &c, scalarField &dcdt) const
Calculate the derivatives in dydx.
messageStream Info
const scalar RR
Universal gas constant (default in [J/(kmol K)])
virtual tmp< volScalarField > Qdot() const
Return the heat release rate [kg/m/s3].
virtual void omega(const scalarField &c, const scalar T, const scalar p, scalarField &dcdt) const
dc/dt = omega, rate of change in concentration, for each species
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
A class for managing temporary objects.
Definition: PtrList.H:53
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
const List< specieCoeffs > & rhs() const
Definition: ReactionI.H:59