thermophysicalModels::reaction: New concrete base class for Reaction

to provide reaction specie coefficients without the need for a thermodynamics
model.
This commit is contained in:
Henry Weller
2020-10-29 22:21:58 +00:00
parent ab12f38c2a
commit e8fba9844a
15 changed files with 455 additions and 497 deletions

View File

@ -65,36 +65,36 @@ Foam::radiationModels::sootModels::mixtureFraction<ThermoType>::mixtureFraction
combustionModel::combustionPropertiesName
);
const multiComponentMixture<ThermoType>& mixture = combustion.mixture();
const basicSpecieMixture& mixture = combustion.mixture();
const Reaction<ThermoType>& reaction = combustion.reaction();
const reaction& singleReaction = combustion.singleReaction();
scalar totalMol = 0;
forAll(reaction.rhs(), i)
forAll(singleReaction.rhs(), i)
{
const scalar stoichCoeff = reaction.rhs()[i].stoichCoeff;
const scalar stoichCoeff = singleReaction.rhs()[i].stoichCoeff;
totalMol += mag(stoichCoeff);
}
totalMol += nuSoot_;
scalarList Xi(reaction.rhs().size());
scalarList Xi(singleReaction.rhs().size());
scalar Wm = 0;
forAll(reaction.rhs(), i)
forAll(singleReaction.rhs(), i)
{
const label speciei = reaction.rhs()[i].index;
const scalar stoichCoeff = reaction.rhs()[i].stoichCoeff;
const label speciei = singleReaction.rhs()[i].index;
const scalar stoichCoeff = singleReaction.rhs()[i].stoichCoeff;
Xi[i] = mag(stoichCoeff)/totalMol;
Wm += Xi[i]*mixture.specieThermos()[speciei].W();
Wm += Xi[i]*mixture.Wi(speciei);
}
scalarList Yprod0(mixture.species().size(), 0.0);
forAll(reaction.rhs(), i)
forAll(singleReaction.rhs(), i)
{
const label speciei = reaction.rhs()[i].index;
Yprod0[speciei] = mixture.specieThermos()[speciei].W()/Wm*Xi[i];
const label speciei = singleReaction.rhs()[i].index;
Yprod0[speciei] = mixture.Wi(speciei)/Wm*Xi[i];
}
const scalar XSoot = nuSoot_/totalMol;
@ -106,7 +106,7 @@ Foam::radiationModels::sootModels::mixtureFraction<ThermoType>::mixtureFraction
if (mappingFieldName_ == "none")
{
const label index = reaction.rhs()[0].index;
const label index = singleReaction.rhs()[0].index;
mappingFieldName_ = mixture.Y(index).name();
}

View File

