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 cavitationModels
35 {
38  (
42  );
43 }
44 }
45 
46 
47 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
48 
50 (
51  const dictionary& dict,
52  const incompressibleTwoPhases& phases
53 )
54 :
55  cavitationModel(dict, phases),
56 
57  n_("n", dimless/dimVolume, dict),
58  dNuc_("dNuc", dimLength, dict),
59  Cv_("Cv", dimless, dict),
60  Cc_("Cc", dimless, dict),
61 
62  p0_("0", pSat().dimensions(), 0.0)
63 {
64  correct();
65 }
66 
67 
68 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
69 
71 Foam::cavitationModels::SchnerrSauer::rRb
72 (
73  const volScalarField::Internal& limitedAlphal
74 ) const
75 {
76  return pow
77  (
79  *limitedAlphal/(1.0 + alphaNuc() - limitedAlphal),
80  1.0/3.0
81  );
82 }
83 
84 
86 Foam::cavitationModels::SchnerrSauer::alphaNuc() const
87 {
89  return Vnuc/(1 + Vnuc);
90 }
91 
92 
94 Foam::cavitationModels::SchnerrSauer::pCoeff
95 (
97 ) const
98 {
99  const volScalarField::Internal limitedAlphal
100  (
101  min(max(alphal(), scalar(0)), scalar(1))
102  );
103 
105  (
106  limitedAlphal*rhol()
107  + (scalar(1) - limitedAlphal)*rhov()
108  );
109 
110  return
111  (3*rhol()*rhov())*sqrt(2/(3*rhol()))
112  *rRb(limitedAlphal)/(rho*sqrt(mag(p - pSat()) + 0.01*pSat()));
113 }
114 
115 
118 {
119  const volScalarField::Internal& p =
120  phases_.mesh().lookupObject<volScalarField>("p");
121 
122  const volScalarField::Internal pCoeff(this->pCoeff(p));
123 
124  const volScalarField::Internal limitedAlphal
125  (
126  min(max(alphal(), scalar(0)), scalar(1))
127  );
128 
130  (
131  Cc_*limitedAlphal*pCoeff*max(p - pSat(), p0_),
132  -Cv_*(1.0 + alphaNuc() - limitedAlphal)*pCoeff*min(p - pSat(), p0_)
133  );
134 }
135 
136 
139 {
140  const volScalarField::Internal& p =
141  phases_.mesh().lookupObject<volScalarField>("p");
142 
143  const volScalarField::Internal pCoeff(this->pCoeff(p));
144 
145  const volScalarField::Internal limitedAlphal
146  (
147  min(max(alphal(), scalar(0)), scalar(1))
148  );
149 
150  const volScalarField::Internal apCoeff(limitedAlphal*pCoeff);
151 
153  (
154  Cc_*(1.0 - limitedAlphal)*pos0(p - pSat())*apCoeff,
155  (-Cv_)*(1.0 + alphaNuc() - limitedAlphal)*neg(p - pSat())*apCoeff
156  );
157 }
158 
159 
161 {}
162 
163 
165 {
167  {
168  dict.lookup("n") >> n_;
169  dict.lookup("dNuc") >> dNuc_;
170  dict.lookup("Cv") >> Cv_;
171  dict.lookup("Cc") >> Cc_;
172 
173  return true;
174  }
175  else
176  {
177  return false;
178  }
179 }
180 
181 
182 // ************************************************************************* //
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.
DimensionedField< Type, GeoMesh > Internal
Type of the internal field from which this GeometricField is derived.
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: Pair.H:65
Abstract base class for cavitation models.
virtual bool read(const dictionary &dict)=0
Read the dictionary and update.
SchnerrSauer cavitation model.
Definition: SchnerrSauer.H:127
virtual void correct()
Correct the SchnerrSauer phaseChange model.
Definition: SchnerrSauer.C:160
virtual bool read(const dictionary &dict)
Read the dictionary and update.
Definition: SchnerrSauer.C:164
virtual Pair< tmp< volScalarField::Internal > > mDotcvAlpha() const
Return the mass condensation and vaporisation rates as a.
Definition: SchnerrSauer.C:117
SchnerrSauer(const dictionary &dict, const incompressibleTwoPhases &phases)
Construct for phases.
Definition: SchnerrSauer.C:50
virtual Pair< tmp< volScalarField::Internal > > mDotcvP() const
Return the mass condensation and vaporisation rates as coefficients.
Definition: SchnerrSauer.C:138
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Interface to two constant density phases.
A class for managing temporary objects.
Definition: tmp.H:55
defineTypeNameAndDebug(Kunz, 0)
addToRunTimeSelectionTable(cavitationModel, Kunz, dictionary)
Namespace for OpenFOAM.
dimensionedScalar pos0(const dimensionedScalar &ds)
const HashTable< dimensionSet > & dimensions()
Get the table of dimension sets.
Definition: dimensionSets.C:96
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