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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 2019-2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
31#include "surfaceFields.H"
32#include "fvc.H"
33
34// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
35
36void Foam::incompressibleThreePhaseMixture::calcNu()
37{
38 nuModel1_->correct();
39 nuModel2_->correct();
40 nuModel3_->correct();
41
42 // Average kinematic viscosity calculated from dynamic viscosity
43 nu_ = mu()/(alpha1_*rho1_ + alpha2_*rho2_ + alpha3_*rho3_);
44}
45
46
47// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
48
50(
51 const volVectorField& U,
53)
54:
55 IOdictionary
56 (
57 IOobject
58 (
59 "transportProperties",
60 U.time().constant(),
61 U.db(),
62 IOobject::MUST_READ_IF_MODIFIED,
63 IOobject::NO_WRITE
64 )
65 ),
66
67 phase1Name_(get<wordList>("phases")[0]),
68 phase2Name_(get<wordList>("phases")[1]),
69 phase3Name_(get<wordList>("phases")[2]),
70
71 alpha1_
72 (
73 IOobject
74 (
75 IOobject::groupName("alpha", phase1Name_),
76 U.time().timeName(),
77 U.mesh(),
78 IOobject::MUST_READ,
79 IOobject::AUTO_WRITE
80 ),
81 U.mesh()
82 ),
83
84 alpha2_
85 (
86 IOobject
87 (
88 IOobject::groupName("alpha", phase2Name_),
89 U.time().timeName(),
90 U.mesh(),
91 IOobject::MUST_READ,
92 IOobject::AUTO_WRITE
93 ),
94 U.mesh()
95 ),
96
97 alpha3_
98 (
99 IOobject
100 (
101 IOobject::groupName("alpha", phase3Name_),
102 U.time().timeName(),
103 U.mesh(),
104 IOobject::MUST_READ,
105 IOobject::AUTO_WRITE
106 ),
107 U.mesh()
108 ),
109
110 U_(U),
111 phi_(phi),
112
113 nu_
114 (
115 IOobject
116 (
117 "nu",
118 U.time().timeName(),
119 U.db()
120 ),
121 U.mesh(),
122 dimensionedScalar(dimensionSet(0, 2, -1, 0, 0), Zero),
123 calculatedFvPatchScalarField::typeName
124 ),
125
126 nuModel1_
127 (
128 viscosityModel::New
129 (
130 "nu1",
131 subDict(phase1Name_),
132 U,
133 phi
134 )
135 ),
136 nuModel2_
137 (
138 viscosityModel::New
139 (
140 "nu2",
141 subDict(phase2Name_),
142 U,
143 phi
144 )
145 ),
146 nuModel3_
147 (
148 viscosityModel::New
149 (
150 "nu3",
151 subDict(phase3Name_),
152 U,
153 phi
154 )
155 ),
156
157 rho1_("rho", dimDensity, nuModel1_->viscosityProperties()),
158 rho2_("rho", dimDensity, nuModel2_->viscosityProperties()),
159 rho3_("rho", dimDensity, nuModel3_->viscosityProperties())
160{
161 alpha3_ == 1.0 - alpha1_ - alpha2_;
162 calcNu();
163}
164
165
166// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
167
170{
171 return tmp<volScalarField>
172 (
174 (
175 "mu",
176 alpha1_*rho1_*nuModel1_->nu()
177 + alpha2_*rho2_*nuModel2_->nu()
178 + alpha3_*rho3_*nuModel3_->nu()
179 )
180 );
181}
182
183
186{
187 surfaceScalarField alpha1f(fvc::interpolate(alpha1_));
188 surfaceScalarField alpha2f(fvc::interpolate(alpha2_));
189 surfaceScalarField alpha3f(fvc::interpolate(alpha3_));
190
191 return tmp<surfaceScalarField>
192 (
194 (
195 "mu",
196 alpha1f*rho1_*fvc::interpolate(nuModel1_->nu())
197 + alpha2f*rho2_*fvc::interpolate(nuModel2_->nu())
198 + alpha3f*rho3_*fvc::interpolate(nuModel3_->nu())
199 )
200 );
201}
202
203
206{
207 surfaceScalarField alpha1f(fvc::interpolate(alpha1_));
208 surfaceScalarField alpha2f(fvc::interpolate(alpha2_));
209 surfaceScalarField alpha3f(fvc::interpolate(alpha3_));
210
211 return tmp<surfaceScalarField>
212 (
214 (
215 "nu",
216 (
217 alpha1f*rho1_*fvc::interpolate(nuModel1_->nu())
218 + alpha2f*rho2_*fvc::interpolate(nuModel2_->nu())
219 + alpha3f*rho3_*fvc::interpolate(nuModel3_->nu())
220 )/(alpha1f*rho1_ + alpha2f*rho2_ + alpha3f*rho3_)
221 )
222 );
223}
224
225
227{
229 {
230 if
231 (
232 nuModel1_().read(*this)
233 && nuModel2_().read(*this)
234 && nuModel3_().read(*this)
235 )
236 {
237 nuModel1_->viscosityProperties().readEntry("rho", rho1_);
238 nuModel2_->viscosityProperties().readEntry("rho", rho2_);
239 nuModel3_->viscosityProperties().readEntry("rho", rho3_);
240
241 return true;
242 }
243 }
244
245 return false;
246}
247
248
249// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
surfaceScalarField & phi
tmp< surfaceScalarField > muf() const
Return the face-interpolated dynamic laminar viscosity.
tmp< volScalarField > mu() const
Return the dynamic laminar viscosity.
bool read()
Read base transportProperties dictionary.
tmp< surfaceScalarField > nuf() const
Return the face-interpolated dynamic laminar viscosity.
constant condensation/saturation model.
A class for managing temporary objects.
Definition: tmp.H:65
virtual bool read()=0
Read transportProperties dictionary.
U
Definition: pEqn.H:72
dynamicFvMesh & mesh
word timeName
Definition: getTimeIndex.H:3
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
List< word > wordList
A List of words.
Definition: fileName.H:63
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:83
bool read(const char *buf, int32_t &val)
Same as readInt32.
Definition: int32.H:108
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
const dimensionSet dimDensity
Foam::surfaceFields.