@ -38,83 +38,79 @@ namespace combustionModels
template<class ThermoType>
void singleStepCombustion<ThermoType>::calculateqFuel()
{
const Reaction<ThermoType>& reaction = reaction_();
const scalar Wu = mixture_.specieThermos()[fuelIndex_].W();
const scalar Wu = mixture_.Wi(fuelIndex_);
forAll(reaction.lhs(), i)
forAll(reaction_.lhs(), i)
{
const label speciei = reaction.lhs()[i].index;
const scalar stoichCoeff = reaction.lhs()[i].stoichCoeff;
const label speciei = reaction_.lhs()[i].index;
const scalar stoichCoeff = reaction_.lhs()[i].stoichCoeff;
specieStoichCoeffs_[speciei] = -stoichCoeff;
qFuel_.value() += mixture_.specieThermos()[speciei].hc()*stoichCoeff/Wu;
qFuel_.value() +=
mixture_.Hf(speciei)*mixture_.Wi(speciei)*stoichCoeff/Wu;
}
forAll(reaction.rhs(), i)
forAll(reaction_.rhs(), i)
{
const label speciei = reaction.rhs()[i].index;
const scalar stoichCoeff = reaction.rhs()[i].stoichCoeff;
const label speciei = reaction_.rhs()[i].index;
const scalar stoichCoeff = reaction_.rhs()[i].stoichCoeff;
specieStoichCoeffs_[speciei] = stoichCoeff;
qFuel_.value() -= mixture_.specieThermos()[speciei].hc()*stoichCoeff/Wu;
qFuel_.value() -=
mixture_.Hf(speciei)*mixture_.Wi(speciei)*stoichCoeff/Wu;
specieProd_[speciei] = -1;
}
Info << "Fuel heat of combustion :" << qFuel_.value() << endl;
Info << "Fuel heat of combustion: " << qFuel_.value() << endl;
}
template<class ThermoType>
void singleStepCombustion<ThermoType>:: massAndAirStoichRatios()
void singleStepCombustion<ThermoType>::massAndAirStoichRatios()
{
const label O2Index = mixture_.species()["O2"];
const scalar Wu = mixture_.specieThermos()[fuelIndex_].W();
const scalar Wu = mixture_.Wi(fuelIndex_);
stoicRatio_ =
(mixture_.specieThermos()[mixture_.defaultSpecie()].W()
* specieStoichCoeffs_[mixture_.defaultSpecie()]
+ mixture_.specieThermos()[O2Index].W()
* mag(specieStoichCoeffs_[O2Index]))
/ (Wu*mag(specieStoichCoeffs_[fuelIndex_]));
(
mixture_.Wi(mixture_.defaultSpecie())
*specieStoichCoeffs_[mixture_.defaultSpecie()]
+ mixture_.Wi(O2Index)*mag(specieStoichCoeffs_[O2Index])
)/(Wu*mag(specieStoichCoeffs_[fuelIndex_]));
s_ =
(mixture_.specieThermos()[O2Index].W()
* mag(specieStoichCoeffs_[O2Index]))
/ (Wu*mag(specieStoichCoeffs_[fuelIndex_]));
s_ = mixture_.Wi(O2Index)*mag(specieStoichCoeffs_[O2Index])
/(Wu*mag(specieStoichCoeffs_[fuelIndex_]));
Info << "stoichiometric air-fuel ratio :" << stoicRatio_.value() << endl;
Info << "stoichiometric oxygen-fuel ratio :" << s_.value() << endl;
Info << "stoichiometric air-fuel ratio: " << stoicRatio_.value() << endl;
Info << "stoichiometric oxygen-fuel ratio: " << s_.value() << endl;
}
template<class ThermoType>
void singleStepCombustion<ThermoType>:: calculateMaxProducts()
void singleStepCombustion<ThermoType>::calculateMaxProducts()
{
const Reaction<ThermoType>& reaction = reaction_();
scalar Wm = 0.0;
scalar totalMol = 0.0;
forAll(reaction.rhs(), i)
forAll(reaction_.rhs(), i)
{
label speciei = reaction.rhs()[i].index;
label speciei = reaction_.rhs()[i].index;
totalMol += mag(specieStoichCoeffs_[speciei]);
}
scalarList Xi(reaction.rhs().size());
scalarList Xi(reaction_.rhs().size());
forAll(reaction.rhs(), i)
forAll(reaction_.rhs(), i)
{
const label speciei = reaction.rhs()[i].index;
const label speciei = reaction_.rhs()[i].index;
Xi[i] = mag(specieStoichCoeffs_[speciei])/totalMol;
Wm += Xi[i]*mixture_.specieThermos()[speciei].W();
Wm += Xi[i]*mixture_.Wi(speciei);
}
forAll(reaction.rhs(), i)
forAll(reaction_.rhs(), i)
{
const label speciei = reaction.rhs()[i].index;
Yprod0_[speciei] = mixture_.specieThermos()[speciei].W()/Wm*Xi[i];
const label speciei = reaction_.rhs()[i].index;
Yprod0_[speciei] = mixture_.Wi(speciei)/Wm*Xi[i];
}
Info << "Maximum products mass concentrations:" << nl;
Info << "Maximum products mass concentrations: " << nl;
forAll(Yprod0_, i)
{
if (Yprod0_[i] > 0)
@ -127,10 +123,8 @@ void singleStepCombustion<ThermoType>:: calculateMaxProducts()
forAll(specieStoichCoeffs_, i)
{
specieStoichCoeffs_[i] =
specieStoichCoeffs_[i]
* mixture_.specieThermos()[i].W()
/ (mixture_.specieThermos()[fuelIndex_].W()
* mag(specieStoichCoeffs_[fuelIndex_]));
specieStoichCoeffs_[i]*mixture_.Wi(i)
/(mixture_.Wi(fuelIndex_)*mag(specieStoichCoeffs_[fuelIndex_]));
}
}
@ -147,19 +141,8 @@ singleStepCombustion<ThermoType>::singleStepCombustion
)
:
combustionModel(modelType, thermo, turb, combustionProperties),
mixture_
(
dynamic_cast<const multiComponentMixture<ThermoType>&>(this->thermo())
),
reaction_
(
Reaction<ThermoType>::New
(
mixture_.species(),
mixture_.specieThermos(),
this->subDict("reaction")
)
),
mixture_(dynamic_cast<const basicSpecieMixture&>(this->thermo())),
reaction_(mixture_.species(), this->subDict("reaction")),
stoicRatio_(dimensionedScalar("stoicRatio", dimless, 0)),
s_(dimensionedScalar("s", dimless, 0)),
qFuel_(dimensionedScalar("qFuel", sqr(dimVelocity), 0)),
@ -287,16 +270,14 @@ bool singleStepCombustion<ThermoType>::read()
template<class ThermoType>
void singleStepCombustion<ThermoType>::fresCorrect()
{
const Reaction<ThermoType>& reaction = reaction_();
const label O2Index = mixture_.species()["O2"];
const volScalarField& YFuel = mixture_.Y()[fuelIndex_];
const volScalarField& YO2 = mixture_.Y()[O2Index];
// reactants
forAll(reaction.lhs(), i)
forAll(reaction_.lhs(), i)
{
const label speciei = reaction.lhs()[i].index;
const label speciei = reaction_.lhs()[i].index;
if (speciei == fuelIndex_)
{
fres_[speciei] = max(YFuel - YO2/s_, scalar(0));
@ -308,9 +289,9 @@ void singleStepCombustion<ThermoType>::fresCorrect()
}
// products
forAll(reaction.rhs(), i)
forAll(reaction_.rhs(), i)
{
const label speciei = reaction.rhs()[i].index;
const label speciei = reaction_.rhs()[i].index;
if (speciei != mixture_.defaultSpecie())
{
forAll(fres_[speciei], celli)

View File

@ -25,7 +25,7 @@ Class
Foam::combustionModels::singleStepCombustion
Description
Base class for combustion models using multiComponentMixture.
Base class for combustion models using basicSpecieMixture.
SourceFiles
singleStepCombustion.C
@ -36,8 +36,8 @@ SourceFiles
#define singleStepCombustion_H
#include "combustionModel.H"
#include "multiComponentMixture.H"
#include "ReactionList.H"
#include "basicSpecieMixture.H"
#include "reaction.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -60,10 +60,10 @@ protected:
// Protected data
//- Reference to the mixture
const multiComponentMixture<ThermoType>& mixture_;
const basicSpecieMixture& mixture_;
//- The single-step reaction
autoPtr<Reaction<ThermoType>> reaction_;
reaction reaction_;
//- Stoichiometric air-fuel mass ratio
dimensionedScalar stoicRatio_;
@ -134,10 +134,10 @@ public:
// Access functions
//- Return the mixture
const multiComponentMixture<ThermoType>& mixture() const;
const basicSpecieMixture& mixture() const;
//- Return the single step reaction
inline const Reaction<ThermoType>& reaction() const;
inline const reaction& singleReaction() const;
//- Return the stoichiometric air-fuel mass ratio
inline const dimensionedScalar& stoicRatio() const;

View File

@ -35,7 +35,7 @@ namespace combustionModels
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
template<class ThermoType>
const multiComponentMixture<ThermoType>&
const Foam::basicSpecieMixture&
singleStepCombustion<ThermoType>::mixture() const
{
return mixture_;
@ -43,10 +43,10 @@ singleStepCombustion<ThermoType>::mixture() const
template<class ThermoType>
inline const Foam::Reaction<ThermoType>&
singleStepCombustion<ThermoType>::reaction() const
inline const Foam::reaction&
singleStepCombustion<ThermoType>::singleReaction() const
{
return reaction_();
return reaction_;
}

View File

@ -26,7 +26,6 @@ License
#include "makeReaction.H"
#include "ArrheniusReactionRate.H"
#include "infiniteReactionRate.H"
#include "LandauTellerReactionRate.H"
#include "thirdBodyArrheniusReactionRate.H"
@ -71,10 +70,6 @@ namespace Foam
forCommonLiquids(makeIRNReactions, ArrheniusReactionRate);
forPolynomials(makeIRNReactions, ArrheniusReactionRate);
forCommonGases(makeIRNReactions, infiniteReactionRate);
forCommonLiquids(makeIRNReactions, infiniteReactionRate);
forPolynomials(makeIRNReactions, infiniteReactionRate);
forCommonGases(makeIRNReactions, LandauTellerReactionRate);
forCommonLiquids(makeIRNReactions, LandauTellerReactionRate);
forPolynomials(makeIRNReactions, LandauTellerReactionRate);

View File

@ -1,6 +1,7 @@
atomicWeights/atomicWeights.C
specie/specie.C
reaction/specieCoeffs/specieCoeffs.C
reaction/reaction/reaction.C
thermophysicalFunctions/thermophysicalFunction/thermophysicalFunction.C

View File

@ -25,10 +25,8 @@ License
#include "Reaction.H"
// * * * * * * * * * * * * * * * * Static Data * * * * * * * * * * * * * * * //
template<class ReactionThermo>
Foam::label Foam::Reaction<ReactionThermo>::nUnNamedReactions(0);
// * * * * * * * * * * * * * * * * Static Data * * * * * * * * * * * * * * * //
template<class ReactionThermo>
Foam::scalar Foam::Reaction<ReactionThermo>::TlowDefault(0);
@ -39,13 +37,6 @@ Foam::scalar Foam::Reaction<ReactionThermo>::ThighDefault(great);
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
template<class ReactionThermo>
Foam::label Foam::Reaction<ReactionThermo>::getNewReactionID()
{
return nUnNamedReactions++;
}
template<class ReactionThermo>
void Foam::Reaction<ReactionThermo>::setThermo
(
@ -54,32 +45,32 @@ void Foam::Reaction<ReactionThermo>::setThermo
{
typename ReactionThermo::thermoType rhsThermo
(
rhs_[0].stoichCoeff
*(*thermoDatabase[species_[rhs_[0].index]]).W()
*(*thermoDatabase[species_[rhs_[0].index]])
rhs()[0].stoichCoeff
*(*thermoDatabase[species()[rhs()[0].index]]).W()
*(*thermoDatabase[species()[rhs()[0].index]])
);
for (label i=1; i<rhs_.size(); ++i)
for (label i=1; i<rhs().size(); ++i)
{
rhsThermo +=
rhs_[i].stoichCoeff
*(*thermoDatabase[species_[rhs_[i].index]]).W()
*(*thermoDatabase[species_[rhs_[i].index]]);
rhs()[i].stoichCoeff
*(*thermoDatabase[species()[rhs()[i].index]]).W()
*(*thermoDatabase[species()[rhs()[i].index]]);
}
typename ReactionThermo::thermoType lhsThermo
(
lhs_[0].stoichCoeff
*(*thermoDatabase[species_[lhs_[0].index]]).W()
*(*thermoDatabase[species_[lhs_[0].index]])
lhs()[0].stoichCoeff
*(*thermoDatabase[species()[lhs()[0].index]]).W()
*(*thermoDatabase[species()[lhs()[0].index]])
);
for (label i=1; i<lhs_.size(); ++i)
for (label i=1; i<lhs().size(); ++i)
{
lhsThermo +=
lhs_[i].stoichCoeff
*(*thermoDatabase[species_[lhs_[i].index]]).W()
*(*thermoDatabase[species_[lhs_[i].index]]);
lhs()[i].stoichCoeff
*(*thermoDatabase[species()[lhs()[i].index]]).W()
*(*thermoDatabase[species()[lhs()[i].index]]);
}
// Check for mass imbalance in the reaction
@ -108,13 +99,10 @@ Foam::Reaction<ReactionThermo>::Reaction
const HashPtrTable<ReactionThermo>& thermoDatabase
)
:
reaction(species, lhs, rhs),
ReactionThermo::thermoType(*thermoDatabase[species[0]]),
name_("un-named-reaction-" + Foam::name(getNewReactionID())),
species_(species),
Tlow_(TlowDefault),
Thigh_(ThighDefault),
lhs_(lhs),
rhs_(rhs)
Thigh_(ThighDefault)
{
setThermo(thermoDatabase);
}
@ -127,13 +115,10 @@ Foam::Reaction<ReactionThermo>::Reaction
const speciesTable& species
)
:
reaction(r, species),
ReactionThermo::thermoType(r),
name_(r.name() + "Copy"),
species_(species),
Tlow_(r.Tlow()),
Thigh_(r.Thigh()),
lhs_(r.lhs_),
rhs_(r.rhs_)
Thigh_(r.Thigh())
{}
@ -145,19 +130,11 @@ Foam::Reaction<ReactionThermo>::Reaction
const dictionary& dict
)
:
reaction(species, dict),
ReactionThermo::thermoType(*thermoDatabase[species[0]]),
name_(dict.dictName()),
species_(species),
Tlow_(dict.lookupOrDefault<scalar>("Tlow", TlowDefault)),
Thigh_(dict.lookupOrDefault<scalar>("Thigh", ThighDefault))
{
specieCoeffs::setLRhs
(
IStringStream(dict.lookup("reaction"))(),
species_,
lhs_,
rhs_
);
setThermo(thermoDatabase);
}
@ -304,13 +281,7 @@ Foam::Reaction<ReactionThermo>::New
template<class ReactionThermo>
void Foam::Reaction<ReactionThermo>::write(Ostream& os) const
{
OStringStream reaction;
writeEntry
(
os,
"reaction",
specieCoeffs::reactionStr(reaction, species_, lhs_, rhs_)
);
reaction::write(os);
}
@ -358,16 +329,16 @@ void Foam::Reaction<ReactionThermo>::omega
p, T, c, li, pf, cf, lRef, pr, cr, rRef
);
forAll(lhs_, i)
forAll(lhs(), i)
{
const label si = lhs_[i].index;
const scalar sl = lhs_[i].stoichCoeff;
const label si = lhs()[i].index;
const scalar sl = lhs()[i].stoichCoeff;
dcdt[si] -= sl*omegaI;
}
forAll(rhs_, i)
forAll(rhs(), i)
{
const label si = rhs_[i].index;
const scalar sr = rhs_[i].stoichCoeff;
const label si = rhs()[i].index;
const scalar sr = rhs()[i].stoichCoeff;
dcdt[si] += sr*omegaI;
}
}
@ -397,34 +368,34 @@ Foam::scalar Foam::Reaction<ReactionThermo>::omega
pf = 1;
pr = 1;
const label Nl = lhs_.size();
const label Nr = rhs_.size();
const label Nl = lhs().size();
const label Nr = rhs().size();
label slRef = 0;
lRef = lhs_[slRef].index;
lRef = lhs()[slRef].index;
pf = kf;
for (label s = 1; s < Nl; s++)
{
const label si = lhs_[s].index;
const label si = lhs()[s].index;
if (c[si] < c[lRef])
{
const scalar exp = lhs_[slRef].exponent;
const scalar exp = lhs()[slRef].exponent;
pf *= pow(max(c[lRef], 0), exp);
lRef = si;
slRef = s;
}
else
{
const scalar exp = lhs_[s].exponent;
const scalar exp = lhs()[s].exponent;
pf *= pow(max(c[si], 0), exp);
}
}
cf = max(c[lRef], 0);
{
const scalar exp = lhs_[slRef].exponent;
const scalar exp = lhs()[slRef].exponent;
if (exp < 1)
{
if (cf > small)
@ -443,30 +414,30 @@ Foam::scalar Foam::Reaction<ReactionThermo>::omega
}
label srRef = 0;
rRef = rhs_[srRef].index;
rRef = rhs()[srRef].index;
// Find the matrix element and element position for the rhs
pr = kr;
for (label s = 1; s < Nr; s++)
{
const label si = rhs_[s].index;
const label si = rhs()[s].index;
if (c[si] < c[rRef])
{
const scalar exp = rhs_[srRef].exponent;
const scalar exp = rhs()[srRef].exponent;
pr *= pow(max(c[rRef], 0), exp);
rRef = si;
srRef = s;
}
else
{
const scalar exp = rhs_[s].exponent;
const scalar exp = rhs()[s].exponent;
pr *= pow(max(c[si], 0), exp);
}
}
cr = max(c[rRef], 0);
{
const scalar exp = rhs_[srRef].exponent;
const scalar exp = rhs()[srRef].exponent;
if (exp < 1)
{
if (cr > small)
@ -509,30 +480,30 @@ void Foam::Reaction<ReactionThermo>::dwdc
omegaI = omega(p, T, c, li, pf, cf, lRef, pr, cr, rRef);
forAll(lhs_, i)
forAll(lhs(), i)
{
const label si = reduced ? c2s[lhs_[i].index] : lhs_[i].index;
const scalar sl = lhs_[i].stoichCoeff;
const label si = reduced ? c2s[lhs()[i].index] : lhs()[i].index;
const scalar sl = lhs()[i].stoichCoeff;
dcdt[si] -= sl*omegaI;
}
forAll(rhs_, i)
forAll(rhs(), i)
{
const label si = reduced ? c2s[rhs_[i].index] : rhs_[i].index;
const scalar sr = rhs_[i].stoichCoeff;
const label si = reduced ? c2s[rhs()[i].index] : rhs()[i].index;
const scalar sr = rhs()[i].stoichCoeff;
dcdt[si] += sr*omegaI;
}
kfwd = this->kf(p, T, c, li);
kbwd = this->kr(kfwd, p, T, c, li);
forAll(lhs_, j)
forAll(lhs(), j)
{
const label sj = reduced ? c2s[lhs_[j].index] : lhs_[j].index;
const label sj = reduced ? c2s[lhs()[j].index] : lhs()[j].index;
scalar kf = kfwd;
forAll(lhs_, i)
forAll(lhs(), i)
{
const label si = lhs_[i].index;
const scalar el = lhs_[i].exponent;
const label si = lhs()[i].index;
const scalar el = lhs()[i].exponent;
if (i == j)
{
if (el < 1)
@ -557,28 +528,28 @@ void Foam::Reaction<ReactionThermo>::dwdc
}
}
forAll(lhs_, i)
forAll(lhs(), i)
{
const label si = reduced ? c2s[lhs_[i].index] : lhs_[i].index;
const scalar sl = lhs_[i].stoichCoeff;
const label si = reduced ? c2s[lhs()[i].index] : lhs()[i].index;
const scalar sl = lhs()[i].stoichCoeff;
J(si, sj) -= sl*kf;
}
forAll(rhs_, i)
forAll(rhs(), i)
{
const label si = reduced ? c2s[rhs_[i].index] : rhs_[i].index;
const scalar sr = rhs_[i].stoichCoeff;
const label si = reduced ? c2s[rhs()[i].index] : rhs()[i].index;
const scalar sr = rhs()[i].stoichCoeff;
J(si, sj) += sr*kf;
}
}
forAll(rhs_, j)
forAll(rhs(), j)
{
const label sj = reduced ? c2s[rhs_[j].index] : rhs_[j].index;
const label sj = reduced ? c2s[rhs()[j].index] : rhs()[j].index;
scalar kr = kbwd;
forAll(rhs_, i)
forAll(rhs(), i)
{
const label si = rhs_[i].index;
const scalar er = rhs_[i].exponent;
const label si = rhs()[i].index;
const scalar er = rhs()[i].exponent;
if (i == j)
{
if (er < 1)
@ -603,16 +574,16 @@ void Foam::Reaction<ReactionThermo>::dwdc
}
}
forAll(lhs_, i)
forAll(lhs(), i)
{
const label si = reduced ? c2s[lhs_[i].index] : lhs_[i].index;
const scalar sl = lhs_[i].stoichCoeff;
const label si = reduced ? c2s[lhs()[i].index] : lhs()[i].index;
const scalar sl = lhs()[i].stoichCoeff;
J(si, sj) += sl*kr;
}
forAll(rhs_, i)
forAll(rhs(), i)
{
const label si = reduced ? c2s[rhs_[i].index] : rhs_[i].index;
const scalar sr = rhs_[i].stoichCoeff;
const label si = reduced ? c2s[rhs()[i].index] : rhs()[i].index;
const scalar sr = rhs()[i].stoichCoeff;
J(si, sj) -= sr*kr;
}
}
@ -632,18 +603,18 @@ void Foam::Reaction<ReactionThermo>::dwdc
sj = reduced ? c2s[sj] : sj;
if (sj != -1)
{
forAll(lhs_, i)
forAll(lhs(), i)
{
const label si =
reduced ? c2s[lhs_[i].index] : lhs_[i].index;
const scalar sl = lhs_[i].stoichCoeff;
reduced ? c2s[lhs()[i].index] : lhs()[i].index;
const scalar sl = lhs()[i].stoichCoeff;
J(si, sj) -= sl*dcidc[j]*omegaI;
}
forAll(rhs_, i)
forAll(rhs(), i)
{
const label si =
reduced ? c2s[rhs_[i].index] : rhs_[i].index;
const scalar sr = rhs_[i].stoichCoeff;
reduced ? c2s[rhs()[i].index] : rhs()[i].index;
const scalar sr = rhs()[i].stoichCoeff;
J(si, sj) += sr*dcidc[j]*omegaI;
}
}
@ -675,10 +646,10 @@ void Foam::Reaction<ReactionThermo>::dwdT
scalar dkrdT = this->dkrdT(p, T, c, li, dkfdT, kr);
scalar sumExp = 0.0;
forAll(lhs_, i)
forAll(lhs(), i)
{
const label si = lhs_[i].index;
const scalar el = lhs_[i].exponent;
const label si = lhs()[i].index;
const scalar el = lhs()[i].exponent;
const scalar cExp = pow(c[si], el);
dkfdT *= cExp;
kf *= cExp;
@ -687,10 +658,10 @@ void Foam::Reaction<ReactionThermo>::dwdT
kf *= -sumExp/T;
sumExp = 0.0;
forAll(rhs_, i)
forAll(rhs(), i)
{
const label si = rhs_[i].index;
const scalar er = rhs_[i].exponent;
const label si = rhs()[i].index;
const scalar er = rhs()[i].exponent;
const scalar cExp = pow(c[si], er);
dkrdT *= cExp;
kr *= cExp;
@ -707,26 +678,19 @@ void Foam::Reaction<ReactionThermo>::dwdT
dcidT *= omegaI;
// J(i, indexT) = sum_reactions nu_i dqdT
forAll(lhs_, i)
forAll(lhs(), i)
{
const label si = reduced ? c2s[lhs_[i].index] : lhs_[i].index;
const scalar sl = lhs_[i].stoichCoeff;
const label si = reduced ? c2s[lhs()[i].index] : lhs()[i].index;
const scalar sl = lhs()[i].stoichCoeff;
J(si, indexT) -= sl*(dqidT + dcidT);
}
forAll(rhs_, i)
forAll(rhs(), i)
{
const label si = reduced ? c2s[rhs_[i].index] : rhs_[i].index;
const scalar sr = rhs_[i].stoichCoeff;
const label si = reduced ? c2s[rhs()[i].index] : rhs()[i].index;
const scalar sr = rhs()[i].stoichCoeff;
J(si, indexT) += sr*(dqidT + dcidT);
}
}
template<class ReactionThermo>
const Foam::speciesTable& Foam::Reaction<ReactionThermo>::species() const
{
return species_;
}
// ************************************************************************* //

View File

@ -37,8 +37,7 @@ SourceFiles
#ifndef Reaction_H
#define Reaction_H
#include "speciesTable.H"
#include "specieCoeffs.H"
#include "reaction.H"
#include "HashPtrTable.H"
#include "scalarField.H"
#include "simpleMatrix.H"
@ -61,6 +60,7 @@ inline Ostream& operator<<(Ostream&, const Reaction<ReactionThermo>&);
class objectRegistry;
/*---------------------------------------------------------------------------*\
Class Reaction Declaration
\*---------------------------------------------------------------------------*/
@ -68,6 +68,7 @@ class objectRegistry;
template<class ReactionThermo>
class Reaction
:
public reaction,
public ReactionThermo::thermoType
{
@ -75,40 +76,23 @@ public:
// Static data
//- Number of un-named reactions
static label nUnNamedReactions;
//- Default temperature limits of applicability of reaction rates
static scalar TlowDefault, ThighDefault;
private:
// Private Data
//- Name of reaction
const word name_;
//- List of specie names present in reaction system
const speciesTable& species_;
//- Temperature limits of applicability of reaction rates
scalar Tlow_, Thigh_;
//- Specie info for the left-hand-side of the reaction
List<specieCoeffs> lhs_;
//- Specie info for the right-hand-side of the reaction
List<specieCoeffs> rhs_;
// Private Member Functions
//- Construct reaction thermo
void setThermo(const HashPtrTable<ReactionThermo>& thermoDatabase);
//- Return new reaction ID for un-named reactions
label getNewReactionID();
public:
@ -219,24 +203,12 @@ public:
// Access
//- Return the name of the reaction
inline const word& name() const;
//- Return the lower temperature limit for the reaction
inline scalar Tlow() const;
//- Return the upper temperature limit for the reaction
inline scalar Thigh() const;
//- Return the components of the left hand side
inline const List<specieCoeffs>& lhs() const;
//- Return the components of the right hand side
inline const List<specieCoeffs>& rhs() const;
//- Return the specie list
const speciesTable& species() const;
// Reaction rate coefficients

View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -25,60 +25,34 @@ License
#include "Reaction.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class ReactionThermo>
inline const word& Reaction<ReactionThermo>::name() const
{
return name_;
}
template<class ReactionThermo>
inline scalar Reaction<ReactionThermo>::Tlow() const
inline Foam::scalar Foam::Reaction<ReactionThermo>::Tlow() const
{
return Tlow_;
}
template<class ReactionThermo>
inline scalar Reaction<ReactionThermo>::Thigh() const
inline Foam::scalar Foam::Reaction<ReactionThermo>::Thigh() const
{
return Thigh_;
}
template<class ReactionThermo>
inline const List<specieCoeffs>& Reaction<ReactionThermo>::lhs() const
{
return lhs_;
}
template<class ReactionThermo>
inline const List<specieCoeffs>& Reaction<ReactionThermo>::rhs() const
{
return rhs_;
}
// * * * * * * * * * * * * * * * Ostream Operator * * * * * * * * * * * * * //
template<class ReactionThermo>
inline Ostream& operator<<(Ostream& os, const Reaction<ReactionThermo>& r)
inline Foam::Ostream& Foam::operator<<
(
Ostream& os,
const Reaction<ReactionThermo>& r
)
{
r.write(os);
return os;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// ************************************************************************* //

View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2020 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -23,87 +23,82 @@ License
\*---------------------------------------------------------------------------*/
#include "reaction.H"
#include "dictionary.H"
// * * * * * * * * * * * * * * * * Static Data * * * * * * * * * * * * * * * //
Foam::label Foam::reaction::nUnNamedReactions(0);
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
Foam::label Foam::reaction::getNewReactionID()
{
return nUnNamedReactions++;
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
inline Foam::infiniteReactionRate::infiniteReactionRate()
{}
inline Foam::infiniteReactionRate::infiniteReactionRate
Foam::reaction::reaction
(
const speciesTable&,
const dictionary&
const speciesTable& species,
const List<specieCoeffs>& lhs,
const List<specieCoeffs>& rhs
)
:
name_("un-named-reaction-" + Foam::name(getNewReactionID())),
species_(species),
lhs_(lhs),
rhs_(rhs)
{}
inline void Foam::infiniteReactionRate::write(Ostream& os) const
Foam::reaction::reaction
(
const reaction& r,
const speciesTable& species
)
:
name_(r.name() + "Copy"),
species_(species),
lhs_(r.lhs_),
rhs_(r.rhs_)
{}
Foam::reaction::reaction
(
const speciesTable& species,
const dictionary& dict
)
:
name_(dict.dictName()),
species_(species)
{
specieCoeffs::setLRhs
(
IStringStream(dict.lookup("reaction"))(),
species_,
lhs_,
rhs_
);
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
inline Foam::scalar Foam::infiniteReactionRate::operator()
(
const scalar p,
const scalar,
const scalarField&,
const label
) const
void Foam::reaction::write(Ostream& os) const
{
return (1);
}
inline Foam::scalar Foam::infiniteReactionRate::ddT
(
const scalar p,
const scalar,
const scalarField&,
const label
) const
{
return (0);
}
inline const Foam::List<Foam::Tuple2<Foam::label, Foam::scalar>>&
Foam::infiniteReactionRate::beta() const
{
return NullObjectRef<List<Tuple2<label, scalar>>>();
}
inline void Foam::infiniteReactionRate::dcidc
(
const scalar p,
const scalar T,
const scalarField& c,
const label li,
scalarField& dcidc
) const
{}
inline Foam::scalar Foam::infiniteReactionRate::dcidT
(
const scalar p,
const scalar T,
const scalarField& c,
const label li
) const
{
return 0;
}
inline Foam::Ostream& Foam::operator<<
(
Ostream& os,
const infiniteReactionRate& rr
)
{
rr.write(os);
return os;
OStringStream reaction;
writeEntry
(
os,
"reaction",
specieCoeffs::reactionStr(reaction, species_, lhs_, rhs_)
);
}

View File

@ -0,0 +1,164 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2020 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::reaction
Description
Reaction base-class holding the specie names and coefficients
SourceFiles
reactionI.H
reaction.C
\*---------------------------------------------------------------------------*/
#ifndef reaction_H
#define reaction_H
#include "speciesTable.H"
#include "specieCoeffs.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// Forward declaration of friend functions and operators
class reaction;
inline Ostream& operator<<(Ostream&, const reaction&);
/*---------------------------------------------------------------------------*\
Class reaction Declaration
\*---------------------------------------------------------------------------*/
class reaction
{
public:
// Static data
//- Number of un-named reactions
static label nUnNamedReactions;
private:
// Private Data
//- Name of reaction
const word name_;
//- List of specie names present in reaction system
const speciesTable& species_;
//- Specie info for the left-hand-side of the reaction
List<specieCoeffs> lhs_;
//- Specie info for the right-hand-side of the reaction
List<specieCoeffs> rhs_;
// Private Member Functions
//- Return new reaction ID for un-named reactions
label getNewReactionID();
public:
// Constructors
//- Construct from components
reaction
(
const speciesTable& species,
const List<specieCoeffs>& lhs,
const List<specieCoeffs>& rhs
);
//- Construct as copy given new speciesTable
reaction(const reaction&, const speciesTable& species);
//- Construct from dictionary
reaction
(
const speciesTable& species,
const dictionary& dict
);
//- Destructor
~reaction()
{}
// Member Functions
// Access
//- Return the name of the reaction
inline const word& name() const;
//- Return the components of the left hand side
inline const List<specieCoeffs>& lhs() const;
//- Return the components of the right hand side
inline const List<specieCoeffs>& rhs() const;
//- Return the specie list
inline const speciesTable& species() const;
//- Write
void write(Ostream&) const;
// Member Operators
//- Disallow default bitwise assignment
void operator=(const reaction&) = delete;
// Ostream Operator
friend Ostream& operator<<(Ostream&, const reaction&);
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "reactionI.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,63 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2020 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "reaction.H"
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
inline const Foam::word& Foam::reaction::name() const
{
return name_;
}
inline const Foam::List<Foam::specieCoeffs>& Foam::reaction::lhs() const
{
return lhs_;
}
inline const Foam::List<Foam::specieCoeffs>& Foam::reaction::rhs() const
{
return rhs_;
}
inline const Foam::speciesTable& Foam::reaction::species() const
{
return species_;
}
// * * * * * * * * * * * * * * * Ostream Operator * * * * * * * * * * * * * //
inline Foam::Ostream& Foam::operator<<(Foam::Ostream& os, const reaction& r)
{
r.write(os);
return os;
}
// ************************************************************************* //

View File

@ -1,149 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::infiniteReactionRate
Description
infinite reaction rate.
SourceFiles
infiniteReactionRateI.H
\*---------------------------------------------------------------------------*/
#ifndef infiniteReactionRate_H
#define infiniteReactionRate_H
#include "scalarField.H"
#include "typeInfo.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// Forward declaration of friend functions and operators
class infiniteReactionRate;
Ostream& operator<<(Ostream&, const infiniteReactionRate&);
/*---------------------------------------------------------------------------*\
Class infiniteReactionRate Declaration
\*---------------------------------------------------------------------------*/
class infiniteReactionRate
{
public:
// Constructors
//- Null constructor
inline infiniteReactionRate
();
//- Construct from dictionary
inline infiniteReactionRate
(
const speciesTable& species,
const dictionary& dict
);
// Member Functions
//- Return the type name
static word type()
{
return "infinite";
}
inline scalar operator()
(
const scalar p,
const scalar T,
const scalarField& c,
const label li
) const;
inline scalar ddT
(
const scalar p,
const scalar T,
const scalarField& c,
const label li
) const;
//- Third-body efficiencies (beta = 1-alpha)
// non-empty only for third-body reactions
// with enhanced molecularity (alpha != 1)
inline const List<Tuple2<label, scalar>>& beta() const;
//- Species concentration derivative of the pressure dependent term
inline void dcidc
(
const scalar p,
const scalar T,
const scalarField& c,
const label li,
scalarField& dcidc
) const;
//- Temperature derivative of the pressure dependent term
inline scalar dcidT
(
const scalar p,
const scalar T,
const scalarField& c,
const label li
) const;
//- Write to stream
inline void write(Ostream& os) const;
// Ostream Operator
inline friend Ostream& operator<<
(
Ostream&,
const infiniteReactionRate&
);
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "infiniteReactionRateI.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -8,7 +8,6 @@
reaction
{
type irreversibleInfinite;
reaction "CH4 + 2O2 + 7.5N2 = CO2 + 2H2O + 7.5N2";
}

View File

@ -8,7 +8,6 @@
reaction
{
type irreversibleInfinite;
reaction "CH4 + 2O2 + 7.5N2 = CO2 + 2H2O + 7.5N2";
}