phasePressureModel.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) 2013-2020 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::RASModels::phasePressureModel
26 
27 Description
28  Particle-particle phase-pressure RAS model
29 
30  The derivative of the phase-pressure with respect to the phase-fraction
31  is evaluated as
32 
33  g0*min(exp(preAlphaExp*(alpha - alphaMax)), expMax)
34 
35  The default model coefficients correspond to the following:
36  \verbatim
37  phasePressureCoeffs
38  {
39  preAlphaExp 500;
40  expMax 1000;
41  alphaMax 0.62;
42  g0 1000;
43  }
44  \endverbatim
45 
46 SourceFiles
47  phasePressureModel.C
48 
49 \*---------------------------------------------------------------------------*/
50 
51 #ifndef phasePressureModel_H
52 #define phasePressureModel_H
53 
54 #include "RASModel.H"
55 #include "eddyViscosity.H"
57 #include "phaseModel.H"
58 
59 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
60 
61 namespace Foam
62 {
63 namespace RASModels
64 {
65 
66 /*---------------------------------------------------------------------------*\
67  Class phasePressureModel Declaration
68 \*---------------------------------------------------------------------------*/
69 
71 :
72  public eddyViscosity<RASModel<phaseCompressibleMomentumTransportModel>>
73 {
74  // Private Data
75 
76  // Kinetic Theory Model coefficients
77 
78  //- Maximum packing phase-fraction
79  scalar alphaMax_;
80 
81  //- Pre-exponential factor
82  scalar preAlphaExp_;
83 
84  //- Maximum limit of the exponential
85  scalar expMax_;
86 
87  //- g0
89 
90 
91  // Private Member Functions
92 
93  void correctNut()
94  {}
95 
96 
97 public:
98 
99  //- Runtime type information
100  TypeName("phasePressure");
101 
102 
103  // Constructors
104 
105  //- Construct from components
107  (
108  const volScalarField& alpha,
109  const volScalarField& rho,
110  const volVectorField& U,
112  const surfaceScalarField& phi,
113  const phaseModel& transport,
114  const word& type = typeName
115  );
116 
117  //- Disallow default bitwise copy construction
118  phasePressureModel(const phasePressureModel&) = delete;
119 
120 
121  //- Destructor
122  virtual ~phasePressureModel();
123 
124 
125  // Member Functions
126 
127  //- Re-read model coefficients if they have changed
128  virtual bool read();
129 
130  //- Return the effective viscosity
131  virtual tmp<volScalarField> nuEff() const
132  {
133  return this->nut();
134  }
135 
136  //- Return the effective viscosity on patch
137  virtual tmp<scalarField> nuEff(const label patchi) const
138  {
139  return this->nut(patchi);
140  }
141 
142  //- Return the turbulence kinetic energy
143  virtual tmp<volScalarField> k() const;
144 
145  //- Return the turbulence kinetic energy dissipation rate
146  virtual tmp<volScalarField> epsilon() const;
147 
148  //- Return the stress tensor [m^2/s^2]
149  virtual tmp<volSymmTensorField> sigma() const;
150 
151  //- Return the phase-pressure'
152  // (derivative of phase-pressure w.r.t. phase-fraction)
153  virtual tmp<volScalarField> pPrime() const;
154 
155  //- Return the face-phase-pressure'
156  // (derivative of phase-pressure w.r.t. phase-fraction)
157  virtual tmp<surfaceScalarField> pPrimef() const;
158 
159  //- Return the effective stress tensor
160  virtual tmp<volSymmTensorField> devTau() const;
161 
162  //- Return the source term for the momentum equation
163  virtual tmp<fvVectorMatrix> divDevTau(volVectorField& U) const;
164 
165  //- Solve the kinetic theory equations and correct the viscosity
166  virtual void correct();
167 
168 
169  // Member Operators
170 
171  //- Disallow default bitwise assignment
172  void operator=(const phasePressureModel&) = delete;
173 };
174 
175 
176 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
177 
178 } // End namespace RASModels
179 } // End namespace Foam
180 
181 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
182 
183 #endif
184 
185 // ************************************************************************* //
const transportModel & transport() const
Access function to incompressible transport model.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
virtual tmp< volSymmTensorField > devTau() const
Return the effective stress tensor.
phasePressureModel(const volScalarField &alpha, const volScalarField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const phaseModel &transport, const word &type=typeName)
Construct from components.
virtual tmp< volScalarField > nut() const
Return the turbulence viscosity.
virtual tmp< volSymmTensorField > sigma() const
Return the stress tensor [m^2/s^2].
virtual tmp< fvVectorMatrix > divDevTau(volVectorField &U) const
Return the source term for the momentum equation.
const volScalarField & rho() const
Return the density field.
Particle-particle phase-pressure RAS model.
Eddy viscosity turbulence model base class.
Definition: eddyViscosity.H:49
virtual tmp< volScalarField > pPrime() const
Return the phase-pressure&#39;.
A class for handling words, derived from string.
Definition: word.H:59
virtual tmp< surfaceScalarField > phi() const
Return the volumetric flux field.
virtual void correct()
Solve the kinetic theory equations and correct the viscosity.
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
virtual bool read()
Re-read model coefficients if they have changed.
const volVectorField & U() const
Access function to velocity field.
const surfaceScalarField & alphaRhoPhi() const
Access function to phase flux field.
virtual ~phasePressureModel()
Destructor.
label patchi
void operator=(const phasePressureModel &)=delete
Disallow default bitwise assignment.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
A class for managing temporary objects.
Definition: PtrList.H:53
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:53
TypeName("phasePressure")
Runtime type information.
virtual tmp< volScalarField > nuEff() const
Return the effective viscosity.
Namespace for OpenFOAM.
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate.
virtual tmp< surfaceScalarField > pPrimef() const
Return the face-phase-pressure&#39;.