Go to the documentation of this file.
35 template<
class ChemistryModel>
38 typename ChemistryModel::reactionThermo&
thermo
42 coeffsDict_(this->subDict(
"EulerImplicitCoeffs")),
43 cTauChem_(coeffsDict_.get<scalar>(
"cTauChem")),
44 eqRateLimiter_(coeffsDict_.lookup(
"equilibriumRateLimiter")),
51 template<
class ChemistryModel>
58 template<
class ChemistryModel>
73 this->reactions_[index];
77 const label si =
R.lhs()[
s].index;
78 const scalar sl =
R.lhs()[
s].stoichCoeff;
79 RR[si][rRef] -= sl*pr*corr;
80 RR[si][lRef] += sl*pf*corr;
85 const label si =
R.rhs()[
s].index;
86 const scalar sr =
R.rhs()[
s].stoichCoeff;
87 RR[si][lRef] -= sr*pf*corr;
88 RR[si][rRef] += sr*pr*corr;
93 template<
class ChemistryModel>
112 const scalar cTot =
sum(
c);
113 typename ChemistryModel::thermoType
mixture
115 (this->specieThermo_[0].W()*
c[0])*this->specieThermo_[0]
119 mixture += (this->specieThermo_[i].W()*
c[i])*this->specieThermo_[i];
122 const scalar deltaTEst =
min(deltaT, subDeltaT);
124 forAll(this->reactions(), i)
126 scalar pf, cf, pr, cr;
129 const scalar omegai =
130 this->omegaI(i,
c,
T,
p, pf, cf, lRef, pr, cr, rRef);
137 corr = 1/(1 + pr*deltaTEst);
141 corr = 1/(1 + pf*deltaTEst);
145 updateRRInReactionI(i, pr, pf, corr, lRef, rRef,
p,
T,
RR);
161 tMin =
min(tMin, -(
c[i] + SMALL)/d);
166 const scalar cm =
max(cTot -
c[i], 1
e-5);
167 tMin =
min(tMin, cm/d);
171 subDeltaT = cTauChem_*tMin;
172 deltaT =
min(deltaT, subDeltaT);
177 RR(i, i) += 1/deltaT;
178 RR.source()[i] =
c[i]/deltaT;
191 mixture = (this->specieThermo_[0].W()*
c[0])*this->specieThermo_[0];
194 mixture += (this->specieThermo_[i].W()*
c[i])*this->specieThermo_[i];
const scalar RR
Universal gas constant: default in [J/(kmol K)].
A simple square matrix solver with scalar coefficients.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.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))
An abstract base class for solving chemistry.
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
EulerImplicit(typename ChemistryModel::reactionThermo &thermo)
Construct from thermo.
#define forAll(list, i)
Loop across all elements in list.
virtual void solve(scalarField &c, scalar &T, scalar &p, scalar &deltaT, scalar &subDeltaT) const
Update the concentrations and return the chemical time.
#define R(A, B, C, D, E, F, K, M)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Macros for easy insertion into run-time selection tables.
virtual ~EulerImplicit()
Destructor.
const dimensionedScalar e
Elementary charge.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
const dimensionedScalar c
Speed of light in a vacuum.
Simple extension of ReactionThermo to handle reaction kinetics in addition to the equilibrium thermod...
Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n"<< endl;autoPtr< temperaturePhaseChangeTwoPhaseMixture > mixture
void updateRRInReactionI(const label index, const scalar pr, const scalar pf, const scalar corr, const label lRef, const label rRef, const scalar p, const scalar T, simpleMatrix< scalar > &RR) const