Saito.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) 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 "Saito.H"
28 
29 #include "mathematicalConstants.H"
31 
34 
35 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39 namespace compressible
40 {
41 namespace cavitationModels
42 {
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  Ca_("Ca", dimless/dimLength, dict),
61  Cv_("Cv", dimless, dict),
62  Cc_("Cc", dimless, dict),
63  alphaNuc_("alphaNuc", dimless, dict),
64 
65  p0_("0", dimPressure, 0)
66 {
67  correct();
68 }
69 
70 
71 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
72 
74 Foam::compressible::cavitationModels::Saito::fT(const rhoThermo& thermo) const
75 {
76  return sqrt(2*pi*(RR/thermo.W()()())*thermo.T()());
77 }
78 
79 
82 {
83  const volScalarField::Internal& p = thermol().p();
84 
85  const volScalarField::Internal alphav
86  (
87  min(max(this->alphav(), scalar(0)), scalar(1))
88  );
89 
90  const volScalarField::Internal alphal
91  (
92  min(max(this->alphal(), scalar(0)), scalar(1))
93  );
94 
95  const volScalarField::Internal alphavNuc(max(alphav, alphaNuc_));
96 
97  const volScalarField::Internal A(Ca_*alphal*alphavNuc);
98 
99  const volScalarField::Internal mvCoeff(Cv_*A*rhol()/(rhov()*fT(thermol())));
100  const volScalarField::Internal mcCoeff(Cc_*A/fT(thermol()));
101 
103  (
104  mcCoeff*alphal*max(p - pSatv(), p0_),
105  -mvCoeff*alphavNuc*min(p - pSatl(), p0_)
106  );
107 }
108 
109 
112 {
113  const volScalarField::Internal& p = thermol().p();
114 
115  const volScalarField::Internal alphav
116  (
117  min(max(this->alphav(), scalar(0)), scalar(1))
118  );
119 
120  const volScalarField::Internal alphal
121  (
122  min(max(this->alphal(), scalar(0)), scalar(1))
123  );
124 
125  const volScalarField::Internal alphavNuc(max(alphav, alphaNuc_));
126 
127  const volScalarField::Internal A(Ca_*alphal*alphavNuc);
128 
129  const volScalarField::Internal mvCoeff(Cv_*A*rhol()/(rhov()*fT(thermol())));
130  const volScalarField::Internal mcCoeff(Cc_*A/fT(thermol()));
131 
133  (
134  mcCoeff*alphal*alphav*pos0(p - pSatv()),
135  -mvCoeff*alphal*alphavNuc*neg(p - pSatl())
136  );
137 }
138 
139 
141 {}
142 
143 
145 (
146  const dictionary& dict
147 )
148 {
150  {
151  dict.lookup("Ca") >> Ca_;
152  dict.lookup("Cv") >> Cv_;
153  dict.lookup("Cc") >> Cc_;
154  dict.lookup("alphaNuc") >> alphaNuc_;
155 
156  return true;
157  }
158  else
159  {
160  return false;
161  }
162 }
163 
164 
165 // ************************************************************************* //
static const Foam::dimensionedScalar A("A", Foam::dimPressure, 611.21)
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...
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.
Saito(const dictionary &dict, const compressibleTwoPhases &phases)
Construct for phases.
Definition: Saito.C:53
virtual Pair< tmp< volScalarField::Internal > > mDotcvAlphal() const
Return the mass condensation and vaporisation rates as a.
Definition: Saito.C:81
virtual void correct()
Correct the Saito phaseChange model.
Definition: Saito.C:140
virtual bool read(const dictionary &dict)
Read the dictionary and update.
Definition: Saito.C:145
virtual Pair< tmp< volScalarField::Internal > > mDotcvP() const
Return the mass condensation and vaporisation rates as coefficients.
Definition: Saito.C:111
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:160
Base-class for fluid thermodynamic properties based on density.
Definition: rhoThermo.H:55
A class for managing temporary objects.
Definition: tmp.H:55
addToRunTimeSelectionTable(cavitationModel, Kunz, dictionary)
const dimensionedScalar RR
Universal gas constant: default SI units: [J/kmol/K].
Namespace for OpenFOAM.
dimensionedScalar pos0(const dimensionedScalar &ds)
const dimensionSet dimPressure
const dimensionSet dimless
const dimensionSet dimLength
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar sqrt(const dimensionedScalar &ds)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar neg(const dimensionedScalar &ds)
dictionary dict
volScalarField & p
fluidMulticomponentThermo & thermo
Definition: createFields.H:31