30 template<
class CompType,
class ThermoType>
38 searchInitSet_(this->coeffsDict_.subDict(
"initialSet").size())
46 if (initSet.
found(chemistry.
Y()[i].member()))
48 searchInitSet_[j++] = i;
52 if (j<searchInitSet_.size())
55 << searchInitSet_.size()-j
56 <<
" species in the initial set is not in the mechanism " 65 template<
class CompType,
class ThermoType>
72 template<
class CompType,
class ThermoType>
80 scalarField& completeC(this->chemistry_.completeC());
83 for (
label i=0; i<this->nSpecie_; i++)
89 c1[this->nSpecie_] =
T;
90 c1[this->nSpecie_+1] =
p;
105 scalar pf, cf, pr, cr;
107 forAll(this->chemistry_.reactions(), i)
111 scalar omegai = this->chemistry_.omega
113 R, c1, T, p, pf, cf, lRef, pr, cr, rRef
126 scalar sl = -R.
lhs()[
s].stoichCoeff;
139 wA.append(sl*omegai);
146 scalar sl = R.
rhs()[
s].stoichCoeff;
159 wA.append(sl*omegai);
167 label curID = wAID[id];
168 scalar curwA = wA[id];
184 deltaBi[curID] =
false;
186 while(!usedIndex.empty())
190 if (deltaBi[curIndex])
192 deltaBi[curIndex] =
false;
193 if (rABPos(curID, curIndex)==-1)
195 rABPos(curID, curIndex) = NbrABInit[curID];
196 rABOtherSpec(curID, NbrABInit[curID]) = curIndex;
199 PAB(curID, NbrABInit[curID]) = curwA;
203 CAB(curID, NbrABInit[curID]) = -curwA;
211 PAB(curID, rABPos(curID, curIndex)) += curwA;
215 CAB(curID, rABPos(curID, curIndex)) += -curwA;
225 if (PA[curID] == 0.0)
236 if (CA[curID] == 0.0)
268 this->nSpecie_, this->nSpecie_, -1
273 for (
int i=0; i<NbrABInit[A]; i++)
275 label ri = rABOtherSpec(A, i);
276 scalar maxPACA =
max(PA[ri],CA[ri]);
277 if (maxPACA > vSmall)
279 for (
int j=0; j<NbrABInit[ri]; j++)
281 label B = rABOtherSpec(ri, j);
284 if (rABPos2nd(A, B)==-1)
286 rABPos2nd(A, B) = NbrABInit2nd[A]++;
287 rABOtherSpec2nd(A, rABPos2nd(A, B)) = B;
288 PAB2nd(A, rABPos2nd(A, B)) =
289 PAB(A, i)*PAB(ri, j)/maxPACA;
290 CAB2nd(A, rABPos2nd(A, B)) =
291 CAB(A, i)*CAB(ri, j)/maxPACA;
295 PAB2nd(A, rABPos2nd(A, B)) +=
296 PAB(A, i)*PAB(ri, j)/maxPACA;
297 CAB2nd(A, rABPos2nd(A, B)) +=
298 CAB(A, i)*CAB(ri, j)/maxPACA;
307 label speciesNumber = 0;
311 for (
label i=0; i<this->nSpecie_; i++)
313 this->activeSpecies_[i] =
false;
317 const labelList& SIS(this->searchInitSet_);
323 this->activeSpecies_[q] =
true;
332 scalar Den =
max(PA[u],CA[u]);
337 for (
label v=0; v<NbrABInit[u]; v++)
339 label otherSpec = rABOtherSpec(u, v);
340 scalar rAB = (PAB(u, v)+CAB(u, v))/Den;
341 label id2nd = rABPos2nd(u, otherSpec);
344 rAB += (PAB2nd(u, id2nd)+CAB2nd(u, id2nd))/Den;
349 rAB >= this->tolerance()
350 && !this->activeSpecies_[otherSpec]
354 this->activeSpecies_[otherSpec] =
true;
360 for (
label v=0; v<NbrABInit2nd[u]; v++)
362 label otherSpec = rABOtherSpec2nd(u, v);
363 scalar rAB = (PAB2nd(u, v)+CAB2nd(u, v))/Den;
367 rAB >= this->tolerance()
368 && !this->activeSpecies_[otherSpec]
372 this->activeSpecies_[otherSpec] =
true;
380 forAll(this->chemistry_.reactions(), i)
383 this->chemistry_.reactionsDisabled()[i] =
false;
388 if (!this->activeSpecies_[ss])
390 this->chemistry_.reactionsDisabled()[i] =
true;
394 if (!this->chemistry_.reactionsDisabled()[i])
399 if (!this->activeSpecies_[ss])
401 this->chemistry_.reactionsDisabled()[i] =
true;
408 this->NsSimp_ = speciesNumber;
409 scalarField& simplifiedC(this->chemistry_.simplifiedC());
410 simplifiedC.
setSize(this->NsSimp_+2);
413 Field<label>& c2s(this->chemistry_.completeToSimplifiedIndex());
416 for (
label i=0; i<this->nSpecie_; i++)
418 if (this->activeSpecies_[i])
421 simplifiedC[j] = c[i];
423 if (!this->chemistry_.active(i))
425 this->chemistry_.setActive(i);
433 simplifiedC[this->NsSimp_] =
T;
434 simplifiedC[this->NsSimp_+1] =
p;
435 this->chemistry_.setNsDAC(this->NsSimp_);
438 this->chemistry_.setNSpecie(this->NsSimp_);
A FIFO stack based on a singly-linked list.
virtual label nSpecie() const
The number of species.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
PtrList< volScalarField > & Y()
A list of keyword definitions, which are a keyword followed by any number of values (e...
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
void size(const label)
Override size to be inconsistent with allocated storage.
const List< specieCoeffs > & lhs() const
Return the components of the left hand side.
BasicChemistryModel< rhoReactionThermo > & chemistry
Simple extension of ReactionThermo to handle reaction kinetics in addition to the equilibrium thermod...
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
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))
void setSize(const label)
Alter the addressed list size.
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
errorManip< error > abort(error &err)
virtual ~PFA()
Destructor.
void setSize(const label)
Reset size of List.
#define R(A, B, C, D, E, F, K, M)
void push(const T &a)
Push an element onto the stack.
const dimensionedScalar c1
First radiation constant: default SI units: [W/m^2].
virtual void reduceMechanism(const scalarField &c, const scalar T, const scalar p)
Reduce the mechanism.
T pop()
Pop the bottom element off the stack.
const List< specieCoeffs > & rhs() const
Return the components of the right hand side.
An abstract class for methods of chemical mechanism reduction.
PFA(const IOdictionary &dict, TDACChemistryModel< CompType, ThermoType > &chemistry)
Construct from components.