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-2023 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 compressible
35 {
36 namespace cavitationModels
37 {
40  (
44  );
45 }
46 }
47 }
48 
49 
50 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
51 
53 (
54  const dictionary& dict,
55  const compressibleTwoPhases& phases
56 )
57 :
58  cavitationModel(dict, phases),
59 
60  n_("n", dimless/dimVolume, dict),
61  dNuc_("dNuc", dimLength, dict),
62  Cv_("Cv", dimless, dict),
63  Cc_("Cc", dimless, dict),
64 
65  p0_("0", dimPressure, 0.0)
66 {
67  correct();
68 }
69 
70 
71 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
72 
74 Foam::compressible::cavitationModels::SchnerrSauer::rRb
75 (
76  const volScalarField::Internal& limitedAlphal
77 ) const
78 {
79  return pow
80  (
82  *limitedAlphal/(1 + alphaNuc() - limitedAlphal),
83  1.0/3.0
84  );
85 }
86 
87 
89 Foam::compressible::cavitationModels::SchnerrSauer::alphaNuc() const
90 {
92  return Vnuc/(1 + Vnuc);
93 }
94 
95 
97 Foam::compressible::cavitationModels::SchnerrSauer::pCoeff
98 (
100  const volScalarField::Internal& pSat
101 ) const
102 {
103  const volScalarField::Internal limitedAlphal
104  (
105  min(max(alphal(), scalar(0)), scalar(1))
106  );
107 
109  (
110  limitedAlphal*rhol() + (1 - limitedAlphal)*rhov()
111  );
112 
113  return
114  (3*rhol()*rhov())*sqrt((2.0/3.0)/rhol())
115  *rRb(limitedAlphal)/(rho*sqrt(mag(p - pSat) + 0.01*pSat));
116 }
117 
118 
121 {
122  const volScalarField::Internal& p =
123  phases_.mesh().lookupObject<volScalarField>("p");
124 
125  const volScalarField::Internal limitedAlphal
126  (
127  min(max(alphal(), scalar(0)), scalar(1))
128  );
129 
130  const volScalarField::Internal pSatv(this->pSatv());
131  const volScalarField::Internal pSatl(this->pSatl());
132 
134  (
135  Cc_*limitedAlphal*pCoeff(p, pSatv)*max(p - pSatv, p0_),
136  -Cv_
137  *(1 + alphaNuc() - limitedAlphal)
138  *pCoeff(p, pSatl)
139  *min(p - pSatl, p0_)
140  );
141 }
142 
143 
146 {
147  const volScalarField::Internal& p =
148  phases_.mesh().lookupObject<volScalarField>("p");
149 
150  const volScalarField::Internal limitedAlphal
151  (
152  min(max(alphal(), scalar(0)), scalar(1))
153  );
154 
155  const volScalarField::Internal pSatv(this->pSatv());
156  const volScalarField::Internal pSatl(this->pSatl());
157 
159  (
160  Cc_*(1 - limitedAlphal)*pos0(p - pSatv)*limitedAlphal*pCoeff(p, pSatv),
161  -Cv_
162  *(1 + alphaNuc() - limitedAlphal)
163  *neg(p - pSatl)
164  *limitedAlphal
165  *pCoeff(p, pSatl)
166  );
167 }
168 
169 
171 {}
172 
173 
175 (
176  const dictionary& dict
177 )
178 {
180  {
181  dict.lookup("n") >> n_;
182  dict.lookup("dNuc") >> dNuc_;
183  dict.lookup("Cv") >> Cv_;
184  dict.lookup("Cc") >> Cc_;
185 
186  return true;
187  }
188  else
189  {
190  return false;
191  }
192 }
193 
194 
195 // ************************************************************************* //
Macros for easy insertion into run-time selection tables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Generic GeometricField class.
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: Pair.H:65
virtual bool read(const dictionary &dict)=0
Read the dictionary and update.
virtual Pair< tmp< volScalarField::Internal > > mDotcvAlphal() const
Return the mass condensation and vaporisation rates as a.
Definition: SchnerrSauer.C:120
virtual void correct()
Correct the SchnerrSauer phaseChange model.
Definition: SchnerrSauer.C:170
SchnerrSauer(const dictionary &dict, const compressibleTwoPhases &phases)
Construct for phases.
Definition: SchnerrSauer.C:53
virtual bool read(const dictionary &dict)
Read the dictionary and update.
Definition: SchnerrSauer.C:175
virtual Pair< tmp< volScalarField::Internal > > mDotcvP() const
Return the mass condensation and vaporisation rates as coefficients.
Definition: SchnerrSauer.C:145
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
A class for managing temporary objects.
Definition: tmp.H:55
addToRunTimeSelectionTable(cavitationModel, Kunz, dictionary)
Namespace for OpenFOAM.
dimensionedScalar pos0(const dimensionedScalar &ds)
const dimensionSet dimPressure
dimensionedScalar pow3(const dimensionedScalar &ds)
const dimensionSet dimless
const dimensionSet dimLength
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
const dimensionSet dimVolume
dimensioned< scalar > mag(const dimensioned< Type > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar neg(const dimensionedScalar &ds)
dictionary dict
volScalarField & p