incompressibleTwoPhaseInteractingMixture.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) 2014-2018 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 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35  defineTypeNameAndDebug(incompressibleTwoPhaseInteractingMixture, 0);
36 }
37 
38 
39 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
40 
43 (
44  const volVectorField& U,
45  const surfaceScalarField& phi
46 )
47 :
48  IOdictionary
49  (
50  IOobject
51  (
52  "transportProperties",
53  U.time().constant(),
54  U.db(),
55  IOobject::MUST_READ_IF_MODIFIED,
56  IOobject::NO_WRITE
57  )
58  ),
59  twoPhaseMixture(U.mesh(), *this),
60 
61  muModel_
62  (
63  mixtureViscosityModel::New
64  (
65  "mu",
66  subDict(phase1Name_),
67  U,
68  phi
69  )
70  ),
71 
72  nucModel_
73  (
74  viscosityModel::New
75  (
76  "nuc",
77  subDict(phase2Name_),
78  U,
79  phi
80  )
81  ),
82 
83  rhod_("rho", dimDensity, muModel_->viscosityProperties()),
84  rhoc_("rho", dimDensity, nucModel_->viscosityProperties()),
85  dd_
86  (
87  "d",
88  dimLength,
89  muModel_->viscosityProperties().lookupOrDefault("d", 0.0)
90  ),
91  alphaMax_(muModel_->viscosityProperties().lookupOrDefault("alphaMax", 1.0)),
92 
93  U_(U),
94  phi_(phi),
95 
96  mu_
97  (
98  IOobject
99  (
100  "mu",
101  U_.time().timeName(),
102  U_.db()
103  ),
104  U_.mesh(),
105  dimensionedScalar("mu", dimensionSet(1, -1, -1, 0, 0), 0),
106  calculatedFvPatchScalarField::typeName
107  )
108 {
109  correct();
110 }
111 
112 
113 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
114 
116 {
117  if (regIOobject::read())
118  {
119  if
120  (
123  )
124  {
125  muModel_->viscosityProperties().lookup("rho") >> rhod_;
126  nucModel_->viscosityProperties().lookup("rho") >> rhoc_;
127 
129  (
130  "d",
131  dimLength,
132  muModel_->viscosityProperties().lookupOrDefault("d", 0)
133  );
134 
135  alphaMax_ =
136  muModel_->viscosityProperties().lookupOrDefault
137  (
138  "alphaMax",
139  1.0
140  );
141 
142  return true;
143  }
144  else
145  {
146  return false;
147  }
148  }
149  else
150  {
151  return false;
152  }
153 }
154 
155 
156 // ************************************************************************* //
Foam::surfaceFields.
surfaceScalarField & phi
virtual bool read()
Read base transportProperties dictionary.
virtual bool read()
Read object.
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
Macros for easy insertion into run-time selection tables.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:692
dynamicFvMesh & mesh
scalar alphaMax_
Optional maximum dispersed phase-fraction (e.g. packing limit)
word timeName
Definition: getTimeIndex.H:3
incompressibleTwoPhaseInteractingMixture(const volVectorField &U, const surfaceScalarField &phi)
Construct from components.
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(psi *p+alphal *rhol0+((alphav *psiv+alphal *psil) - psi) *pSat, rhoMin);# 1 "/home/ubuntu/OpenFOAM-6/applications/solvers/multiphase/cavitatingFoam/alphavPsi.H" 1{ alphav=max(min((rho - rholSat)/(rhovSat - rholSat), scalar(1)), scalar(0));alphal=1.0 - alphav;Info<< "max-min alphav: "<< max(alphav).value()<< " "<< min(alphav).value()<< endl;psiModel-> correct()
Definition: pEqn.H:72
defineTypeNameAndDebug(combustionModel, 0)
const dimensionSet dimDensity
U
Definition: pEqn.H:72
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedScalar dd_
Optional diameter of the dispersed phase particles.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Namespace for OpenFOAM.