incompressibleTwoPhaseVoFMixture.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2024 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 
27 #include "surfaceInterpolate.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
35 }
36 
37 
38 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
39 
41 (
42  const fvMesh& mesh
43 )
44 :
45  twoPhaseVoFMixture(mesh),
46 
47  nuModel1_(viscosityModel::New(mesh, phase1Name())),
48  nuModel2_(viscosityModel::New(mesh, phase2Name())),
49 
50  rho1_("rho", dimDensity, nuModel1_()),
51  rho2_("rho", dimDensity, nuModel2_()),
52 
53  rho_
54  (
55  IOobject
56  (
57  "rho",
58  mesh.time().name(),
59  mesh,
60  IOobject::NO_READ,
61  IOobject::AUTO_WRITE
62  ),
63  mesh,
64  dimensionedScalar("rho", dimDensity, 0)
65  ),
66 
67  nu_
68  (
69  IOobject
70  (
71  "nu",
72  mesh.time().name(),
73  mesh
74  ),
75  mesh,
77  calculatedFvPatchScalarField::typeName
78  )
79 {
80  correct();
81 }
82 
83 
84 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
85 
88 {
89  const volScalarField limitedAlpha1
90  (
91  min(max(alpha1(), scalar(0)), scalar(1))
92  );
93 
94  return volScalarField::New
95  (
96  "mu",
97  limitedAlpha1*rho1_*nuModel1_->nu()
98  + (scalar(1) - limitedAlpha1)*rho2_*nuModel2_->nu()
99  );
100 }
101 
102 
105 {
106  const surfaceScalarField alpha1f
107  (
108  min(max(fvc::interpolate(alpha1()), scalar(0)), scalar(1))
109  );
110 
112  (
113  "muf",
114  alpha1f*rho1_*fvc::interpolate(nuModel1_->nu())
115  + (scalar(1) - alpha1f)*rho2_*fvc::interpolate(nuModel2_->nu())
116  );
117 }
118 
119 
122 {
123  const surfaceScalarField alpha1f
124  (
125  min(max(fvc::interpolate(alpha1()), scalar(0)), scalar(1))
126  );
127 
129  (
130  "nuf",
131  (
132  alpha1f*rho1_*fvc::interpolate(nuModel1_->nu())
133  + (scalar(1) - alpha1f)*rho2_*fvc::interpolate(nuModel2_->nu())
134  )/(alpha1f*rho1_ + (scalar(1) - alpha1f)*rho2_)
135  );
136 }
137 
138 
140 {
142  {
143  nuModel1_->lookup("rho") >> rho1_;
144  nuModel2_->lookup("rho") >> rho2_;
145 
146  return true;
147  }
148  else
149  {
150  return false;
151  }
152 }
153 
154 
156 {
157  rho_ = alpha1()*rho1_ + alpha2()*rho2_;
158 
159  nuModel1_->correct();
160  nuModel2_->correct();
161 
162  const volScalarField limitedAlpha1
163  (
164  "limitedAlpha1",
165  min(max(alpha1(), scalar(0)), scalar(1))
166  );
167 
168  // Average kinematic viscosity calculated from dynamic viscosity
169  nu_ = mu()/(limitedAlpha1*rho1_ + (scalar(1) - limitedAlpha1)*rho2_);
170 }
171 
172 
173 // ************************************************************************* //
Macros for easy insertion into run-time selection tables.
Generic GeometricField class.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &, const HashPtrTable< Source > &=HashPtrTable< Source >())
Return a temporary field constructed from name,.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
Class to represent a mixture of two constant density phases.
tmp< surfaceScalarField > muf() const
Return the face-interpolated dynamic laminar viscosity.
tmp< surfaceScalarField > nuf() const
Return the face-interpolated kinematic laminar viscosity.
virtual void correct()
Correct the mixture density and laminar viscosity.
tmp< volScalarField > mu() const
Return the dynamic laminar viscosity.
incompressibleTwoPhaseVoFMixture(const fvMesh &mesh)
Construct from a mesh.
virtual bool read()
Read base phaseProperties dictionary.
A class for managing temporary objects.
Definition: tmp.H:55
Class to represent a VoF mixture.
virtual bool read()=0
Read base phaseProperties dictionary.
An abstract base class for Newtonian viscosity models.
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
const dimensionedScalar mu
Atomic mass unit.
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
Namespace for OpenFOAM.
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
const dimensionSet dimKinematicViscosity
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
const dimensionSet dimDensity
defineTypeNameAndDebug(combustionModel, 0)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)