All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
DAC.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2016-2020 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 "DAC.H"
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template<class CompType, class ThermoType>
32 (
33  const IOdictionary& dict,
35 )
36 :
38  searchInitSet_(this->coeffsDict_.subDict("initialSet").size()),
39  zprime_(0),
40  nbCLarge_(3),
41  sC_(this->nSpecie_,0),
42  sH_(this->nSpecie_,0),
43  sO_(this->nSpecie_,0),
44  sN_(this->nSpecie_,0),
45  CO2Id_(-1),
46  COId_(-1),
47  HO2Id_(-1),
48  H2OId_(-1),
49  NOId_(-1),
50  automaticSIS_(true),
51  phiTol_(this->tolerance()),
52  NOxThreshold_(1800),
53  CO2Name_
54  (
55  dict.subDict("reduction").lookupOrDefault<word>
56  (
57  "CO2Name","CO2"
58  )
59  ),
60  COName_
61  (
62  dict.subDict("reduction").lookupOrDefault<word>
63  (
64  "COName","CO"
65  )
66  ),
67  HO2Name_
68  (
69  dict.subDict("reduction").lookupOrDefault<word>
70  (
71  "HO2Name","HO2"
72  )
73  ),
74  H2OName_
75  (
76  dict.subDict("reduction").lookupOrDefault<word>
77  (
78  "H2OName","H2O"
79  )
80  ),
81  NOName_
82  (
83  dict.subDict("reduction").lookupOrDefault<word>
84  (
85  "NOName","NO"
86  )
87  ),
88  forceFuelInclusion_(false)
89 {
90  label j=0;
91  dictionary initSet = this->coeffsDict_.subDict("initialSet");
92  for (label i=0; i<chemistry.nSpecie(); i++)
93  {
94  if (initSet.found(chemistry.Y()[i].member()))
95  {
96  searchInitSet_[j++] = i;
97  }
98  }
99  if (j<searchInitSet_.size())
100  {
102  << searchInitSet_.size()-j
103  << " species in the initial set is not in the mechanism "
104  << initSet
105  << exit(FatalError);
106  }
107 
108  if (this->coeffsDict_.found("automaticSIS"))
109  {
110  automaticSIS_.readIfPresent("automaticSIS", this->coeffsDict_);
111  }
112 
113  if (this->coeffsDict_.found("forceFuelInclusion"))
114  {
115  forceFuelInclusion_.readIfPresent
116  (
117  "forceFuelInclusion",this->coeffsDict_
118  );
119  }
120 
121  if (this->coeffsDict_.found("phiTol"))
122  {
123  phiTol_ = this->coeffsDict_.template lookup<scalar>("phiTol");
124  }
125 
126  if (this->coeffsDict_.found("NOxThreshold"))
127  {
128  NOxThreshold_ =
129  this->coeffsDict_.template lookup<scalar>("NOxThreshold");
130  }
131  const List<List<specieElement>>& specieComposition =
132  chemistry.specieComp();
133 
134  for (label i=0; i<this->nSpecie_; i++)
135  {
136  const List<specieElement>& curSpecieComposition =
137  specieComposition[i];
138  // For all elements in the current species
139  forAll(curSpecieComposition, j)
140  {
141  const specieElement& curElement =
142  curSpecieComposition[j];
143  if (curElement.name() == "C")
144  {
145  sC_[i] = curElement.nAtoms();
146  }
147  else if (curElement.name() == "H")
148  {
149  sH_[i] = curElement.nAtoms();
150  }
151  else if (curElement.name() == "O")
152  {
153  sO_[i] = curElement.nAtoms();
154  }
155  else if (curElement.name() == "N")
156  {
157  sN_[i] = curElement.nAtoms();
158  }
159  else
160  {
161  Info<< "element not considered"<<endl;
162  }
163  }
164  if (this->chemistry_.Y()[i].member() == CO2Name_)
165  {
166  CO2Id_ = i;
167  }
168  else if (this->chemistry_.Y()[i].member() == COName_)
169  {
170  COId_ = i;
171  }
172  else if (this->chemistry_.Y()[i].member() == HO2Name_)
173  {
174  HO2Id_ = i;
175  }
176  else if (this->chemistry_.Y()[i].member() == H2OName_)
177  {
178  H2OId_ = i;
179  }
180  else if (this->chemistry_.Y()[i].member() == NOName_)
181  {
182  NOId_ = i;
183  }
184  }
185 
186  if ((CO2Id_==-1 || COId_==-1 || HO2Id_==-1 || H2OId_==-1) && automaticSIS_)
187  {
189  << "The name of the species used in automatic SIS are not found in "
190  << " the mechanism. You should either set the name for CO2, CO, H2O"
191  << " and HO2 properly or set automaticSIS to off "
192  << exit(FatalError);
193  }
194 
195  // To compute zprime, the fuel species should be specified.
196  // According to the given mass fraction, an equivalent O/C ratio is computed
197  if (automaticSIS_)
198  {
199  dictionary fuelDict;
200  if (this->coeffsDict_.found("fuelSpecies"))
201  {
202  fuelDict = this->coeffsDict_.subDict("fuelSpecies");
203  fuelSpecies_ = fuelDict.toc();
204  if (fuelSpecies_.size() == 0)
205  {
207  << "With automatic SIS, the fuel species should be "
208  << "specified in the fuelSpecies subDict"
209  << exit(FatalError);
210  }
211  }
212  else
213  {
215  << "With automatic SIS, the fuel species should be "
216  << "specified in the fuelSpecies subDict"
217  << exit(FatalError);
218  }
219 
220  if (this->coeffsDict_.found("nbCLarge"))
221  {
222  nbCLarge_ = fuelDict.lookup<label>("nbCLarge");
223  }
224 
225  fuelSpeciesID_.setSize(fuelSpecies_.size());
226  fuelSpeciesProp_.setSize(fuelSpecies_.size());
227  scalar Mmtot(0.0);
228 
229  forAll(fuelSpecies_, i)
230  {
231  fuelSpeciesProp_[i] = fuelDict.lookup<scalar>(fuelSpecies_[i]);
232  for (label j=0; j<this->nSpecie_; j++)
233  {
234  if (this->chemistry_.Y()[j].member() == fuelSpecies_[i])
235  {
236  fuelSpeciesID_[i] = j;
237  break;
238  }
239  }
240  scalar curMm =
241  this->chemistry_.specieThermos()[fuelSpeciesID_[i]].W();
242  Mmtot += fuelSpeciesProp_[i]/curMm;
243  }
244 
245  Mmtot = 1.0/Mmtot;
246  scalar nbC(0.0);
247  scalar nbO(0.0);
248  forAll(fuelSpecies_, i)
249  {
250  label curID = fuelSpeciesID_[i];
251  scalar curMm = this->chemistry_.specieThermos()[curID].W();
252 
253  nbC += fuelSpeciesProp_[i]*Mmtot/curMm*sC_[curID];
254  nbO += fuelSpeciesProp_[i]*Mmtot/curMm*sO_[curID];
255  }
256  zprime_ = nbO/nbC;
257  }
258 }
259 
260 
261 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
262 
263 template<class CompType, class ThermoType>
265 {}
266 
267 
268 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
269 
270 
271 template<class CompType, class ThermoType>
273 (
274  const scalar p,
275  const scalar T,
276  const scalarField& c,
277  const label li
278 )
279 {
280  scalarField& completeC(this->chemistry_.completeC());
281  scalarField c1(this->chemistry_.nEqns(), 0.0);
282  for(label i=0; i<this->nSpecie_; i++)
283  {
284  c1[i] = c[i];
285  completeC[i] = c[i];
286  }
287 
288  c1[this->nSpecie_] = T;
289  c1[this->nSpecie_+1] = p;
290 
291  // Compute the rAB matrix
292  RectangularMatrix<scalar> rABNum(this->nSpecie_,this->nSpecie_,0.0);
293  scalarField PA(this->nSpecie_,0.0);
294  scalarField CA(this->nSpecie_,0.0);
295 
296  // Number of initialized rAB for each lines
297  Field<label> NbrABInit(this->nSpecie_,0);
298  // Position of the initialized rAB, -1 when not initialized
299  RectangularMatrix<label> rABPos(this->nSpecie_, this->nSpecie_, -1);
300  // Index of the other species involved in the rABNum
301  RectangularMatrix<label> rABOtherSpec(this->nSpecie_, this->nSpecie_, -1);
302 
303  scalar pf, cf, pr, cr;
304  label lRef, rRef;
305  forAll(this->chemistry_.reactions(), i)
306  {
307  const Reaction<ThermoType>& R = this->chemistry_.reactions()[i];
308 
309  // for each reaction compute omegai
310  scalar omegai = this->chemistry_.omega
311  (
312  R, p, T, c1, li, pf, cf, lRef, pr, cr, rRef
313  );
314 
315  // Then for each pair of species composing this reaction,
316  // compute the rAB matrix (separate the numerator and
317  // denominator)
318 
319  // While computing the rAB for all the species involved in the reaction
320  // we should consider that one can write a reaction A+B=2C as A+B=C+C
321  // In this case, the following algorithm only take once the effect
322  // of the species. It stores the species encountered in the reaction but
323  // use another list to see if this species has already been used
324 
325  DynamicList<scalar> wA(R.lhs().size()+R.rhs().size());
326  DynamicList<label> wAID(R.lhs().size()+R.rhs().size());
327 
328  forAll(R.lhs(), s) // Compute rAB for all species in the left hand side
329  {
330  label ss = R.lhs()[s].index;
331  scalar sl = -R.lhs()[s].stoichCoeff; // vAi = v''-v' => here -v'
332  List<bool> deltaBi(this->nSpecie_, false);
333  FIFOStack<label> usedIndex;
334  forAll(R.lhs(), j)
335  {
336  label sj = R.lhs()[j].index;
337  usedIndex.push(sj);
338  deltaBi[sj] = true;
339  }
340  forAll(R.rhs(), j)
341  {
342  label sj = R.rhs()[j].index;
343  usedIndex.push(sj);
344  deltaBi[sj] = true;
345  }
346 
347  // Disable for self reference (by definition rAA=0)
348  deltaBi[ss] = false;
349 
350  while(!usedIndex.empty())
351  {
352  label curIndex = usedIndex.pop();
353  if (deltaBi[curIndex])
354  {
355  // Disable to avoid counting it more than once
356  deltaBi[curIndex] = false;
357  // Test if this rAB is not initialized
358  if (rABPos(ss, curIndex)==-1)
359  {
360  // It starts at rABPos(ss, sj)=0
361  rABPos(ss, curIndex) = NbrABInit[ss];
362  NbrABInit[ss]++;
363  // to avoid overflow
364  rABNum(ss, rABPos(ss, curIndex)) = sl*omegai;
365  // store the other specie involved
366  rABOtherSpec(ss, rABPos(ss, curIndex)) = curIndex;
367  }
368  else
369  {
370  rABNum(ss, rABPos(ss, curIndex)) += sl*omegai;
371  }
372  }
373  }
374 
375  bool found(false);
376  forAll(wAID, id)
377  {
378  if (ss==wAID[id])
379  {
380  wA[id] += sl*omegai;
381  found = true;
382  }
383  }
384  if (!found)
385  {
386  wA.append(sl*omegai);
387  wAID.append(ss);
388  }
389  }
390 
391  forAll(R.rhs(), s) // Compute rAB for all species in the right hand side
392  {
393  label ss = R.rhs()[s].index;
394  scalar sl = R.rhs()[s].stoichCoeff; // vAi = v''-v' => here v''
395  List<bool> deltaBi(this->nSpecie_, false);
396  FIFOStack<label> usedIndex;
397  forAll(R.lhs(), j)
398  {
399  label sj = R.lhs()[j].index;
400  usedIndex.push(sj);
401  deltaBi[sj] = true;
402  }
403  forAll(R.rhs(), j)
404  {
405  label sj = R.rhs()[j].index;
406  usedIndex.push(sj);
407  deltaBi[sj] = true;
408  }
409 
410  // Disable for self reference (by definition rAA=0)
411  deltaBi[ss] = false;
412 
413  while(!usedIndex.empty())
414  {
415  label curIndex = usedIndex.pop();
416  if (deltaBi[curIndex])
417  {
418  // Disable to avoid counting it more than once
419  deltaBi[curIndex] = false;
420 
421  // Test if this rAB is not initialized
422  if (rABPos(ss, curIndex) == -1)
423  {
424  // it starts at rABPos(ss, sj)=0
425  rABPos(ss, curIndex) = NbrABInit[ss];
426  NbrABInit[ss]++;
427  rABNum(ss, rABPos(ss, curIndex)) = sl*omegai;
428  rABOtherSpec(ss, rABPos(ss, curIndex)) = curIndex;
429  }
430  else
431  {
432  rABNum(ss, rABPos(ss, curIndex)) += sl*omegai;
433  }
434  }
435  }
436 
437  bool found(false);
438  forAll(wAID, id)
439  {
440  if (ss==wAID[id])
441  {
442  wA[id] += sl*omegai;
443  found = true;
444  }
445  }
446  if (!found)
447  {
448  wA.append(sl*omegai);
449  wAID.append(ss);
450  }
451  }
452  wAID.shrink();
453 
454  // Now that every species of the reactions has been visited, we can
455  // compute the production and consumption rate. This way, it avoids
456  // getting wrong results when species are present in both lhs and rhs
457  forAll(wAID, id)
458  {
459  if (wA[id] > 0.0)
460  {
461  if (PA[wAID[id]] == 0.0)
462  {
463  PA[wAID[id]] = wA[id];
464  }
465  else
466  {
467  PA[wAID[id]] += wA[id];
468  }
469  }
470  else
471  {
472  if (CA[wAID[id]] == 0.0)
473  {
474  CA[wAID[id]] = -wA[id];
475  }
476  else
477  {
478  CA[wAID[id]] += -wA[id];
479  }
480  }
481  }
482  }
483  // rii = 0.0 by definition
484 
485  scalar phiLarge(0.0);
486  scalar phiProgress(0.0);
487  if (automaticSIS_)
488  {
489  // Compute the progress equivalence ratio
490  // and the equivalence ratio for fuel decomposition
491  label nElements = 4; // 4 main elements (C, H, O, N)
492 
493  // Total number of C, H and O (in this order)
494  scalarList Na(nElements,0.0);
495  scalarList Nal(nElements,0.0); // for large hydrocarbons
496 
497  for (label i=0; i<this->nSpecie_; i++)
498  {
499  // Complete combustion products are not considered
500  if
501  (
502  this->chemistry_.Y()[i].member() == "CO2"
503  || this->chemistry_.Y()[i].member() == "H2O"
504  )
505  {
506  continue;
507  }
508  Na[0] += sC_[i]*c[i];
509  Na[1] += sH_[i]*c[i];
510  Na[2] += sO_[i]*c[i];
511  if (sC_[i]>nbCLarge_ || this->chemistry_.Y()[i].member() == "O2")
512  {
513  Nal[0] += sC_[i]*c[i];
514  Nal[1] += sH_[i]*c[i];
515  Nal[2] += sO_[i]*c[i];
516  }
517  }
518 
519  // 2C(-CO2) + H(-H2O)/2 - z'C(-CO2)
520  // Progress equivalence ratio = ----------------------------------
521  // O(-CO2-H2O) - z' C(-CO2)
522  // where minus means that this species is not considered for the number
523  // of atoms and z' is the ratio of the number of O and C in the fuel(s)
524  phiProgress = (2*Na[0]+Na[1]/2-zprime_*Na[0])/(Na[2]-zprime_*Na[0]);
525 
526  // 2Cl + Hl/2
527  // Equivalence ratio for fuel decomposition = ----------
528  // Ol(+O2)
529  phiLarge = (2*Nal[0]+Nal[1]/2)/Nal[2];
530  }
531 
532  // Using the rAB matrix (numerator and denominator separated)
533  // compute the R value according to the search initiating set
534  scalarField Rvalue(this->nSpecie_,0.0);
535  label speciesNumber = 0;
536 
537  // Set all species to inactive and activate them according
538  // to rAB and initial set
539  for (label i=0; i<this->nSpecie_; i++)
540  {
541  this->activeSpecies_[i] = false;
542  }
543 
544  // Initialize the FIFOStack for search set
546 
547  const labelList& SIS(searchInitSet_);
548 
549  // If automaticSIS is on, the search initiating set is selected according to
550  // phiProgress and phiLarge
551  if (automaticSIS_)
552  {
553  if (phiLarge >= phiTol_ && phiProgress >= phiTol_)
554  {
555  // When phiLarge and phiProgress >= phiTol then
556  // CO, HO2 and fuel are in the SIS
557  Q.push(COId_);
558  speciesNumber++;
559  this->activeSpecies_[COId_] = true;
560  Rvalue[COId_] = 1.0;
561  Q.push(HO2Id_);
562  speciesNumber++;
563  this->activeSpecies_[HO2Id_] = true;
564  Rvalue[HO2Id_] = 1.0;
565  forAll(fuelSpeciesID_,i)
566  {
567  Q.push(fuelSpeciesID_[i]);
568  speciesNumber++;
569  this->activeSpecies_[fuelSpeciesID_[i]] = true;
570  Rvalue[fuelSpeciesID_[i]] = 1.0;
571  }
572 
573  }
574  else if (phiLarge < phiTol_ && phiProgress >= phiTol_)
575  {
576  // When phiLarge < phiTol and phiProgress >= phiTol then
577  // CO, HO2 are in the SIS
578  Q.push(COId_);
579  speciesNumber++;
580  this->activeSpecies_[COId_] = true;
581  Rvalue[COId_] = 1.0;
582  Q.push(HO2Id_);
583  speciesNumber++;
584  this->activeSpecies_[HO2Id_] = true;
585  Rvalue[HO2Id_] = 1.0;
586 
587  if (forceFuelInclusion_)
588  {
589  forAll(fuelSpeciesID_,i)
590  {
591  Q.push(fuelSpeciesID_[i]);
592  speciesNumber++;
593  this->activeSpecies_[fuelSpeciesID_[i]] = true;
594  Rvalue[fuelSpeciesID_[i]] = 1.0;
595  }
596  }
597  }
598  else
599  {
600  // When phiLarge and phiProgress< phiTol then
601  // CO2, H2O are in the SIS
602  Q.push(CO2Id_);
603  speciesNumber++;
604  this->activeSpecies_[CO2Id_] = true;
605  Rvalue[CO2Id_] = 1.0;
606 
607  Q.push(H2OId_);
608  speciesNumber++;
609  this->activeSpecies_[H2OId_] = true;
610  Rvalue[H2OId_] = 1.0;
611  if (forceFuelInclusion_)
612  {
613  forAll(fuelSpeciesID_,i)
614  {
615  Q.push(fuelSpeciesID_[i]);
616  speciesNumber++;
617  this->activeSpecies_[fuelSpeciesID_[i]] = true;
618  Rvalue[fuelSpeciesID_[i]] = 1.0;
619  }
620  }
621  }
622 
623  if (T>NOxThreshold_ && NOId_!=-1)
624  {
625  Q.push(NOId_);
626  speciesNumber++;
627  this->activeSpecies_[NOId_] = true;
628  Rvalue[NOId_] = 1.0;
629  }
630  }
631  else // No automaticSIS => all species of the SIS are added
632  {
633  for (label i=0; i<SIS.size(); i++)
634  {
635  label q = SIS[i];
636  this->activeSpecies_[q] = true;
637  speciesNumber++;
638  Q.push(q);
639  Rvalue[q] = 1.0;
640  }
641  }
642 
643  // Execute the main loop for R-value
644  while (!Q.empty())
645  {
646  label u = Q.pop();
647  scalar Den = max(PA[u],CA[u]);
648  if (Den != 0)
649  {
650  for (label v=0; v<NbrABInit[u]; v++)
651  {
652  label otherSpec = rABOtherSpec(u, v);
653  scalar rAB = mag(rABNum(u, v))/Den;
654  if (rAB > 1)
655  {
656  rAB = 1;
657  }
658  // The direct link is weaker than the user-defined tolerance
659  if (rAB >= this->tolerance())
660  {
661  scalar Rtemp = Rvalue[u]*rAB;
662  // a link analysed previously is stronger
663  // the (composed) link is stronger than the user-defined
664  // tolerance
665  if ((Rvalue[otherSpec]<Rtemp) && (Rtemp>=this->tolerance()))
666  {
667  Q.push(otherSpec);
668  Rvalue[otherSpec] = Rtemp;
669  if (!this->activeSpecies_[otherSpec])
670  {
671  this->activeSpecies_[otherSpec] = true;
672  speciesNumber++;
673  }
674  }
675  }
676  }
677  }
678  }
679 
680  // Put a flag on the reactions containing at least one removed species
681  forAll(this->chemistry_.reactions(), i)
682  {
683  const Reaction<ThermoType>& R = this->chemistry_.reactions()[i];
684  this->chemistry_.reactionsDisabled()[i] = false;
685 
686  forAll(R.lhs(), s)
687  {
688  label ss = R.lhs()[s].index;
689  if (!this->activeSpecies_[ss])
690  {
691  // Flag the reaction to disable it
692  this->chemistry_.reactionsDisabled()[i] = true;
693  break;
694  }
695  }
696  if (!this->chemistry_.reactionsDisabled()[i])
697  {
698  forAll(R.rhs(), s)
699  {
700  label ss = R.rhs()[s].index;
701  if (!this->activeSpecies_[ss])
702  {
703  // Flag the reaction to disable it
704  this->chemistry_.reactionsDisabled()[i] = true;
705  break;
706  }
707  }
708  }
709  }
710 
711  this->NsSimp_ = speciesNumber;
712  scalarField& simplifiedC(this->chemistry_.simplifiedC());
713  simplifiedC.setSize(this->NsSimp_+2);
714  DynamicList<label>& s2c(this->chemistry_.simplifiedToCompleteIndex());
715  s2c.setSize(this->NsSimp_);
716  Field<label>& c2s(this->chemistry_.completeToSimplifiedIndex());
717 
718  label j = 0;
719  for (label i=0; i<this->nSpecie_; i++)
720  {
721  if (this->activeSpecies_[i])
722  {
723  s2c[j] = i;
724  simplifiedC[j] = c[i];
725  c2s[i] = j++;
726  if (!this->chemistry_.active(i))
727  {
728  this->chemistry_.setActive(i);
729  }
730  }
731  else
732  {
733  c2s[i] = -1;
734  }
735  }
736 
737  simplifiedC[this->NsSimp_] = T;
738  simplifiedC[this->NsSimp_+1] = p;
739  this->chemistry_.setNsDAC(this->NsSimp_);
740 
741  // Change temporary Ns in chemistryModel
742  // to make the function nEqns working
743  this->chemistry_.setNSpecie(this->NsSimp_);
744 }
745 
746 
747 // ************************************************************************* //
A FIFO stack based on a singly-linked list.
Definition: FIFOStack.H:51
virtual label nSpecie() const
The number of species.
dictionary dict
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:667
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
virtual void reduceMechanism(const scalar p, const scalar T, const scalarField &c, const label li)
Reduce the mechanism.
Definition: DAC.C:273
const dimensionedScalar & c1
First radiation constant: default SI units: [W/m^2].
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
const PtrList< volScalarField > & Y()
const List< specieCoeffs > & lhs() const
Return the components of the left hand side.
Definition: ReactionI.H:57
BasicChemistryModel< rhoReactionThermo > & chemistry
wordList toc() const
Return the table of contents.
Definition: dictionary.C:1018
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:52
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:934
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.
Definition: DynamicListI.H:175
A class for handling words, derived from string.
Definition: word.H:59
DAC(const IOdictionary &dict, TDACChemistryModel< CompType, ThermoType > &chemistry)
Construct from components.
Definition: DAC.C:32
const word & name() const
Return the name of the element.
label nAtoms() const
Return the number of atoms of this element in the specie.
List< List< specieElement > > & specieComp()
const volScalarField & T
virtual ~DAC()
Destructor.
Definition: DAC.C:264
void setSize(const label)
Reset size of List.
Definition: List.C:281
#define R(A, B, C, D, E, F, K, M)
void push(const T &a)
Push an element onto the stack.
Definition: FIFOStack.H:84
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
T pop()
Pop the bottom element off the stack.
Definition: FIFOStack.H:90
const List< specieCoeffs > & rhs() const
Return the components of the right hand side.
Definition: ReactionI.H:64
An abstract class for methods of chemical mechanism reduction.
volScalarField & p
bool found
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:812