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 | Copyright (C) 2015 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 "NonRandomTwoLiquid.H"
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template<class Thermo, class OtherThermo>
33 (
34  const dictionary& dict,
35  const phasePair& pair
36 )
37 :
38  InterfaceCompositionModel<Thermo, OtherThermo>(dict, pair),
39  gamma1_
40  (
41  IOobject
42  (
43  IOobject::groupName("gamma1", pair.name()),
44  pair.phase1().mesh().time().timeName(),
45  pair.phase1().mesh()
46  ),
47  pair.phase1().mesh(),
48  dimensionedScalar("one", dimless, 1)
49  ),
50  gamma2_
51  (
52  IOobject
53  (
54  IOobject::groupName("gamma2", pair.name()),
55  pair.phase1().mesh().time().timeName(),
56  pair.phase1().mesh()
57  ),
58  pair.phase1().mesh(),
59  dimensionedScalar("one", dimless, 1)
60  ),
61  beta12_("", dimless/dimTemperature, 0),
62  beta21_("", dimless/dimTemperature, 0)
63 {
64  if (this->speciesNames_.size() != 2)
65  {
67  << "NonRandomTwoLiquid model is suitable for two species only."
68  << exit(FatalError);
69  }
70 
71  species1Name_ = this->speciesNames_[0];
72  species2Name_ = this->speciesNames_[1];
73 
74  species1Index_ = this->thermo_.composition().species()[species1Name_];
75  species2Index_ = this->thermo_.composition().species()[species2Name_];
76 
77  alpha12_ = dimensionedScalar
78  (
79  "alpha12",
80  dimless,
81  dict.subDict(species1Name_).lookup("alpha")
82  );
83  alpha21_ = dimensionedScalar
84  (
85  "alpha21",
86  dimless,
87  dict.subDict(species2Name_).lookup("alpha")
88  );
89 
90  beta12_ = dimensionedScalar
91  (
92  "beta12",
94  dict.subDict(species1Name_).lookup("beta")
95  );
96  beta21_ = dimensionedScalar
97  (
98  "beta21",
100  dict.subDict(species2Name_).lookup("beta")
101  );
102 
103  saturationModel12_.reset
104  (
106  (
107  dict.subDict(species1Name_).subDict("interaction")
108  ).ptr()
109  );
110  saturationModel21_.reset
111  (
113  (
114  dict.subDict(species2Name_).subDict("interaction")
115  ).ptr()
116  );
117 
118  speciesModel1_.reset
119  (
121  (
122  dict.subDict(species1Name_),
123  pair
124  ).ptr()
125  );
126 
127  speciesModel2_.reset
128  (
130  (
131  dict.subDict(species2Name_),
132  pair
133  ).ptr()
134  );
135 }
136 
137 
138 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
139 
140 template<class Thermo, class OtherThermo>
143 {}
144 
145 
146 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
147 
148 template<class Thermo, class OtherThermo>
149 void
151 update
152 (
153  const volScalarField& Tf
154 )
155 {
156  volScalarField W(this->thermo_.composition().W());
157 
158  volScalarField X1
159  (
160  this->thermo_.composition().Y(species1Index_)
161  *W
162  /this->thermo_.composition().W(species1Index_)
163  );
164 
165  volScalarField X2
166  (
167  this->thermo_.composition().Y(species2Index_)
168  *W
169  /this->thermo_.composition().W(species2Index_)
170  );
171 
172  volScalarField alpha12(alpha12_ + Tf*beta12_);
173  volScalarField alpha21(alpha21_ + Tf*beta21_);
174 
175  volScalarField tau12(saturationModel12_->lnPSat(Tf));
176  volScalarField tau21(saturationModel21_->lnPSat(Tf));
177 
178  volScalarField G12(exp(- alpha12*tau12));
179  volScalarField G21(exp(- alpha21*tau21));
180 
181  gamma1_ =
182  exp
183  (
184  sqr(X2)
185  *(
186  tau21*sqr(G21)/max(sqr(X1 + X2*G21), SMALL)
187  + tau12*G12/max(sqr(X2 + X1*G12), SMALL)
188  )
189  );
190  gamma2_ =
191  exp
192  (
193  sqr(X1)
194  *(
195  tau12*sqr(G12)/max(sqr(X2 + X1*G12), SMALL)
196  + tau21*G21/max(sqr(X1 + X2*G21), SMALL)
197  )
198  );
199 }
200 
201 
202 template<class Thermo, class OtherThermo>
205 (
206  const word& speciesName,
207  const volScalarField& Tf
208 ) const
209 {
210  if (speciesName == species1Name_)
211  {
212  return
213  this->otherThermo_.composition().Y(speciesName)
214  *speciesModel1_->Yf(speciesName, Tf)
215  *gamma1_;
216  }
217  else if(speciesName == species2Name_)
218  {
219  return
220  this->otherThermo_.composition().Y(speciesName)
221  *speciesModel2_->Yf(speciesName, Tf)
222  *gamma2_;
223  }
224  else
225  {
226  return
227  this->thermo_.composition().Y(speciesName)
228  *(scalar(1) - Yf(species1Name_, Tf) - Yf(species2Name_, Tf));
229  }
230 }
231 
232 
233 template<class Thermo, class OtherThermo>
236 YfPrime
237 (
238  const word& speciesName,
239  const volScalarField& Tf
240 ) const
241 {
242  if (speciesName == species1Name_)
243  {
244  return
245  this->otherThermo_.composition().Y(speciesName)
246  *speciesModel1_->YfPrime(speciesName, Tf)
247  *gamma1_;
248  }
249  else if(speciesName == species2Name_)
250  {
251  return
252  this->otherThermo_.composition().Y(speciesName)
253  *speciesModel2_->YfPrime(speciesName, Tf)
254  *gamma2_;
255  }
256  else
257  {
258  return
259  - this->thermo_.composition().Y(speciesName)
260  *(YfPrime(species1Name_, Tf) + YfPrime(species2Name_, Tf));
261  }
262 }
263 
264 
265 // ************************************************************************* //
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
dimensionedSymmTensor sqr(const dimensionedVector &dv)
static autoPtr< saturationModel > New(const dictionary &dict)
Select null constructed.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
static autoPtr< interfaceCompositionModel > New(const dictionary &dict, const phasePair &pair)
virtual tmp< volScalarField > Yf(const word &speciesName, const volScalarField &Tf) const
The interface species fraction.
virtual tmp< volScalarField > YfPrime(const word &speciesName, const volScalarField &Tf) const
The interface species fraction derivative w.r.t. temperature.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
dynamicFvMesh & mesh
dimensionedScalar exp(const dimensionedScalar &ds)
virtual tmp< volScalarField > Yf(const word &speciesName, const volScalarField &Tf) const =0
Interface mass fraction.
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
Definition: dimensionSets.H:52
NonRandomTwoLiquid(const dictionary &dict, const phasePair &pair)
Construct from components.
word timeName
Definition: getTimeIndex.H:3
phaseModel & phase1
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
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual tmp< volScalarField > YfPrime(const word &speciesName, const volScalarField &Tf) const =0
The interface mass fraction derivative w.r.t. temperature.
A class for managing temporary objects.
Definition: PtrList.H:53
virtual void update(const volScalarField &Tf)
Update the composition.