readInitialConditions.H
Go to the documentation of this file.
1  word constProp(initialConditions.lookup("constantProperty"));
2  if ((constProp != "pressure") && (constProp != "volume"))
3  {
4  FatalError << "in initialConditions, unknown constantProperty type "
5  << constProp << nl << " Valid types are: pressure volume."
6  << abort(FatalError);
7  }
8 
9  word fractionBasis(initialConditions.lookup("fractionBasis"));
10  if ((fractionBasis != "mass") && (fractionBasis != "mole"))
11  {
12  FatalError << "in initialConditions, unknown fractionBasis type " << nl
13  << "Valid types are: mass or mole."
15  }
16 
17  label nSpecie = Y.size();
18  PtrList<gasHThermoPhysics> specieData(Y.size());
20  {
21  specieData.set
22  (
23  i,
25  (
26  dynamic_cast<const reactingMixture<gasHThermoPhysics>&>
27  (thermo).speciesData()[i]
28  )
29  );
30  }
31 
32  scalarList Y0(nSpecie, 0.0);
33  scalarList X0(nSpecie, 0.0);
34 
35  dictionary fractions(initialConditions.subDict("fractions"));
36  if (fractionBasis == "mole")
37  {
38  forAll(Y, i)
39  {
40  const word& name = Y[i].name();
41  if (fractions.found(name))
42  {
43  X0[i] = readScalar(fractions.lookup(name));
44  }
45  }
46 
47  scalar mw = 0.0;
48  const scalar mTot = sum(X0);
49  forAll(Y, i)
50  {
51  X0[i] /= mTot;
52  mw += specieData[i].W()*X0[i];
53  }
54 
55  forAll(Y, i)
56  {
57  Y0[i] = X0[i]*specieData[i].W()/mw;
58  }
59  }
60  else // mass fraction
61  {
62  forAll(Y, i)
63  {
64  const word& name = Y[i].name();
65  if (fractions.found(name))
66  {
67  Y0[i] = readScalar(fractions.lookup(name));
68  }
69  }
70 
71  scalar invW = 0.0;
72  const scalar mTot = sum(Y0);
73  forAll(Y, i)
74  {
75  Y0[i] /= mTot;
76  invW += Y0[i]/specieData[i].W();
77  }
78  const scalar mw = 1.0/invW;
79 
80  forAll(Y, i)
81  {
82  X0[i] = Y0[i]*mw/specieData[i].W();
83  }
84  }
85 
86  scalar h0 = 0.0;
87  forAll(Y, i)
88  {
89  Y[i] = Y0[i];
90  h0 += Y0[i]*specieData[i].Hs(p[0], T0);
91  }
92 
93  thermo.he() = dimensionedScalar("h", dimEnergy/dimMass, h0);
94  thermo.correct();
95 
96  rho = thermo.rho();
97  scalar rho0 = rho[0];
98  scalar u0 = h0 - p0/rho0;
99  scalar R0 = p0/(rho0*T0);
100  Rspecific[0] = R0;
101 
102  scalar integratedHeat = 0.0;
103 
104  Info << constProp << " will be held constant." << nl
105  << " p = " << p[0] << " [Pa]" << nl
106  << " T = " << thermo.T()[0] << " [K] " << nl
107  << " rho = " << rho[0] << " [kg/m3]" << nl
108  << endl;
#define readScalar
Definition: doubleScalar.C:38
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 scalar mTot
error FatalError
forAll(specieData, i)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
scalar rho0
const scalar mw
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
scalarList Y0(nSpecie, 0.0)
scalarList X0(nSpecie, 0.0)
psiReactionThermo & thermo
Definition: createFields.H:31
word constProp(initialConditions.lookup("constantProperty"))
List< scalar > scalarList
A List of scalars.
Definition: scalarList.H:50
errorManip< error > abort(error &err)
Definition: errorManip.H:131
word fractionBasis(initialConditions.lookup("fractionBasis"))
static const char nl
Definition: Ostream.H:262
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const dimensionSet dimEnergy
PtrList< gasHThermoPhysics > specieData(Y.size())
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
PtrList< volScalarField > & Y
scalar R0
messageStream Info
dictionary fractions(initialConditions.subDict("fractions"))
Rspecific[0]
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
volScalarField & p
sutherlandTransport< species::thermo< janafThermo< perfectGas< specie > >, sensibleEnthalpy > > gasHThermoPhysics