MaxwellStefan.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) 2021-2024 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 "MaxwellStefan.H"
27 #include "fvcDiv.H"
28 #include "fvcLaplacian.H"
29 #include "fvcSnGrad.H"
30 #include "fvmSup.H"
31 #include "surfaceInterpolate.H"
32 #include "Function2Evaluate.H"
33 
34 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38 
39 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
40 
41 template<class BasicThermophysicalTransportModel>
42 void MaxwellStefan<BasicThermophysicalTransportModel>::
43 transformDiffusionCoefficient() const
44 {
45  const label d = this->thermo().defaultSpecie();
46 
47  // Calculate the molecular weight of the mixture and the mole fractions
48  scalar Wm = 0;
49 
50  forAll(W, i)
51  {
52  X[i] = Y[i]/W[i];
53  Wm += X[i];
54  }
55 
56  Wm = 1/Wm;
57  X *= Wm;
58 
59  // i counter for the specie sub-system without the default specie
60  label is = 0;
61 
62  // Calculate the A and B matrices from the binary mass diffusion
63  // coefficients and specie mole fractions
64  forAll(X, i)
65  {
66  if (i != d)
67  {
68  A(is, is) = -X[i]*Wm/(DD(i, d)*W[d]);
69  B(is, is) = -(X[i]*Wm/W[d] + (1 - X[i])*Wm/W[i]);
70 
71  // j counter for the specie sub-system without the default specie
72  label js = 0;
73 
74  forAll(X, j)
75  {
76  if (j != i)
77  {
78  A(is, is) -= X[j]*Wm/(DD(i, j)*W[i]);
79 
80  if (j != d)
81  {
82  A(is, js) =
83  X[i]*(Wm/(DD(i, j)*W[j]) - Wm/(DD(i, d)*W[d]));
84 
85  B(is, js) = X[i]*(Wm/W[j] - Wm/W[d]);
86  }
87  }
88 
89  if (j != d)
90  {
91  js++;
92  }
93  }
94 
95  is++;
96  }
97  }
98 
99  // LU decompose A and invert
100  A.decompose();
101  A.inv(invA);
102 
103  // Calculate the generalised Fick's law diffusion coefficients
104  multiply(D, invA, B);
105 }
106 
107 
108 template<class BasicThermophysicalTransportModel>
109 void MaxwellStefan<BasicThermophysicalTransportModel>::
110 transformDiffusionCoefficientFields() const
111 {
112  const label d = this->thermo().defaultSpecie();
113 
114  // For each cell or patch face
115  forAll(*(YPtrs[0]), pi)
116  {
117  forAll(W, i)
118  {
119  // Map YPtrs -> Y
120  Y[i] = (*YPtrs[i])[pi];
121 
122  // Map DijPtrs -> DD
123  forAll(W, j)
124  {
125  DD(i, j) = (*DijPtrs[i][j])[pi];
126  }
127  }
128 
129  // Transform DD -> D
130  transformDiffusionCoefficient();
131 
132  // i counter for the specie sub-system without the default specie
133  label is = 0;
134 
135  forAll(W, i)
136  {
137  if (i != d)
138  {
139  // j counter for the specie sub-system
140  // without the default specie
141  label js = 0;
142 
143  // Map D -> DijPtrs
144  forAll(W, j)
145  {
146  if (j != d)
147  {
148  (*DijPtrs[i][j])[pi] = D(is, js);
149 
150  js++;
151  }
152  }
153 
154  is++;
155  }
156  }
157  }
158 }
159 
160 
161 template<class BasicThermophysicalTransportModel>
162 void MaxwellStefan<BasicThermophysicalTransportModel>::transform
163 (
164  List<PtrList<volScalarField>>& Dij
165 ) const
166 {
167  const PtrList<volScalarField>& Y = this->thermo().Y();
168  const volScalarField& Y0 = Y[0];
169 
170  forAll(W, i)
171  {
172  // Map this->thermo().Y() internal fields -> YPtrs
173  YPtrs[i] = &Y[i].primitiveField();
174 
175  // Map Dii_ internal fields -> DijPtrs
176  DijPtrs[i][i] = &Dii_[i].primitiveFieldRef();
177 
178  // Map Dij internal fields -> DijPtrs
179  forAll(W, j)
180  {
181  if (j != i)
182  {
183  DijPtrs[i][j] = &Dij[i][j].primitiveFieldRef();
184  }
185  }
186  }
187 
188  // Transform binary mass diffusion coefficients internal field DijPtrs ->
189  // generalised Fick's law diffusion coefficients DijPtrs
190  transformDiffusionCoefficientFields();
191 
192  forAll(Y0.boundaryField(), patchi)
193  {
194  forAll(W, i)
195  {
196  // Map this->thermo().Y() patch fields -> YPtrs
197  YPtrs[i] = &Y[i].boundaryField()[patchi];
198 
199  // Map Dii_ patch fields -> DijPtrs
200  DijPtrs[i][i] = &Dii_[i].boundaryFieldRef()[patchi];
201 
202  // Map Dij patch fields -> DijPtrs
203  forAll(W, j)
204  {
205  if (j != i)
206  {
207  DijPtrs[i][j] = &Dij[i][j].boundaryFieldRef()[patchi];
208  }
209  }
210  }
211 
212  // Transform binary mass diffusion coefficients patch field DijPtrs ->
213  // generalised Fick's law diffusion coefficients DijPtrs
214  transformDiffusionCoefficientFields();
215  }
216 }
217 
218 
219 template<class BasicThermophysicalTransportModel>
221 {
222  const label d = this->thermo().defaultSpecie();
223 
224  const PtrList<volScalarField>& Y = this->thermo().Y();
225  const volScalarField& p = this->thermo().p();
226  const volScalarField& T = this->thermo().T();
227  const volScalarField& rho = this->momentumTransport().rho();
228 
229  Dii_.setSize(Y.size());
230  jexp_.setSize(Y.size());
231 
232  List<PtrList<volScalarField>> Dij(Y.size());
233 
234  // Evaluate the specie binary mass diffusion coefficient functions
235  // and initialise the explicit part of the specie mass flux fields
236  forAll(Y, i)
237  {
238  if (i != d)
239  {
240  if (jexp_.set(i))
241  {
242  jexp_[i] = Zero;
243  }
244  else
245  {
246  jexp_.set
247  (
248  i,
250  (
251  "jexp" + Y[i].name(),
252  Y[i].mesh(),
253  dimensionedScalar(dimensionSet(1, -2, -1, 0, 0), 0)
254  )
255  );
256  }
257  }
258 
259  Dii_.set(i, evaluate(DFuncs_[i][i], dimKinematicViscosity, p, T));
260 
261  Dij[i].setSize(Y.size());
262 
263  forAll(Y, j)
264  {
265  if (j > i)
266  {
267  Dij[i].set
268  (
269  j,
270  evaluate(DFuncs_[i][j], dimKinematicViscosity, p, T)
271  );
272  }
273  else if (j < i)
274  {
275  Dij[i].set(j, Dij[j][i].clone());
276  }
277  }
278  }
279 
280  //- Transform the binary mass diffusion coefficients into the
281  // the generalised Fick's law diffusion coefficients
282  transform(Dij);
283 
284  // Accumulate the explicit part of the specie mass flux fields
285  forAll(Y, j)
286  {
287  if (j != d)
288  {
289  const surfaceScalarField snGradYj(fvc::snGrad(Y[j]));
290 
291  forAll(Y, i)
292  {
293  if (i != d && i != j)
294  {
295  jexp_[i] -= fvc::interpolate(rho*Dij[i][j])*snGradYj;
296  }
297  }
298  }
299  }
300 
301  // Optionally add the Soret thermal diffusion contribution to the
302  // explicit part of the specie mass flux fields
303  if (DTFuncs_.size())
304  {
306 
307  forAll(Y, i)
308  {
309  if (i != d)
310  {
311  jexp_[i] -= fvc::interpolate
312  (
313  evaluate(DTFuncs_[i], dimDynamicViscosity, p, T)
314  )*gradTbyT;
315  }
316  }
317  }
318 }
319 
320 
321 template<class BasicThermophysicalTransportModel>
322 const PtrList<volScalarField>&
323 MaxwellStefan<BasicThermophysicalTransportModel>::Dii() const
324 {
325  if (!Dii_.size())
326  {
327  updateDii();
328  }
329 
330  return Dii_;
331 }
332 
333 
334 template<class BasicThermophysicalTransportModel>
335 const PtrList<surfaceScalarField>&
336 MaxwellStefan<BasicThermophysicalTransportModel>::jexp() const
337 {
338  if (!jexp_.size())
339  {
340  updateDii();
341  }
342 
343  return jexp_;
344 }
345 
346 
347 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
348 
349 template<class BasicThermophysicalTransportModel>
351 (
352  const word& type,
353  const momentumTransportModel& momentumTransport,
354  const thermoModel& thermo
355 )
356 :
357  BasicThermophysicalTransportModel
358  (
359  type,
360  momentumTransport,
361  thermo
362  ),
363 
365 
366  DFuncs_(this->thermo().species().size()),
367 
368  DTFuncs_
369  (
370  this->coeffDict_.found("DT")
371  ? this->thermo().species().size()
372  : 0
373  ),
374 
375  W(this->thermo().species().size()),
376 
377  YPtrs(W.size()),
378  DijPtrs(W.size()),
379 
380  Y(W.size()),
381  X(W.size()),
382  DD(W.size()),
383  A(W.size() - 1),
384  B(A.m()),
385  invA(A.m()),
386  D(W.size())
387 {
388  // Set the molecular weights of the species
389  forAll(W, i)
390  {
391  W[i] = this->thermo().Wi(i).value();
392  }
393 }
394 
395 
396 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
397 
398 template<class BasicThermophysicalTransportModel>
400 {
401  if
402  (
404  )
405  {
406  const speciesTable& species = this->thermo().species();
407 
408  const dictionary& Ddict = this->coeffDict_.subDict("D");
409 
410  // Read the array of specie binary mass diffusion coefficient functions
411  forAll(species, i)
412  {
413  DFuncs_[i].setSize(species.size());
414 
415  forAll(species, j)
416  {
417  if (j >= i)
418  {
419  const word nameij(species[i] + '-' + species[j]);
420  const word nameji(species[j] + '-' + species[i]);
421 
422  word Dname;
423 
424  if (Ddict.found(nameij) && Ddict.found(nameji))
425  {
426  if (i != j)
427  {
429  << "Binary mass diffusion coefficients "
430  "for both " << nameij << " and " << nameji
431  << " provided, using " << nameij << endl;
432  }
433 
434  Dname = nameij;
435  }
436  else if (Ddict.found(nameij))
437  {
438  Dname = nameij;
439  }
440  else if (Ddict.found(nameji))
441  {
442  Dname = nameji;
443  }
444  else
445  {
447  << "Binary mass diffusion coefficients for pair "
448  << nameij << " or " << nameji << " not provided"
449  << exit(FatalIOError);
450  }
451 
452  DFuncs_[i].set
453  (
454  j,
456  (
457  Dname,
458  dimPressure,
461  Ddict
462  ).ptr()
463  );
464  }
465  }
466  }
467 
468  // Optionally read the List of specie Soret thermal diffusion
469  // coefficient functions
470  if (this->coeffDict_.found("DT"))
471  {
472  const dictionary& DTdict = this->coeffDict_.subDict("DT");
473 
474  forAll(species, i)
475  {
476  DTFuncs_.set
477  (
478  i,
480  (
481  species[i],
482  dimPressure,
485  DTdict
486  ).ptr()
487  );
488  }
489  }
490 
491  return true;
492  }
493  else
494  {
495  return false;
496  }
497 }
498 
499 
500 template<class BasicThermophysicalTransportModel>
502 (
503  const volScalarField& Yi
504 ) const
505 {
506  return volScalarField::New
507  (
508  "DEff",
509  this->momentumTransport().rho()*Dii()[this->thermo().specieIndex(Yi)]
510  );
511 }
512 
513 
514 template<class BasicThermophysicalTransportModel>
516 (
517  const volScalarField& Yi,
518  const label patchi
519 ) const
520 {
521  return
522  this->momentumTransport().rho().boundaryField()[patchi]
523  *Dii()[this->thermo().specieIndex(Yi)].boundaryField()[patchi];
524 }
525 
526 
527 template<class BasicThermophysicalTransportModel>
530 {
532  (
534  (
536  (
537  "q",
538  this->momentumTransport().alphaRhoPhi().group()
539  ),
540  -fvc::interpolate(this->alpha()*this->kappaEff())
541  *fvc::snGrad(this->thermo().T())
542  )
543  );
544 
545  const label d = this->thermo().defaultSpecie();
546 
547  const PtrList<volScalarField>& Y = this->thermo().Y();
548  const volScalarField& p = this->thermo().p();
549  const volScalarField& T = this->thermo().T();
550 
551  if (Y.size())
552  {
553  surfaceScalarField sumJ
554  (
556  (
557  "sumJ",
558  Y[0].mesh(),
560  )
561  );
562 
563  surfaceScalarField sumJh
564  (
566  (
567  "sumJh",
568  Y[0].mesh(),
570  )
571  );
572 
573  forAll(Y, i)
574  {
575  if (i != d)
576  {
577  const volScalarField hi(this->thermo().hsi(i, p, T));
578 
579  const surfaceScalarField ji(this->j(Y[i]));
580  sumJ += ji;
581 
582  sumJh += ji*fvc::interpolate(hi);
583  }
584  }
585 
586  {
587  const label i = d;
588 
589  const volScalarField hi(this->thermo().hsi(i, p, T));
590 
591  sumJh -= sumJ*fvc::interpolate(hi);
592  }
593 
594  tmpq.ref() += sumJh;
595  }
596 
597  return tmpq;
598 }
599 
600 
601 template<class BasicThermophysicalTransportModel>
603 (
605 ) const
606 {
607  tmp<fvScalarMatrix> tmpDivq
608  (
609  fvm::Su
610  (
611  -fvc::laplacian(this->alpha()*this->kappaEff(), this->thermo().T()),
612  he
613  )
614  );
615 
616  const label d = this->thermo().defaultSpecie();
617 
618  const PtrList<volScalarField>& Y = this->thermo().Y();
619  const volScalarField& p = this->thermo().p();
620  const volScalarField& T = this->thermo().T();
621 
622  tmpDivq.ref() -=
623  fvm::laplacianCorrection(this->alpha()*this->alphaEff(), he);
624 
625  surfaceScalarField sumJ
626  (
628  (
629  "sumJ",
630  he.mesh(),
632  )
633  );
634 
635  surfaceScalarField sumJh
636  (
638  (
639  "sumJh",
640  he.mesh(),
641  dimensionedScalar(sumJ.dimensions()*he.dimensions(), 0)
642  )
643  );
644 
645  forAll(Y, i)
646  {
647  if (i != d)
648  {
649  const volScalarField hi(this->thermo().hsi(i, p, T));
650 
651  const surfaceScalarField ji(this->j(Y[i]));
652  sumJ += ji;
653 
654  sumJh += ji*fvc::interpolate(hi);
655  }
656  }
657 
658  {
659  const label i = d;
660 
661  const volScalarField hi(this->thermo().hsi(i, p, T));
662 
663  sumJh -= sumJ*fvc::interpolate(hi);
664  }
665 
666  tmpDivq.ref() += fvc::div(sumJh*he.mesh().magSf());
667 
668  return tmpDivq;
669 }
670 
671 
672 template<class BasicThermophysicalTransportModel>
674 (
675  const volScalarField& Yi
676 ) const
677 {
678  const label d = this->thermo().defaultSpecie();
679 
680  if (this->thermo().specieIndex(Yi) == d)
681  {
682  const PtrList<volScalarField>& Y = this->thermo().Y();
683 
685  (
687  (
689  (
690  "j" + name(d),
691  this->momentumTransport().alphaRhoPhi().group()
692  ),
693  Yi.mesh(),
695  )
696  );
697 
698  surfaceScalarField& jd = tjd.ref();
699 
700  forAll(Y, i)
701  {
702  if (i != d)
703  {
704  jd -= this->j(Y[i]);
705  }
706  }
707 
708  return tjd;
709  }
710  else
711  {
712  return
713  BasicThermophysicalTransportModel::j(Yi)
714  + jexp()[this->thermo().specieIndex(Yi)];
715  }
716 }
717 
718 
719 template<class BasicThermophysicalTransportModel>
721 (
722  volScalarField& Yi
723 ) const
724 {
725  return
726  BasicThermophysicalTransportModel::divj(Yi)
727  + fvc::div(jexp()[this->thermo().specieIndex(Yi)]*Yi.mesh().magSf());
728 }
729 
730 
731 template<class BasicThermophysicalTransportModel>
733 {
734  BasicThermophysicalTransportModel::predict();
735  updateDii();
736 }
737 
738 
739 template<class BasicThermophysicalTransportModel>
741 {
742  return true;
743 }
744 
745 
746 template<class BasicThermophysicalTransportModel>
748 (
749  const polyTopoChangeMap& map
750 )
751 {
752  // Delete the cached Dii and jexp, will be re-created in predict
753  Dii_.clear();
754  jexp_.clear();
755 }
756 
757 
758 template<class BasicThermophysicalTransportModel>
760 (
761  const polyMeshMap& map
762 )
763 {
764  // Delete the cached Dii and jexp, will be re-created in predict
765  Dii_.clear();
766  jexp_.clear();
767 }
768 
769 
770 template<class BasicThermophysicalTransportModel>
772 (
773  const polyDistributionMap& map
774 )
775 {
776  // Delete the cached Dii and jexp, will be re-created in predict
777  Dii_.clear();
778  jexp_.clear();
779 }
780 
781 
782 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
783 
784 } // End namespace Foam
785 
786 // ************************************************************************* //
static const Foam::dimensionedScalar A("A", Foam::dimPressure, 611.21)
static const Foam::dimensionedScalar B("B", Foam::dimless, 18.678)
static const Foam::dimensionedScalar D("D", Foam::dimTemperature, 257.14)
bool found
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
const dimensionSet & dimensions() const
Return dimensions.
const Mesh & mesh() const
Return mesh.
Run-time selectable function of two variables.
Definition: Function2.H:98
Generic GeometricField class.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &, const HashPtrTable< Source > &=HashPtrTable< Source >())
Return a temporary field constructed from name,.
static word groupName(Name name, const word &group)
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:91
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
Base class for multi-component Maxwell Stefan generalised Fick's law diffusion coefficients based tem...
Definition: MaxwellStefan.H:76
virtual bool movePoints()
Update for mesh motion.
virtual tmp< volScalarField > DEff(const volScalarField &Yi) const
Effective mass diffusion coefficient.
virtual tmp< fvScalarMatrix > divq(volScalarField &he) const
Return the source term for the energy equation.
virtual tmp< fvScalarMatrix > divj(volScalarField &Yi) const
Return the source term for the given specie mass-fraction equation.
virtual void predict()
Update the diffusion coefficients and flux corrections.
virtual void mapMesh(const polyMeshMap &map)
Update from another mesh using the given map.
virtual tmp< surfaceScalarField > j(const volScalarField &Yi) const
Return the specie flux for the given specie mass-fraction [kg/m^2/s].
virtual tmp< surfaceScalarField > q() const
Return the heat flux [W/m^2].
virtual void distribute(const polyDistributionMap &map)
Redistribute or update using the given distribution map.
MaxwellStefan(const word &type, const momentumTransportModel &momentumTransport, const thermoModel &thermo)
Construct from a momentum transport model and a thermo model.
BasicThermophysicalTransportModel::thermoModel thermoModel
BasicThermophysicalTransportModel::momentumTransportModel momentumTransportModel
virtual void topoChange(const polyTopoChangeMap &map)
Update topology using the given map.
virtual bool read()
Read thermophysicalTransport dictionary.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: PtrList.H:75
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
virtual const volScalarField & T() const =0
Temperature [K].
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:843
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:509
Dimension set for the base types.
Definition: dimensionSet.H:125
virtual const volScalarField & p() const =0
Pressure [Pa].
A wordList with hashed indices for faster lookup by name.
virtual const speciesTable & species() const =0
The table of species.
virtual PtrList< volScalarField > & Y()=0
Access the mass-fraction fields.
virtual label defaultSpecie() const =0
The index of the default specie.
label specieIndex(const volScalarField &Yi) const
Access the specie index of the given mass-fraction field.
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
Class containing mesh-to-mesh mapping information.
Definition: polyMeshMap.H:51
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
A class for managing temporary objects.
Definition: tmp.H:55
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
A class for handling words, derived from string.
Definition: word.H:62
volSymmTensorField invA(inv(I *UEqn.A()+drag->Dcu()))
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:346
Calculate the divergence of the given field.
Calculate the laplacian of the given field.
Calculate the snGrad of the given volField.
Calculate the matrix for implicit and explicit sources.
label patchi
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
#define WarningInFunction
Report a warning using Foam::Warning.
void read(Istream &, label &, const dictionary &)
In-place read with dictionary lookup.
const char *const group
Group name for atomic constants.
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< VolField< Type > > laplacian(const VolField< Type > &vf, const word &name)
Definition: fvcLaplacian.C:45
tmp< VolField< Type > > div(const SurfaceField< Type > &ssf)
Definition: fvcDiv.C:47
tmp< SurfaceField< Type > > snGrad(const VolField< Type > &vf, const word &name)
Definition: fvcSnGrad.C:45
tmp< fvMatrix< Type > > Su(const DimensionedField< Type, volMesh > &, const VolField< Type > &)
tmp< fvMatrix< Type > > laplacianCorrection(const VolField< scalar > &gamma, const VolField< Type > &vf)
Definition: fvmLaplacian.C:340
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
static const zero Zero
Definition: zero.H:97
const dimensionSet dimDynamicViscosity
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
const dimensionSet dimEnergy
const dimensionSet dimPressure
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
const dimensionSet dimKinematicViscosity
SurfaceField< scalar > surfaceScalarField
const dimensionSet dimTemperature
const dimensionSet dimTime
scalarList W(const fluidMulticomponentThermo &thermo)
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:64
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:504
T clone(const T &t)
Definition: List.H:55
IOerror FatalIOError
const dimensionSet dimMass
void evaluate(GeometricField< Type, PatchField, GeoMesh > &result, const Function1< Type > &func, const GeometricField< Type, PatchField, GeoMesh > &x)
const dimensionSet dimArea
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
void multiply(FieldField< Field, Type > &f, const FieldField< Field, Type > &f1, const FieldField< Field, scalar > &f2)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
thermo he()
scalarList Y0(nSpecie, 0.0)
volScalarField & p
PtrList< volScalarField > & Y
fluidMulticomponentThermo & thermo
Definition: createFields.H:31