cavitationModel.H
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) 2021-2022 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 Class
25  Foam::cavitationModel
26 
27 Description
28  Abstract base class for cavitation models
29 
30 SourceFiles
31  cavitationModel.C
32 
33 \*---------------------------------------------------------------------------*/
34 
35 #ifndef cavitationModel_H
36 #define cavitationModel_H
37 
39 #include "fvMatricesFwd.H"
40 #include "Pair.H"
41 
42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43 
44 namespace Foam
45 {
46 
47 /*---------------------------------------------------------------------------*\
48  Class cavitationModel
49 \*---------------------------------------------------------------------------*/
50 
51 class cavitationModel
52 {
53 protected:
54 
55  // Protected data
56 
57  //- Mixture properties
59 
60  //- Index of the liquid
61  const bool liquidIndex_;
62 
63  //- Saturation vapour pressure
65 
66 
67  // Protected Member Functions
68 
69  //- Return the liquid density
70  inline const volScalarField::Internal& alphal() const
71  {
72  return phases_.alpha(liquidIndex_);
73  }
74 
75  //- Return the vapour density
76  inline const volScalarField::Internal& alphav() const
77  {
78  return phases_.alpha(!liquidIndex_);
79  }
80 
81  //- Return the liquid density
82  inline const dimensionedScalar& rhol() const
83  {
84  return phases_.rho(liquidIndex_);
85  }
86 
87  //- Return the vapour density
88  inline const dimensionedScalar& rhov() const
89  {
90  return phases_.rho(!liquidIndex_);
91  }
92 
93 
94 public:
95 
96  //- Runtime type information
97  TypeName("cavitationModel");
98 
99 
100  // Declare run-time constructor selection table
101 
103  (
104  autoPtr,
106  dictionary,
107  (
108  const dictionary& dict,
109  const incompressibleTwoPhases& phases
110  ),
111  (dict, phases)
112  );
113 
114 
115  // Constructors
116 
117  //- Construct for phases
119  (
120  const dictionary& dict,
121  const incompressibleTwoPhases& phases
122  );
123 
124 
125  // Selector
127  (
128  const dictionary& dict,
129  const incompressibleTwoPhases& phases
130  );
131 
132 
133  //- Destructor
134  virtual ~cavitationModel()
135  {}
136 
137 
138  // Member Functions
139 
140  //- Return the saturation vapour pressure
141  inline const dimensionedScalar& pSat() const
142  {
143  return pSat_;
144  }
145 
146  //- Return the mass condensation and vaporisation rates as a
147  // coefficient to multiply alphav for the condensation rate and a
148  // coefficient to multiply alphal for the vaporisation rate
149  virtual Pair<tmp<volScalarField::Internal>> mDotcvAlpha() const = 0;
150 
151  //- Return the mass condensation and vaporisation rates as coefficients
152  // to multiply (p - pSat)
153  virtual Pair<tmp<volScalarField::Internal>> mDotcvP() const = 0;
154 
155  //- Return the mass transfer rates of the two phases as coefficients to
156  // multiply the volume fraction of the other phase
158  {
160  }
161 
162  //- Return the mass transfer rates of the two phases as coefficients to
163  // multiply (p - pSat)
165  {
166  return liquidIndex_ ? reverse(mDotcvP()) : mDotcvP();
167  }
168 
169  //- Correct the cavitation model
170  virtual void correct() = 0;
171 
172  //- Read the dictionary and update
173  virtual bool read(const dictionary& dict) = 0;
174 };
175 
176 
177 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
178 
179 } // End namespace Foam
180 
181 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
182 
183 #endif
184 
185 // ************************************************************************* //
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
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
Abstract base class for cavitation models.
Pair< tmp< volScalarField::Internal > > mDot12P() const
Return the mass transfer rates of the two phases as coefficients to.
static autoPtr< cavitationModel > New(const dictionary &dict, const incompressibleTwoPhases &phases)
virtual Pair< tmp< volScalarField::Internal > > mDotcvP() const =0
Return the mass condensation and vaporisation rates as coefficients.
const dimensionedScalar & rhol() const
Return the liquid density.
virtual bool read(const dictionary &dict)=0
Read the dictionary and update.
dimensionedScalar pSat_
Saturation vapour pressure.
Pair< tmp< volScalarField::Internal > > mDot12Alpha() const
Return the mass transfer rates of the two phases as coefficients to.
const volScalarField::Internal & alphav() const
Return the vapour density.
cavitationModel(const dictionary &dict, const incompressibleTwoPhases &phases)
Construct for phases.
virtual void correct()=0
Correct the cavitation model.
const dimensionedScalar & pSat() const
Return the saturation vapour pressure.
declareRunTimeSelectionTable(autoPtr, cavitationModel, dictionary,(const dictionary &dict, const incompressibleTwoPhases &phases),(dict, phases))
const bool liquidIndex_
Index of the liquid.
virtual Pair< tmp< volScalarField::Internal > > mDotcvAlpha() const =0
Return the mass condensation and vaporisation rates as a.
const incompressibleTwoPhases & phases_
Mixture properties.
virtual ~cavitationModel()
Destructor.
TypeName("cavitationModel")
Runtime type information.
const volScalarField::Internal & alphal() const
Return the liquid density.
const dimensionedScalar & rhov() const
Return the vapour density.
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.
const dimensionedScalar & rho(const bool index) const
Return the density of a given phase.
const volScalarField & alpha(const bool index) const
Return the volume-fraction of a given phase.
Definition: twoPhases.H:103
Forward declarations of fvMatrix specialisations.
Namespace for OpenFOAM.
void reverse(UList< T > &, const label n)
Definition: UListI.H:334
dictionary dict