NonRandomTwoLiquid.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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2015-2018 OpenFOAM Foundation
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 \*---------------------------------------------------------------------------*/
27 
28 #include "NonRandomTwoLiquid.H"
29 
30 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
31 
32 template<class Thermo, class OtherThermo>
35 (
36  const dictionary& dict,
37  const phasePair& pair
38 )
39 :
41  gamma1_
42  (
43  IOobject
44  (
45  IOobject::groupName("gamma1", pair.name()),
46  pair.phase1().mesh().time().timeName(),
47  pair.phase1().mesh()
48  ),
49  pair.phase1().mesh(),
50  dimensionedScalar("one", dimless, 1)
51  ),
52  gamma2_
53  (
54  IOobject
55  (
56  IOobject::groupName("gamma2", pair.name()),
57  pair.phase1().mesh().time().timeName(),
58  pair.phase1().mesh()
59  ),
60  pair.phase1().mesh(),
61  dimensionedScalar("one", dimless, 1)
62  ),
63  alpha12_("alpha12", dimless, Zero),
64  alpha21_("alpha21", dimless, Zero),
65  beta12_("beta12", dimless/dimTemperature, Zero),
66  beta21_("beta21", dimless/dimTemperature, Zero)
67 {
68  if (this->speciesNames_.size() != 2)
69  {
71  << "NonRandomTwoLiquid model is suitable for two species only."
72  << exit(FatalError);
73  }
74 
75  species1Name_ = this->speciesNames_[0];
76  species2Name_ = this->speciesNames_[1];
77 
78  species1Index_ = this->thermo_.composition().species()[species1Name_];
79  species2Index_ = this->thermo_.composition().species()[species2Name_];
80 
81  alpha12_.read("alpha", dict.subDict(species1Name_));
82  alpha21_.read("alpha", dict.subDict(species2Name_));
83  beta12_.read("beta", dict.subDict(species1Name_));
84  beta21_.read("beta", dict.subDict(species2Name_));
85 
86  saturationModel12_.reset
87  (
89  (
90  dict.subDict(species1Name_).subDict("interaction"),
91  pair.phase1().mesh()
92  ).ptr()
93  );
94  saturationModel21_.reset
95  (
97  (
98  dict.subDict(species2Name_).subDict("interaction"),
99  pair.phase1().mesh()
100  ).ptr()
101  );
102 
103  speciesModel1_.reset
104  (
106  (
107  dict.subDict(species1Name_),
108  pair
109  ).ptr()
110  );
111 
112  speciesModel2_.reset
113  (
115  (
116  dict.subDict(species2Name_),
117  pair
118  ).ptr()
119  );
120 }
121 
122 
123 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
124 
125 template<class Thermo, class OtherThermo>
128 {}
129 
130 
131 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
132 
133 template<class Thermo, class OtherThermo>
134 void
136 update
137 (
138  const volScalarField& Tf
139 )
140 {
141  volScalarField W(this->thermo_.W());
142 
143  volScalarField X1
144  (
145  this->thermo_.composition().Y(species1Index_)
146  *W
148  (
149  "W",
151  this->thermo_.composition().W(species1Index_)
152  )
153  );
154 
155  volScalarField X2
156  (
157  this->thermo_.composition().Y(species2Index_)
158  *W
160  (
161  "W",
163  this->thermo_.composition().W(species2Index_)
164  )
165  );
166 
167  volScalarField alpha12(alpha12_ + Tf*beta12_);
168  volScalarField alpha21(alpha21_ + Tf*beta21_);
169 
170  volScalarField tau12(saturationModel12_->lnPSat(Tf));
171  volScalarField tau21(saturationModel21_->lnPSat(Tf));
172 
173  volScalarField G12(exp(- alpha12*tau12));
174  volScalarField G21(exp(- alpha21*tau21));
175 
176  gamma1_ =
177  exp
178  (
179  sqr(X2)
180  *(
181  tau21*sqr(G21)/max(sqr(X1 + X2*G21), SMALL)
182  + tau12*G12/max(sqr(X2 + X1*G12), SMALL)
183  )
184  );
185  gamma2_ =
186  exp
187  (
188  sqr(X1)
189  *(
190  tau12*sqr(G12)/max(sqr(X2 + X1*G12), SMALL)
191  + tau21*G21/max(sqr(X1 + X2*G21), SMALL)
192  )
193  );
194 }
195 
196 
197 template<class Thermo, class OtherThermo>
200 (
201  const word& speciesName,
202  const volScalarField& Tf
203 ) const
204 {
205  if (speciesName == species1Name_)
206  {
207  return
208  this->otherThermo_.composition().Y(speciesName)
209  *speciesModel1_->Yf(speciesName, Tf)
210  *gamma1_;
211  }
212  else if (speciesName == species2Name_)
213  {
214  return
215  this->otherThermo_.composition().Y(speciesName)
216  *speciesModel2_->Yf(speciesName, Tf)
217  *gamma2_;
218  }
219  else
220  {
221  return
222  this->thermo_.composition().Y(speciesName)
223  *(scalar(1) - Yf(species1Name_, Tf) - Yf(species2Name_, Tf));
224  }
225 }
226 
227 
228 template<class Thermo, class OtherThermo>
231 YfPrime
232 (
233  const word& speciesName,
234  const volScalarField& Tf
235 ) const
236 {
237  if (speciesName == species1Name_)
238  {
239  return
240  this->otherThermo_.composition().Y(speciesName)
241  *speciesModel1_->YfPrime(speciesName, Tf)
242  *gamma1_;
243  }
244  else if (speciesName == species2Name_)
245  {
246  return
247  this->otherThermo_.composition().Y(speciesName)
248  *speciesModel2_->YfPrime(speciesName, Tf)
249  *gamma2_;
250  }
251  else
252  {
253  return
254  - this->thermo_.composition().Y(speciesName)
255  *(YfPrime(species1Name_, Tf) + YfPrime(species2Name_, Tf));
256  }
257 }
258 
259 
260 // ************************************************************************* //
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Foam::phasePair
Description for mass transfer between a pair of phases. The direction of the mass transfer is from th...
Definition: phasePair.H:53
Foam::InterfaceCompositionModel
Base class for interface composition models, templated on the two thermodynamic models either side of...
Definition: InterfaceCompositionModel.H:58
Foam::interfaceCompositionModels::NonRandomTwoLiquid::NonRandomTwoLiquid
NonRandomTwoLiquid(const dictionary &dict, const phasePair &pair)
Construct from components.
Definition: NonRandomTwoLiquid.C:35
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::interfaceCompositionModels::NonRandomTwoLiquid::Yf
virtual tmp< volScalarField > Yf(const word &speciesName, const volScalarField &Tf) const
The interface species fraction.
Definition: NonRandomTwoLiquid.C:200
Foam::dimMoles
const dimensionSet dimMoles(0, 0, 0, 0, 1, 0, 0)
Definition: dimensionSets.H:55
Foam::exp
dimensionedScalar exp(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:261
Foam::dimensionedScalar
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Definition: dimensionedScalarFwd.H:42
Foam::interfaceCompositionModels::NonRandomTwoLiquid::update
virtual void update(const volScalarField &Tf)
Update the composition.
Definition: NonRandomTwoLiquid.C:137
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::FatalError
error FatalError
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
Foam::dimMass
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:51
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Foam::New
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
Definition: DimensionedFieldReuseFunctions.H:105
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::interfaceCompositionModels::NonRandomTwoLiquid::YfPrime
virtual tmp< volScalarField > YfPrime(const word &speciesName, const volScalarField &Tf) const
The interface species fraction derivative w.r.t. temperature.
Definition: NonRandomTwoLiquid.C:232
Foam::sqr
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Definition: dimensionedSymmTensor.C:51
Foam::phasePair::name
virtual word name() const
Pair name.
Definition: phasePair.C:69
Foam::dimTemperature
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
Definition: dimensionSets.H:54
Foam::phasePair::phase1
const phaseModel & phase1() const
Definition: phasePairI.H:30
NonRandomTwoLiquid.H
Foam::GeometricField< scalar, fvPatchField, volMesh >
Foam::dimless
const dimensionSet dimless
Dimensionless.
Definition: dimensionSets.C:189
Foam::interfaceCompositionModels::NonRandomTwoLiquid::~NonRandomTwoLiquid
virtual ~NonRandomTwoLiquid()
Destructor.
Definition: NonRandomTwoLiquid.C:127