moleFractions.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-2021 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 "moleFractions.H"
27 #include "basicThermo.H"
28 
29 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
30 
31 template<class ThermoType>
33 {
34  const ThermoType& thermo =
35  mesh_.lookupObject<ThermoType>
36  (
37  IOobject::groupName(physicalProperties::typeName, phaseName_)
38  );
39 
40  const PtrList<volScalarField>& Y = thermo.composition().Y();
41 
42  const volScalarField W(thermo.W());
43 
44  forAll(Y, i)
45  {
46  const dimensionedScalar Wi
47  (
48  "Wi",
50  thermo.composition().Wi(i)
51  );
52 
53  X_[i] = W*Y[i]/Wi;
54  }
55 }
56 
57 
58 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
59 
60 template<class ThermoType>
62 (
63  const word& name,
64  const Time& runTime,
65  const dictionary& dict
66 )
67 :
68  fvMeshFunctionObject(name, runTime, dict),
69  phaseName_(dict.lookupOrDefault<word>("phase", word::null))
70 {
71  const word dictName
72  (
73  IOobject::groupName(physicalProperties::typeName, phaseName_)
74  );
75 
76  if (mesh_.foundObject<ThermoType>(dictName))
77  {
78  const ThermoType& thermo = mesh_.lookupObject<ThermoType>(dictName);
79 
80  const PtrList<volScalarField>& Y = thermo.composition().Y();
81 
82  X_.setSize(Y.size());
83 
84  forAll(Y, i)
85  {
86  X_.set
87  (
88  i,
89  new volScalarField
90  (
91  IOobject
92  (
93  "X_" + Y[i].name(),
94  mesh_.time().timeName(),
95  mesh_
96  ),
97  mesh_,
99  )
100  );
101  }
102 
103  calculateMoleFractions();
104  }
105  else
106  {
107  if (phaseName_ != word::null)
108  {
110  << "Cannot find thermodynamics model of type "
111  << ThermoType::typeName
112  << " for phase "
113  << phaseName_
114  << exit(FatalError);
115  }
116  else
117  {
119  << "Cannot find thermodynamics model of type "
120  << ThermoType::typeName
121  << exit(FatalError);
122  }
123  }
124 }
125 
126 
127 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
128 
129 template<class ThermoType>
131 {}
132 
133 
134 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
135 
136 template<class ThermoType>
138 (
139  const dictionary& dict
140 )
141 {
142  return true;
143 }
144 
145 
146 template<class ThermoType>
148 {
149  calculateMoleFractions();
150  return true;
151 }
152 
153 
154 template<class ThermoType>
156 {
157  forAll(X_, i)
158  {
159  X_[i].write();
160  }
161 
162  return true;
163 }
164 
165 
166 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
virtual ~moleFractions()
Destructor.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
const dimensionSet dimless
const ObjectType & lookupObject(const word &fieldName) const
Lookup object from the objectRegistry.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:69
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:58
moleFractions(const word &name, const Time &t, const dictionary &dict)
Construct from Time and dictionary.
Definition: moleFractions.C:62
A class for handling words, derived from string.
Definition: word.H:59
virtual bool write()
The mole-fraction fields auto-write.
This function object calculates mole-fraction fields from the mass-fraction fields of the psi/rhoReac...
Definition: moleFractions.H:76
void setSize(const label)
Reset size of PtrList. If extending the PtrList, new entries are.
Definition: PtrList.C:131
const dimensionSet dimMass
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
virtual bool execute()
Calculate the mole-fraction fields.
const dimensionSet dimMoles
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
const scalarList W(::W(thermo))
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:98
virtual bool read(const dictionary &)
Read the moleFractions data.
const word dictName("noiseDict")