SchnerrSauer.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-2021 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 "SchnerrSauer.H"
27 #include "mathematicalConstants.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace twoPhaseChangeModels
35 {
36  defineTypeNameAndDebug(SchnerrSauer, 0);
38  (
39  cavitationModel,
40  SchnerrSauer,
41  dictionary
42  );
43 }
44 }
45 
46 
47 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
48 
50 (
52 )
53 :
54  cavitationModel(typeName, mixture),
55 
56  n_("n", dimless/dimVolume, twoPhaseChangeModelCoeffs_),
57  dNuc_("dNuc", dimLength, twoPhaseChangeModelCoeffs_),
58  Cc_("Cc", dimless, twoPhaseChangeModelCoeffs_),
59  Cv_("Cv", dimless, twoPhaseChangeModelCoeffs_),
60 
61  p0_("0", pSat().dimensions(), 0.0)
62 {
63  correct();
64 }
65 
66 
67 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
68 
70 Foam::twoPhaseChangeModels::SchnerrSauer::rRb
71 (
72  const volScalarField& limitedAlpha1
73 ) const
74 {
75  return pow
76  (
78  *limitedAlpha1/(1.0 + alphaNuc() - limitedAlpha1),
79  1.0/3.0
80  );
81 }
82 
83 
85 Foam::twoPhaseChangeModels::SchnerrSauer::alphaNuc() const
86 {
88  return Vnuc/(1 + Vnuc);
89 }
90 
91 
93 Foam::twoPhaseChangeModels::SchnerrSauer::pCoeff
94 (
95  const volScalarField& p
96 ) const
97 {
98  const volScalarField limitedAlpha1
99  (
100  min(max(mixture_.alpha1(), scalar(0)), scalar(1))
101  );
102 
103  const volScalarField rho
104  (
105  limitedAlpha1*mixture_.rho1()
106  + (scalar(1) - limitedAlpha1)*mixture_.rho2()
107  );
108 
109  return
110  (3*mixture_.rho1()*mixture_.rho2())*sqrt(2/(3*mixture_.rho1()))
111  *rRb(limitedAlpha1)/(rho*sqrt(mag(p - pSat()) + 0.01*pSat()));
112 }
113 
114 
117 {
118  const volScalarField& p =
119  mixture_.U().db().lookupObject<volScalarField>("p");
120 
121  const volScalarField pCoeff(this->pCoeff(p));
122 
123  const volScalarField limitedAlpha1
124  (
125  min(max(mixture_.alpha1(), scalar(0)), scalar(1))
126  );
127 
129  (
130  Cc_*limitedAlpha1*pCoeff*max(p - pSat(), p0_),
131 
132  Cv_*(1.0 + alphaNuc() - limitedAlpha1)*pCoeff*min(p - pSat(), p0_)
133  );
134 }
135 
136 
139 {
140  const volScalarField& p =
141  mixture_.U().db().lookupObject<volScalarField>("p");
142 
143  const volScalarField pCoeff(this->pCoeff(p));
144 
145  const volScalarField limitedAlpha1
146  (
147  min(max(mixture_.alpha1(), scalar(0)), scalar(1))
148  );
149 
150  const volScalarField apCoeff(limitedAlpha1*pCoeff);
151 
153  (
154  Cc_*(1.0 - limitedAlpha1)*pos0(p - pSat())*apCoeff,
155 
156  (-Cv_)*(1.0 + alphaNuc() - limitedAlpha1)*neg(p - pSat())*apCoeff
157  );
158 }
159 
160 
162 {
164 }
165 
166 
168 {
169  if (cavitationModel::read())
170  {
171  twoPhaseChangeModelCoeffs_ = optionalSubDict(type() + "Coeffs");
172 
173  twoPhaseChangeModelCoeffs_.lookup("n") >> n_;
174  twoPhaseChangeModelCoeffs_.lookup("dNuc") >> dNuc_;
175  twoPhaseChangeModelCoeffs_.lookup("Cc") >> Cc_;
176  twoPhaseChangeModelCoeffs_.lookup("Cv") >> Cv_;
177 
178  return true;
179  }
180  else
181  {
182  return false;
183  }
184 }
185 
186 
187 // ************************************************************************* //
defineTypeNameAndDebug(cavitationModel, 0)
virtual bool read()
Read the transportProperties dictionary and update.
Definition: SchnerrSauer.C:167
SchnerrSauer(const immiscibleIncompressibleTwoPhaseMixture &mixture)
Construct for mixture.
Definition: SchnerrSauer.C:50
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
Abstract base class for cavitation models.
dimensionedScalar sqrt(const dimensionedScalar &ds)
const dimensionSet dimless
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
dimensionedScalar neg(const dimensionedScalar &ds)
Macros for easy insertion into run-time selection tables.
const dimensionSet dimLength
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: contiguous.H:49
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(rho0+psi *p, rhoMin);# 1 "/home/ubuntu/OpenFOAM-9/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:68
virtual void correct()=0
Correct the phaseChange model.
dimensionedScalar pos0(const dimensionedScalar &ds)
addToRunTimeSelectionTable(cavitationModel, Kunz, dictionary)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar pow3(const dimensionedScalar &ds)
An immiscible incompressible two-phase mixture transport model.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
const dimensionSet dimVolume
virtual Pair< tmp< volScalarField > > mDotAlphal() const
Return the mass condensation and vaporisation rates as a.
Definition: SchnerrSauer.C:116
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual void correct()
Correct the SchnerrSauer phaseChange model.
Definition: SchnerrSauer.C:161
A class for managing temporary objects.
Definition: PtrList.H:53
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:312
Namespace for OpenFOAM.
virtual bool read()
Read the transportProperties dictionary and update.
virtual Pair< tmp< volScalarField > > mDotP() const
Return the mass condensation and vaporisation rates as coefficients.
Definition: SchnerrSauer.C:138