chemPointISAT.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-2022 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 "chemPointISAT.H"
27 #include "ISAT.H"
28 #include "odeChemistryModel.H"
29 #include "SVD.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 // Defined as static to be able to dynamically change it during simulations
34 // (all chemPoints refer to the same object)
35 Foam::scalar Foam::chemPointISAT::tolerance_;
36 
37 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
38 
39 void Foam::chemPointISAT::qrDecompose
40 (
41  const label nCols,
43 )
44 {
45  scalarField c(nCols);
46  scalarField d(nCols);
47  scalar scale, sigma, sum;
48 
49  for (label k=0; k<nCols-1; k++)
50  {
51  scale = 0;
52  for (label i=k; i<nCols; i++)
53  {
54  scale=max(scale, mag(R(i, k)));
55  }
56  if (scale == 0)
57  {
58  c[k] = d[k] = 0;
59  }
60  else
61  {
62  for (label i=k; i<nCols; i++)
63  {
64  R(i, k) /= scale;
65  }
66  sum = 0;
67  for (label i=k; i<nCols; i++)
68  {
69  sum += sqr(R(i, k));
70  }
71  sigma = sign(R(k, k))*sqrt(sum);
72  R(k, k) += sigma;
73  c[k] = sigma*R(k, k);
74  d[k] = -scale*sigma;
75  for (label j=k+1; j<nCols; j++)
76  {
77  sum = 0;
78  for ( label i=k; i<nCols; i++)
79  {
80  sum += R(i, k)*R(i, j);
81  }
82  scalar tau = sum/c[k];
83  for ( label i=k; i<nCols; i++)
84  {
85  R(i, j) -= tau*R(i, k);
86  }
87  }
88  }
89  }
90 
91  d[nCols-1] = R(nCols-1, nCols-1);
92 
93  // form R
94  for (label i=0; i<nCols; i++)
95  {
96  R(i, i) = d[i];
97  for ( label j=0; j<i; j++)
98  {
99  R(i, j)=0;
100  }
101  }
102 }
103 
104 
105 void Foam::chemPointISAT::qrUpdate
106 (
108  const label n,
109  const scalarField& u,
110  const scalarField& v
111 )
112 {
113  label k;
114 
115  scalarField w(u);
116  for (k=n-1; k>=0; k--)
117  {
118  if (w[k] != 0)
119  {
120  break;
121  }
122  }
123 
124  if (k < 0)
125  {
126  k = 0;
127  }
128 
129  for (label i=k-1; i>=0; i--)
130  {
131  rotate(R, i, w[i], -w[i+1], n);
132  if (w[i] == 0)
133  {
134  w[i] = mag(w[i+1]);
135  }
136  else if (mag(w[i]) > mag(w[i+1]))
137  {
138  w[i] = mag(w[i])*sqrt(1.0 + sqr(w[i+1]/w[i]));
139  }
140  else
141  {
142  w[i] = mag(w[i+1])*sqrt(1.0 + sqr(w[i]/w[i+1]));
143  }
144  }
145 
146  for (label i=0; i<n; i++)
147  {
148  R(0, i) += w[0]*v[i];
149  }
150 
151  for (label i=0; i<k; i++)
152  {
153  rotate(R, i, R(i, i), -R(i+1, i), n);
154  }
155 }
156 
157 
158 void Foam::chemPointISAT::rotate
159 (
161  const label i,
162  const scalar a,
163  const scalar b,
164  label n
165 )
166 {
167  scalar c, fact, s, w, y;
168 
169  if (a == 0)
170  {
171  c = 0;
172  s = (b >= 0 ? 1 : -1);
173  }
174  else if (mag(a) > mag(b))
175  {
176  fact = b/a;
177  c = sign(a)/sqrt(1.0 + sqr(fact));
178  s = fact*c;
179  }
180  else
181  {
182  fact = a/b;
183  s = sign(b)/sqrt(1.0 + sqr(fact));
184  c = fact*s;
185  }
186 
187  for (label j=i;j<n;j++)
188  {
189  y = R(i, j);
190  w = R(i+1, j);
191  R(i, j) = c*y-s*w;
192  R(i+1, j) = s*y+c*w;
193  }
194 }
195 
196 
197 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
198 
200 (
202  const scalarField& phi,
203  const scalarField& Rphi,
204  const scalarSquareMatrix& A,
205  const scalarField& scaleFactor,
206  const scalar tolerance,
207  const label completeSpaceSize,
208  const label nActive,
209  const dictionary& coeffsDict,
210  binaryNode* node
211 )
212 :
213  table_(table),
214  phi_(phi),
215  Rphi_(Rphi),
216  A_(A),
217  scaleFactor_(scaleFactor),
218  node_(node),
219  completeSpaceSize_(completeSpaceSize),
220  nGrowth_(0),
221  nActive_(nActive),
222  timeTag_(table.timeSteps()),
223  lastTimeUsed_(table.timeSteps()),
224  toRemove_(false),
225  maxNumNewDim_(coeffsDict.lookupOrDefault("maxNumNewDim",0)),
226  printProportion_(coeffsDict.lookupOrDefault("printProportion",false)),
227  numRetrieve_(0),
228  nLifeTime_(0)
229 {
230  tolerance_ = tolerance;
231 
232  iddeltaT_ = completeSpaceSize - 1;
233  scaleFactor_[iddeltaT_] *= phi_[iddeltaT_]/tolerance_;
234 
235  idT_ = completeSpaceSize - 3;
236  idp_ = completeSpaceSize - 2;
237 
238  if (table_.reduction())
239  {
240  simplifiedToCompleteIndex_.setSize(nActive_);
241  completeToSimplifiedIndex_.setSize(completeSpaceSize - 3);
242 
243  forAll(completeToSimplifiedIndex_, i)
244  {
245  completeToSimplifiedIndex_[i] = table_.chemistry().cTos(i);
246  }
247 
248  forAll(simplifiedToCompleteIndex_, i)
249  {
250  simplifiedToCompleteIndex_[i] = table_.chemistry().sToc(i);
251  }
252  }
253 
254  const label reduOrCompDim =
255  table_.reduction() ? nActive_ + 3 : completeSpaceSize;
256 
257  // SVD decomposition A = U*D*V^T
258  SVD svdA(A);
259 
260  scalarDiagonalMatrix D(reduOrCompDim);
261  const scalarDiagonalMatrix& S = svdA.S();
262 
263  // Replace the value of vector D by max(D, 1/2), first ISAT paper
264  for (label i=0; i<reduOrCompDim; i++)
265  {
266  D[i] = max(S[i], 0.5);
267  }
268 
269  // Rebuild A with max length, tol and scale factor before QR decomposition
270  scalarRectangularMatrix Atilde(reduOrCompDim);
271 
272  // Result stored in Atilde
273  multiply(Atilde, svdA.U(), D, svdA.V().T());
274 
275  for (label i=0; i<reduOrCompDim; i++)
276  {
277  for (label j=0; j<reduOrCompDim; j++)
278  {
279  label compi = i;
280 
281  if (table_.reduction())
282  {
283  compi = simplifiedToCompleteIndex(i);
284  }
285 
286  // SF*A/tolerance
287  // (where SF is diagonal with inverse of scale factors)
288  // SF*A is the same as dividing each line by the scale factor
289  // corresponding to the species of this line
290  Atilde(i, j) /= (tolerance*scaleFactor[compi]);
291  }
292  }
293 
294  // The object LT_ (the transpose of the Q) describe the EOA, since we have
295  // A^T B^T B A that should be factorised into L Q^T Q L^T and is set in the
296  // qrDecompose function
297  LT_ = scalarSquareMatrix(Atilde);
298 
299  qrDecompose(reduOrCompDim, LT_);
300 }
301 
302 
304 (
306 )
307 :
308  table_(p.table_),
309  phi_(p.phi()),
310  Rphi_(p.Rphi()),
311  LT_(p.LT()),
312  A_(p.A()),
313  scaleFactor_(p.scaleFactor()),
314  node_(p.node()),
315  completeSpaceSize_(p.completeSpaceSize()),
316  nGrowth_(p.nGrowth()),
317  nActive_(p.nActive()),
318  simplifiedToCompleteIndex_(p.simplifiedToCompleteIndex()),
319  timeTag_(p.timeTag()),
320  lastTimeUsed_(p.lastTimeUsed()),
321  toRemove_(p.toRemove()),
322  maxNumNewDim_(p.maxNumNewDim()),
323  numRetrieve_(0),
324  nLifeTime_(0),
325  completeToSimplifiedIndex_(p.completeToSimplifiedIndex())
326 {
327  tolerance_ = p.tolerance();
328 
329  idT_ = completeSpaceSize() - 3;
330  idp_ = completeSpaceSize() - 2;
331  iddeltaT_ = completeSpaceSize() - 1;
332 }
333 
334 
335 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
336 
338 {
339  const scalarField dphi(phiq - phi());
340 
341  const label dim =
342  table_.reduction() ? nActive_ : completeSpaceSize() - 3;
343 
344  scalar epsTemp = 0;
345  List<scalar> propEps(completeSpaceSize(), scalar(0));
346 
347  for (label i=0; i<completeSpaceSize()-3; i++)
348  {
349  scalar temp = 0;
350 
351  // When mechanism reduction is inactive OR on active species multiply L
352  // by dphi to get the distance in the active species direction else (for
353  // inactive species), just multiply the diagonal element and dphi
354  if
355  (
356  !(table_.reduction())
357  ||(table_.reduction() && completeToSimplifiedIndex_[i] != -1)
358  )
359  {
360  const label si =
361  table_.reduction()
362  ? completeToSimplifiedIndex_[i]
363  : i;
364 
365  for (label j=si; j<dim; j++)// LT is upper triangular
366  {
367  const label sj =
368  table_.reduction()
369  ? simplifiedToCompleteIndex_[j]
370  : j;
371 
372  temp += LT_(si, j)*dphi[sj];
373  }
374 
375  temp += LT_(si, dim)*dphi[idT_];
376  temp += LT_(si, dim+1)*dphi[idp_];
377  temp += LT_(si, dim+2)*dphi[iddeltaT_];
378  }
379  else
380  {
381  temp = dphi[i]/(tolerance_*scaleFactor_[i]);
382  }
383 
384  epsTemp += sqr(temp);
385 
386  if (printProportion_)
387  {
388  propEps[i] = temp;
389  }
390  }
391 
392  // Temperature
393  epsTemp +=
394  sqr
395  (
396  LT_(dim, dim)*dphi[idT_]
397  +LT_(dim, dim+1)*dphi[idp_]
398  +LT_(dim, dim+2)*dphi[iddeltaT_]
399  );
400 
401  // Pressure
402  epsTemp +=
403  sqr
404  (
405  LT_(dim+1, dim+1)*dphi[idp_]
406  +LT_(dim+1, dim+2)*dphi[iddeltaT_]
407  );
408 
409  epsTemp += sqr(LT_[dim+2][dim+2]*dphi[iddeltaT_]);
410 
411  if (printProportion_)
412  {
413  propEps[idT_] = sqr
414  (
415  LT_(dim, dim)*dphi[idT_]
416  + LT_(dim, dim+1)*dphi[idp_]
417  );
418 
419  propEps[idp_] =
420  sqr(LT_(dim+1, dim+1)*dphi[idp_]);
421 
422  propEps[iddeltaT_] =
423  sqr(LT_[dim+2][dim+2]*dphi[iddeltaT_]);
424  }
425 
426  if (sqrt(epsTemp) > 1 + tolerance_)
427  {
428  if (printProportion_)
429  {
430  scalar max = -1;
431  label maxIndex = -1;
432  for (label i=0; i<completeSpaceSize(); i++)
433  {
434  if(max < propEps[i])
435  {
436  max = propEps[i];
437  maxIndex = i;
438  }
439  }
440  word propName;
441  if (maxIndex >= completeSpaceSize() - 3)
442  {
443  if (maxIndex == idT_)
444  {
445  propName = "T";
446  }
447  else if (maxIndex == idp_)
448  {
449  propName = "p";
450  }
451  else if (maxIndex == iddeltaT_)
452  {
453  propName = "deltaT";
454  }
455  }
456  else
457  {
458  propName = table_.chemistry().Y()[maxIndex].member();
459  }
460 
461  Info<< "Direction maximum impact to error in ellipsoid: "
462  << propName << endl;
463 
464  Info<< "Proportion to the total error on the retrieve: "
465  << max/(epsTemp+small) << endl;
466  }
467  return false;
468  }
469  else
470  {
471  return true;
472  }
473 }
474 
475 
477 (
478  const scalarField& phiq,
479  const scalarField& Rphiq
480 )
481 {
482  scalar eps2 = 0;
483  const scalarField dR(Rphiq - Rphi());
484  const scalarField dphi(phiq - phi());
485  const scalarField& scaleFactorV(scaleFactor());
486  const scalarSquareMatrix& Avar(A());
487  scalar dRl = 0;
488 
489  const label dim =
490  table_.reduction() ? nActive_ : completeSpaceSize() - 2;
491 
492  // Since we build only the solution for the species, T and p are not
493  // included
494  for (label i=0; i<completeSpaceSize()-3; i++)
495  {
496  dRl = 0;
497  if (table_.reduction())
498  {
499  const label si = completeToSimplifiedIndex_[i];
500 
501  // If this species is active
502  if (si != -1)
503  {
504  for (label j=0; j<dim; j++)
505  {
506  const label sj = simplifiedToCompleteIndex_[j];
507  dRl += Avar(si, j)*dphi[sj];
508  }
509  dRl += Avar(si, nActive_)*dphi[idT_];
510  dRl += Avar(si, nActive_+1)*dphi[idp_];
511  dRl += Avar(si, nActive_+2)*dphi[iddeltaT_];
512  }
513  else
514  {
515  dRl = dphi[i];
516  }
517  }
518  else
519  {
520  for (label j=0; j<completeSpaceSize(); j++)
521  {
522  dRl += Avar(i, j)*dphi[j];
523  }
524  }
525  eps2 += sqr((dR[i]-dRl)/scaleFactorV[i]);
526  }
527 
528  eps2 = sqrt(eps2);
529  if (eps2 > tolerance())
530  {
531  return false;
532  }
533  else
534  {
535  // if the solution is in the ellipsoid of accuracy
536  return true;
537  }
538 }
539 
540 
542 {
543  const scalarField dphi(phiq - phi());
544  const label initNActiveSpecies(nActive_);
545 
546  if (table_.reduction())
547  {
548  label activeAdded(0);
549  DynamicList<label> dimToAdd(0);
550 
551  // check if the difference of active species is lower than the maximum
552  // number of new dimensions allowed
553  for (label i=0; i<completeSpaceSize()-3; i++)
554  {
555  // first test if the current chemPoint has an inactive species
556  // corresponding to an active one in the query point
557  if
558  (
559  completeToSimplifiedIndex_[i] == -1
560  && table_.chemistry().cTos(i) != -1
561  )
562  {
563  activeAdded++;
564  dimToAdd.append(i);
565  }
566  // then test if an active species in the current chemPoint
567  // corresponds to an inactive on the query side
568  if
569  (
570  completeToSimplifiedIndex_[i] != -1
571  && table_.chemistry().cTos(i) == -1
572  )
573  {
574  activeAdded++;
575  // we don't need to add a new dimension but we count it to have
576  // control on the difference through maxNumNewDim
577  }
578  // finally test if both points have inactive species but
579  // with a dphi!=0
580  if
581  (
582  completeToSimplifiedIndex_[i] == -1
583  && table_.chemistry().cTos(i) == -1
584  && dphi[i] != 0
585  )
586  {
587  activeAdded++;
588  dimToAdd.append(i);
589  }
590  }
591 
592  // if the number of added dimension is too large, growth fail
593  if (activeAdded > maxNumNewDim_)
594  {
595  return false;
596  }
597 
598  // the number of added dimension to the current chemPoint
599  nActive_ += dimToAdd.size();
600  simplifiedToCompleteIndex_.setSize(nActive_);
601  forAll(dimToAdd, i)
602  {
603  label si = nActive_ - dimToAdd.size() + i;
604  // add the new active species
605  simplifiedToCompleteIndex_[si] = dimToAdd[i];
606  completeToSimplifiedIndex_[dimToAdd[i]] = si;
607  }
608 
609  // update LT and A :
610  //-add new column and line for the new active species
611  //-transfer last two lines of the previous matrix (p and T) to the end
612  // (change the diagonal position)
613  //-set all element of the new lines and columns to zero except diagonal
614  // (=1/(tolerance*scaleFactor))
615  if (nActive_ > initNActiveSpecies)
616  {
617  const scalarSquareMatrix LTvar = LT_; // take a copy of LT_
618  const scalarSquareMatrix Avar = A_; // take a copy of A_
619  LT_ = scalarSquareMatrix(nActive_+3, Zero);
620  A_ = scalarSquareMatrix(nActive_+3, Zero);
621 
622  // write the initial active species
623  for (label i=0; i<initNActiveSpecies; i++)
624  {
625  for (label j=0; j<initNActiveSpecies; j++)
626  {
627  LT_(i, j) = LTvar(i, j);
628  A_(i, j) = Avar(i, j);
629  }
630  }
631 
632  // write the columns for temperature and pressure
633  for (label i=0; i<initNActiveSpecies; i++)
634  {
635  for (label j=1; j>=0; j--)
636  {
637  LT_(i, nActive_+j)=LTvar(i, initNActiveSpecies+j);
638  A_(i, nActive_+j)=Avar(i, initNActiveSpecies+j);
639  LT_(nActive_+j, i)=LTvar(initNActiveSpecies+j, i);
640  A_(nActive_+j, i)=Avar(initNActiveSpecies+j, i);
641  }
642  }
643  // end with the diagonal elements for temperature and pressure
644  LT_(nActive_, nActive_)=
645  LTvar(initNActiveSpecies, initNActiveSpecies);
646  A_(nActive_, nActive_)=
647  Avar(initNActiveSpecies, initNActiveSpecies);
648  LT_(nActive_+1, nActive_+1)=
649  LTvar(initNActiveSpecies+1, initNActiveSpecies+1);
650  A_(nActive_+1, nActive_+1)=
651  Avar(initNActiveSpecies+1, initNActiveSpecies+1);
652  LT_(nActive_+2, nActive_+2)=
653  LTvar(initNActiveSpecies+2, initNActiveSpecies+2);
654  A_(nActive_+2, nActive_+2)=
655  Avar(initNActiveSpecies+2, initNActiveSpecies+2);
656 
657  for (label i=initNActiveSpecies; i<nActive_;i++)
658  {
659  LT_(i, i)=
660  1.0
661  /(tolerance_*scaleFactor_[simplifiedToCompleteIndex_[i]]);
662  A_(i, i) = 1;
663  }
664  }
665  }
666 
667  const label dim =
668  table_.reduction() ? nActive_ + 3 : completeSpaceSize();
669 
670  // beginning of grow algorithm
671  scalarField phiTilde(dim, 0);
672  scalar normPhiTilde = 0;
673  // p' = L^T.(p-phi)
674 
675  for (label i=0; i<dim; i++)
676  {
677  for (label j=i; j<dim-3; j++)// LT is upper triangular
678  {
679  const label sj =
680  table_.reduction()
681  ? simplifiedToCompleteIndex_[j]
682  : j;
683 
684  phiTilde[i] += LT_(i, j)*dphi[sj];
685  }
686 
687  phiTilde[i] += LT_(i, dim-3)*dphi[idT_];
688  phiTilde[i] += LT_(i, dim-3+1)*dphi[idp_];
689  phiTilde[i] += LT_(i, dim-3+2)*dphi[iddeltaT_];
690 
691  normPhiTilde += sqr(phiTilde[i]);
692  }
693 
694  const scalar invSqrNormPhiTilde = 1.0/normPhiTilde;
695  normPhiTilde = sqrt(normPhiTilde);
696 
697  // gamma = (1/|p'| - 1)/|p'|^2
698  const scalar gamma = (1/normPhiTilde - 1)*invSqrNormPhiTilde;
699  const scalarField u(gamma*phiTilde);
700  scalarField v(dim, 0);
701 
702  for (label i=0; i<dim; i++)
703  {
704  for (label j=0; j<=i;j++)
705  {
706  v[i] += phiTilde[j]*LT_(j, i);
707  }
708  }
709 
710  qrUpdate(LT_,dim, u, v);
711  nGrowth_++;
712 
713  return true;
714 }
715 
716 
717 // ************************************************************************* //
dimensionedScalar sign(const dimensionedScalar &ds)
const label & timeTag()
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
void multiply(FieldField< Field, Type > &f, const FieldField< Field, Type > &f1, const FieldField< Field, scalar > &f2)
FvWallInfoData< WallInfo, label > label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
label cTos(const label ci) const
Return the index in the simplified set of species.
const List< label > & simplifiedToCompleteIndex() const
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
const scalarRectangularMatrix & U() const
Return U.
Definition: SVDI.H:38
Form T() const
Return the transpose of the matrix.
Definition: Matrix.C:264
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const PtrList< volScalarField > & Y() const
Return a reference to the list of mass fraction fields.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
bool checkSolution(const scalarField &phiq, const scalarField &Rphiq)
If phiq is not in the EOA, then the mapping is computed.
const dimensionedScalar b
Wien displacement law constant: default SI units: [m K].
Definition: createFields.H:27
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Leaf of the binary tree. The chemPoint stores the composition &#39;phi&#39;, the mapping of this composition ...
const scalar & tolerance() const
label nActive() const
label k
Boltzmann constant.
binaryNode *& node()
const dimensionedScalar c
Speed of light in a vacuum.
label sToc(const label si) const
Return the index in the complete set of species.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
const odeChemistryModel & chemistry()
Definition: ISAT.H:246
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m^2/K^4].
const List< label > & completeToSimplifiedIndex() const
scalar y
const scalarDiagonalMatrix & S() const
Return the singular values.
Definition: SVDI.H:50
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))
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
Definition: DynamicList.H:56
A class for handling words, derived from string.
Definition: word.H:59
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
const scalarField & phi() const
tmp< fvMatrix< Type > > S(const Pair< tmp< volScalarField::Internal >> &, const GeometricField< Type, fvPatchField, volMesh > &)
label nGrowth() const
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:296
bool grow(const scalarField &phiq)
More details about the minimum-volume ellipsoid covering an.
Singular value decomposition of a rectangular matrix.
Definition: SVD.H:51
static const zero Zero
Definition: zero.H:97
Node of the binary tree.
Definition: binaryNode.H:48
label timeSteps() const
Return the number of chemistry evaluations.
Definition: ISAT.H:262
chemPointISAT(chemistryTabulationMethods::ISAT &table, const scalarField &phi, const scalarField &Rphi, const scalarSquareMatrix &A, const scalarField &scaleFactor, const scalar tolerance, const label completeSpaceSize, const label nActive, const dictionary &coeffsDict, binaryNode *node=nullptr)
Construct from components.
const scalarRectangularMatrix & V() const
Return the square matrix V.
Definition: SVDI.H:44
const scalarField & Rphi() const
const scalarSquareMatrix & LT() const
void setSize(const label)
Reset size of List.
Definition: List.C:281
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
#define R(A, B, C, D, E, F, K, M)
bool inEOA(const scalarField &phiq)
To RETRIEVE the mapping from the stored chemPoint phi, the query.
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
label n
Implementation of the ISAT (In-situ adaptive tabulation), for chemistry calculation.
Definition: ISAT.H:60
const scalarSquareMatrix & A() const
const scalarField & scaleFactor() const
SquareMatrix< scalar > scalarSquareMatrix
bool reduction() const
Return true if reduction is applied to the state variables.
Definition: ISAT.H:241
label & completeSpaceSize()