Raoult.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2015-2016 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 "Raoult.H"
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template<class Thermo, class OtherThermo>
32 (
33  const dictionary& dict,
34  const phasePair& pair
35 )
36 :
37  InterfaceCompositionModel<Thermo, OtherThermo>(dict, pair),
38  YNonVapour_
39  (
40  IOobject
41  (
42  IOobject::groupName("YNonVapour", pair.name()),
43  pair.phase1().mesh().time().timeName(),
44  pair.phase1().mesh()
45  ),
46  pair.phase1().mesh(),
47  dimensionedScalar("one", dimless, 1)
48  ),
49  YNonVapourPrime_
50  (
51  IOobject
52  (
53  IOobject::groupName("YNonVapourPrime", pair.name()),
54  pair.phase1().mesh().time().timeName(),
55  pair.phase1().mesh()
56  ),
57  pair.phase1().mesh(),
59  )
60 {
61  forAllConstIter(hashedWordList, this->speciesNames_, iter)
62  {
63  speciesModels_.insert
64  (
65  *iter,
66  autoPtr<interfaceCompositionModel>
67  (
69  (
70  dict.subDict(*iter),
71  pair
72  )
73  )
74  );
75  }
76 }
77 
78 
79 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
80 
81 template<class Thermo, class OtherThermo>
83 {}
84 
85 
86 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
87 
88 template<class Thermo, class OtherThermo>
90 (
91  const volScalarField& Tf
92 )
93 {
94  YNonVapour_ = scalar(1);
95 
97  (
98  HashTable<autoPtr<interfaceCompositionModel>>,
99  speciesModels_,
100  iter
101  )
102  {
103  iter()->update(Tf);
104 
105  YNonVapour_ -=
106  this->otherThermo_.composition().Y(iter.key())
107  *iter()->Yf(iter.key(), Tf);
108 
109  YNonVapourPrime_ -=
110  this->otherThermo_.composition().Y(iter.key())
111  *iter()->YfPrime(iter.key(), Tf);
112  }
113 }
114 
115 
116 template<class Thermo, class OtherThermo>
119 (
120  const word& speciesName,
121  const volScalarField& Tf
122 ) const
123 {
124  if (this->speciesNames_.contains(speciesName))
125  {
126  return
127  this->otherThermo_.composition().Y(speciesName)
128  *speciesModels_[speciesName]->Yf(speciesName, Tf);
129  }
130  else
131  {
132  return
133  this->thermo_.composition().Y(speciesName)
134  *YNonVapour_;
135  }
136 }
137 
138 
139 template<class Thermo, class OtherThermo>
142 (
143  const word& speciesName,
144  const volScalarField& Tf
145 ) const
146 {
147  if (this->speciesNames_.contains(speciesName))
148  {
149  return
150  this->otherThermo_.composition().Y(speciesName)
151  *speciesModels_[speciesName]->YfPrime(speciesName, Tf);
152  }
153  else
154  {
155  return
156  this->otherThermo_.composition().Y(speciesName)
157  *YNonVapourPrime_;
158  }
159 }
160 
161 
162 // ************************************************************************* //
const OtherThermo & otherThermo_
Other Thermo.
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:453
static autoPtr< interfaceCompositionModel > New(const dictionary &dict, const phasePair &pair)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
dynamicFvMesh & mesh
virtual void update(const volScalarField &Tf)
Update the composition.
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
Definition: dimensionSets.H:52
word timeName
Definition: getTimeIndex.H:3
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
phaseModel & phase1
virtual tmp< volScalarField > Yf(const word &speciesName, const volScalarField &Tf) const
The interface species fraction.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const hashedWordList speciesNames_
Names of the transferring species.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
Raoult(const dictionary &dict, const phasePair &pair)
Construct from components.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
bool contains(const word &) const
Does the list contain the specified name.
A class for managing temporary objects.
Definition: PtrList.H:53
virtual tmp< volScalarField > YfPrime(const word &speciesName, const volScalarField &Tf) const
The interface species fraction derivative w.r.t. temperature.