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  (
68  "template<class Thermo, class OtherThermo>"
69  "Foam::interfaceCompositionModels::"
70  "NonRandomTwoLiquid<Thermo, OtherThermo>::"
71  "NonRandomTwoLiquid"
72  "( "
73  "const dictionary& dict, "
74  "const phasePair& pair "
75  ")"
76  ) << "NonRandomTwoLiquid model is suitable for two species only."
77  << exit(FatalError);
78  }
79 
80  species1Name_ = this->speciesNames_[0];
81  species2Name_ = this->speciesNames_[1];
82 
83  species1Index_ = this->thermo_.composition().species()[species1Name_];
84  species2Index_ = this->thermo_.composition().species()[species2Name_];
85 
86  alpha12_ = dimensionedScalar
87  (
88  "alpha12",
89  dimless,
90  dict.subDict(species1Name_).lookup("alpha")
91  );
92  alpha21_ = dimensionedScalar
93  (
94  "alpha21",
95  dimless,
96  dict.subDict(species2Name_).lookup("alpha")
97  );
98 
99  beta12_ = dimensionedScalar
100  (
101  "beta12",
103  dict.subDict(species1Name_).lookup("beta")
104  );
105  beta21_ = dimensionedScalar
106  (
107  "beta21",
109  dict.subDict(species2Name_).lookup("beta")
110  );
111 
112  saturationModel12_.reset
113  (
115  (
116  dict.subDict(species1Name_).subDict("interaction")
117  ).ptr()
118  );
119  saturationModel21_.reset
120  (
122  (
123  dict.subDict(species2Name_).subDict("interaction")
124  ).ptr()
125  );
126 
127  speciesModel1_.reset
128  (
130  (
131  dict.subDict(species1Name_),
132  pair
133  ).ptr()
134  );
135 
136  speciesModel2_.reset
137  (
139  (
140  dict.subDict(species2Name_),
141  pair
142  ).ptr()
143  );
144 }
145 
146 
147 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
148 
149 template<class Thermo, class OtherThermo>
152 {}
153 
154 
155 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
156 
157 template<class Thermo, class OtherThermo>
158 void
160 update
161 (
162  const volScalarField& Tf
163 )
164 {
165  volScalarField W(this->thermo_.composition().W());
166 
167  volScalarField X1
168  (
169  this->thermo_.composition().Y(species1Index_)
170  *W
171  /this->thermo_.composition().W(species1Index_)
172  );
173 
174  volScalarField X2
175  (
176  this->thermo_.composition().Y(species2Index_)
177  *W
178  /this->thermo_.composition().W(species2Index_)
179  );
180 
181  volScalarField alpha12(alpha12_ + Tf*beta12_);
182  volScalarField alpha21(alpha21_ + Tf*beta21_);
183 
184  volScalarField tau12(saturationModel12_->lnPSat(Tf));
185  volScalarField tau21(saturationModel21_->lnPSat(Tf));
186 
187  volScalarField G12(exp(- alpha12*tau12));
188  volScalarField G21(exp(- alpha21*tau21));
189 
190  gamma1_ =
191  exp
192  (
193  sqr(X2)
194  *(
195  tau21*sqr(G21)/max(sqr(X1 + X2*G21), SMALL)
196  + tau12*G12/max(sqr(X2 + X1*G12), SMALL)
197  )
198  );
199  gamma2_ =
200  exp
201  (
202  sqr(X1)
203  *(
204  tau12*sqr(G12)/max(sqr(X2 + X1*G12), SMALL)
205  + tau21*G21/max(sqr(X1 + X2*G21), SMALL)
206  )
207  );
208 }
209 
210 
211 template<class Thermo, class OtherThermo>
214 (
215  const word& speciesName,
216  const volScalarField& Tf
217 ) const
218 {
219  if (speciesName == species1Name_)
220  {
221  return
222  this->otherThermo_.composition().Y(speciesName)
223  *speciesModel1_->Yf(speciesName, Tf)
224  *gamma1_;
225  }
226  else if(speciesName == species2Name_)
227  {
228  return
229  this->otherThermo_.composition().Y(speciesName)
230  *speciesModel2_->Yf(speciesName, Tf)
231  *gamma2_;
232  }
233  else
234  {
235  return
236  this->thermo_.composition().Y(speciesName)
237  *(scalar(1) - Yf(species1Name_, Tf) - Yf(species2Name_, Tf));
238  }
239 }
240 
241 
242 template<class Thermo, class OtherThermo>
245 YfPrime
246 (
247  const word& speciesName,
248  const volScalarField& Tf
249 ) const
250 {
251  if (speciesName == species1Name_)
252  {
253  return
254  this->otherThermo_.composition().Y(speciesName)
255  *speciesModel1_->YfPrime(speciesName, Tf)
256  *gamma1_;
257  }
258  else if(speciesName == species2Name_)
259  {
260  return
261  this->otherThermo_.composition().Y(speciesName)
262  *speciesModel2_->YfPrime(speciesName, Tf)
263  *gamma2_;
264  }
265  else
266  {
267  return
268  - this->thermo_.composition().Y(speciesName)
269  *(YfPrime(species1Name_, Tf) + YfPrime(species2Name_, Tf));
270  }
271 }
272 
273 
274 // ************************************************************************* //
phaseModel & phase1
Definition: createFields.H:12
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
virtual void update(const volScalarField &Tf)
Update the composition.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
dimensionedScalar exp(const dimensionedScalar &ds)
dynamicFvMesh & mesh
virtual tmp< volScalarField > YfPrime(const word &speciesName, const volScalarField &Tf) const
The interface species fraction derivative w.r.t. temperature.
autoPtr< BasicCompressibleTurbulenceModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleTurbulenceModel::transportModel &transport, const word &propertiesName)
virtual tmp< volScalarField > Yf(const word &speciesName, const volScalarField &Tf) const
The interface species fraction.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
NonRandomTwoLiquid(const dictionary &dict, const phasePair &pair)
Construct from components.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
error FatalError
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
Definition: dimensionSets.H:52
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A class for managing temporary objects.
Definition: PtrList.H:118
word timeName
Definition: getTimeIndex.H:3