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>
81 scalarField& completeC(this->chemistry_.completeC());
84 for (
label i=0; i<this->nSpecie_; i++)
90 c1[this->nSpecie_] =
T;
91 c1[this->nSpecie_+1] =
p;
106 scalar pf, cf, pr, cr;
108 forAll(this->chemistry_.reactions(), i)
112 scalar omegai = this->chemistry_.omega
114 R, p, T, c1, li, pf, cf, lRef, pr, cr, rRef
127 scalar sl = -R.
lhs()[
s].stoichCoeff;
140 wA.append(sl*omegai);
147 scalar sl = R.
rhs()[
s].stoichCoeff;
160 wA.append(sl*omegai);
168 label curID = wAID[id];
169 scalar curwA = wA[id];
185 deltaBi[curID] =
false;
187 while(!usedIndex.empty())
191 if (deltaBi[curIndex])
193 deltaBi[curIndex] =
false;
194 if (rABPos(curID, curIndex)==-1)
196 rABPos(curID, curIndex) = NbrABInit[curID];
197 rABOtherSpec(curID, NbrABInit[curID]) = curIndex;
200 PAB(curID, NbrABInit[curID]) = curwA;
204 CAB(curID, NbrABInit[curID]) = -curwA;
212 PAB(curID, rABPos(curID, curIndex)) += curwA;
216 CAB(curID, rABPos(curID, curIndex)) += -curwA;
226 if (PA[curID] == 0.0)
237 if (CA[curID] == 0.0)
269 this->nSpecie_, this->nSpecie_, -1
274 for (
int i=0; i<NbrABInit[A]; i++)
276 label ri = rABOtherSpec(A, i);
277 scalar maxPACA =
max(PA[ri],CA[ri]);
278 if (maxPACA > vSmall)
280 for (
int j=0; j<NbrABInit[ri]; j++)
282 label B = rABOtherSpec(ri, j);
285 if (rABPos2nd(A, B)==-1)
287 rABPos2nd(A, B) = NbrABInit2nd[A]++;
288 rABOtherSpec2nd(A, rABPos2nd(A, B)) = B;
289 PAB2nd(A, rABPos2nd(A, B)) =
290 PAB(A, i)*PAB(ri, j)/maxPACA;
291 CAB2nd(A, rABPos2nd(A, B)) =
292 CAB(A, i)*CAB(ri, j)/maxPACA;
296 PAB2nd(A, rABPos2nd(A, B)) +=
297 PAB(A, i)*PAB(ri, j)/maxPACA;
298 CAB2nd(A, rABPos2nd(A, B)) +=
299 CAB(A, i)*CAB(ri, j)/maxPACA;
308 label speciesNumber = 0;
312 for (
label i=0; i<this->nSpecie_; i++)
314 this->activeSpecies_[i] =
false;
318 const labelList& SIS(this->searchInitSet_);
324 this->activeSpecies_[q] =
true;
333 scalar Den =
max(PA[u],CA[u]);
338 for (
label v=0; v<NbrABInit[u]; v++)
340 label otherSpec = rABOtherSpec(u, v);
341 scalar rAB = (PAB(u, v)+CAB(u, v))/Den;
342 label id2nd = rABPos2nd(u, otherSpec);
345 rAB += (PAB2nd(u, id2nd)+CAB2nd(u, id2nd))/Den;
350 rAB >= this->tolerance()
351 && !this->activeSpecies_[otherSpec]
355 this->activeSpecies_[otherSpec] =
true;
361 for (
label v=0; v<NbrABInit2nd[u]; v++)
363 label otherSpec = rABOtherSpec2nd(u, v);
364 scalar rAB = (PAB2nd(u, v)+CAB2nd(u, v))/Den;
368 rAB >= this->tolerance()
369 && !this->activeSpecies_[otherSpec]
373 this->activeSpecies_[otherSpec] =
true;
381 forAll(this->chemistry_.reactions(), i)
384 this->chemistry_.reactionsDisabled()[i] =
false;
389 if (!this->activeSpecies_[ss])
391 this->chemistry_.reactionsDisabled()[i] =
true;
395 if (!this->chemistry_.reactionsDisabled()[i])
400 if (!this->activeSpecies_[ss])
402 this->chemistry_.reactionsDisabled()[i] =
true;
409 this->NsSimp_ = speciesNumber;
410 scalarField& simplifiedC(this->chemistry_.simplifiedC());
411 simplifiedC.
setSize(this->NsSimp_+2);
414 Field<label>& c2s(this->chemistry_.completeToSimplifiedIndex());
417 for (
label i=0; i<this->nSpecie_; i++)
419 if (this->activeSpecies_[i])
422 simplifiedC[j] = c[i];
424 if (!this->chemistry_.active(i))
426 this->chemistry_.setActive(i);
434 simplifiedC[this->NsSimp_] =
T;
435 simplifiedC[this->NsSimp_+1] =
p;
436 this->chemistry_.setNsDAC(this->NsSimp_);
439 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.
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.
const dimensionedScalar & c1
First radiation constant: default SI units: [W/m^2].
void size(const label)
Override size to be inconsistent with allocated storage.
const PtrList< volScalarField > & Y()
const List< specieCoeffs > & lhs() const
Return the components of the left hand side.
BasicChemistryModel< rhoReactionThermo > & chemistry
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.
errorManip< error > abort(error &err)
virtual void reduceMechanism(const scalar p, const scalar T, const scalarField &c, const label li)
Reduce the mechanism.
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.
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.