incompressibleThreePhaseMixture.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) 2011-2014 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 
28 #include "surfaceFields.H"
29 #include "fvc.H"
30 
31 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
32 
33 //- Calculate and return the laminar viscosity
34 void Foam::incompressibleThreePhaseMixture::calcNu()
35 {
36  nuModel1_->correct();
37  nuModel2_->correct();
38  nuModel3_->correct();
39 
40  // Average kinematic viscosity calculated from dynamic viscosity
41  nu_ = mu()/(alpha1_*rho1_ + alpha2_*rho2_ + alpha3_*rho3_);
42 }
43 
44 
45 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
46 
48 (
49  const volVectorField& U,
50  const surfaceScalarField& phi
51 )
52 :
53  IOdictionary
54  (
55  IOobject
56  (
57  "transportProperties",
58  U.time().constant(),
59  U.db(),
60  IOobject::MUST_READ_IF_MODIFIED,
61  IOobject::NO_WRITE
62  )
63  ),
64 
65  phase1Name_(wordList(lookup("phases"))[0]),
66  phase2Name_(wordList(lookup("phases"))[1]),
67  phase3Name_(wordList(lookup("phases"))[2]),
68 
69  alpha1_
70  (
71  IOobject
72  (
73  IOobject::groupName("alpha", phase1Name_),
74  U.time().timeName(),
75  U.mesh(),
76  IOobject::MUST_READ,
77  IOobject::AUTO_WRITE
78  ),
79  U.mesh()
80  ),
81 
82  alpha2_
83  (
84  IOobject
85  (
86  IOobject::groupName("alpha", phase2Name_),
87  U.time().timeName(),
88  U.mesh(),
89  IOobject::MUST_READ,
90  IOobject::AUTO_WRITE
91  ),
92  U.mesh()
93  ),
94 
95  alpha3_
96  (
97  IOobject
98  (
99  IOobject::groupName("alpha", phase3Name_),
100  U.time().timeName(),
101  U.mesh(),
102  IOobject::MUST_READ,
103  IOobject::AUTO_WRITE
104  ),
105  U.mesh()
106  ),
107 
108  U_(U),
109  phi_(phi),
110 
111  nu_
112  (
113  IOobject
114  (
115  "nu",
116  U.time().timeName(),
117  U.db()
118  ),
119  U.mesh(),
120  dimensionedScalar("nu", dimensionSet(0, 2, -1, 0, 0), 0),
121  calculatedFvPatchScalarField::typeName
122  ),
123 
124  nuModel1_
125  (
126  viscosityModel::New
127  (
128  "nu1",
129  subDict(phase1Name_),
130  U,
131  phi
132  )
133  ),
134  nuModel2_
135  (
136  viscosityModel::New
137  (
138  "nu2",
139  subDict(phase2Name_),
140  U,
141  phi
142  )
143  ),
144  nuModel3_
145  (
146  viscosityModel::New
147  (
148  "nu3",
149  subDict(phase3Name_),
150  U,
151  phi
152  )
153  ),
154 
155  rho1_(nuModel1_->viscosityProperties().lookup("rho")),
156  rho2_(nuModel2_->viscosityProperties().lookup("rho")),
157  rho3_(nuModel3_->viscosityProperties().lookup("rho"))
158 {
159  alpha3_ == 1.0 - alpha1_ - alpha2_;
160  calcNu();
161 }
162 
163 
164 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
165 
168 {
169  return tmp<volScalarField>
170  (
171  new volScalarField
172  (
173  "mu",
174  alpha1_*rho1_*nuModel1_->nu()
175  + alpha2_*rho2_*nuModel2_->nu()
176  + alpha3_*rho3_*nuModel3_->nu()
177  )
178  );
179 }
180 
181 
184 {
185  surfaceScalarField alpha1f(fvc::interpolate(alpha1_));
187  surfaceScalarField alpha3f(fvc::interpolate(alpha3_));
188 
189  return tmp<surfaceScalarField>
190  (
192  (
193  "mu",
194  alpha1f*rho1_*fvc::interpolate(nuModel1_->nu())
195  + alpha2f*rho2_*fvc::interpolate(nuModel2_->nu())
196  + alpha3f*rho3_*fvc::interpolate(nuModel3_->nu())
197  )
198  );
199 }
200 
201 
204 {
205  surfaceScalarField alpha1f(fvc::interpolate(alpha1_));
207  surfaceScalarField alpha3f(fvc::interpolate(alpha3_));
208 
209  return tmp<surfaceScalarField>
210  (
212  (
213  "nu",
214  (
215  alpha1f*rho1_*fvc::interpolate(nuModel1_->nu())
216  + alpha2f*rho2_*fvc::interpolate(nuModel2_->nu())
217  + alpha3f*rho3_*fvc::interpolate(nuModel3_->nu())
218  )/(alpha1f*rho1_ + alpha2f*rho2_ + alpha3f*rho3_)
219  )
220  );
221 }
222 
223 
225 {
226  if (transportModel::read())
227  {
228  if
229  (
230  nuModel1_().read(*this)
231  && nuModel2_().read(*this)
232  && nuModel3_().read(*this)
233  )
234  {
235  nuModel1_->viscosityProperties().lookup("rho") >> rho1_;
236  nuModel2_->viscosityProperties().lookup("rho") >> rho2_;
237  nuModel3_->viscosityProperties().lookup("rho") >> rho3_;
238 
239  return true;
240  }
241  else
242  {
243  return false;
244  }
245  }
246  else
247  {
248  return false;
249  }
250 }
251 
252 
253 // ************************************************************************* //
tmp< surfaceScalarField > nuf() const
Return the face-interpolated dynamic laminar viscosity.
dynamicFvMesh & mesh
A class for managing temporary objects.
Definition: PtrList.H:90
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
virtual bool read()=0
Read transportProperties dictionary.
U
Definition: pEqn.H:86
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &vf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
bool read()
Read base transportProperties dictionary.
surfaceScalarField & phi
List< word > wordList
A List of words.
Definition: fileName.H:55
Macros for easy insertion into run-time selection tables.
word timeName
Definition: getTimeIndex.H:3
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
surfaceScalarField alpha2f("alpha2f", scalar(1)-alpha1f)
tmp< volScalarField > mu() const
Return the dynamic laminar viscosity.
tmp< surfaceScalarField > muf() const
Return the face-interpolated dynamic laminar viscosity.
incompressibleThreePhaseMixture(const volVectorField &U, const surfaceScalarField &phi)
Construct from components.
Foam::surfaceFields.