kineticTheoryModel.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) 2011-2021 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::kineticTheoryModel
26 
27 Description
28  Kinetic theory particle phase RAS model
29 
30  Reference:
31  \verbatim
32  van Wachem, B. G. M. (2000).
33  Derivation, implementation, and validation of computer simulation models
34  for gas-solid fluidised beds.
35  PhD Thesis, TU Delft.
36  \endverbatim
37 
38  There are no default model coefficients.
39 
40 SourceFiles
41  kineticTheoryModel.C
42 
43 \*---------------------------------------------------------------------------*/
44 
45 #ifndef kineticTheoryModel_H
46 #define kineticTheoryModel_H
47 
48 #include "RASModel.H"
49 #include "eddyViscosity.H"
51 #include "phaseModel.H"
52 #include "dragModel.H"
54 #include "conductivityModel.H"
55 #include "radialModel.H"
56 #include "granularPressureModel.H"
57 #include "frictionalStressModel.H"
58 
59 
60 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
61 
62 namespace Foam
63 {
64 namespace RASModels
65 {
66 
67 /*---------------------------------------------------------------------------*\
68  Class kineticTheoryModel Declaration
69 \*---------------------------------------------------------------------------*/
70 
72 :
73  public eddyViscosity<RASModel<phaseCompressible::momentumTransportModel>>
74 {
75  // Private Data
76 
77  // Input Fields
78 
79  const phaseModel& phase_;
80 
81  //- Name of the continuous phase
82  word continuousPhaseName_;
83 
84 
85  // Sub-models
86 
87  //- Run-time selected viscosity model
89 
90  //- Run-time selected conductivity model
92 
93  //- Run-time selected radial distribution model
95 
96  //- Run-time selected granular pressure model
98  granularPressureModel_;
99 
100  //- Run-time selected frictional stress model
102  frictionalStressModel_;
103 
104 
105  // Kinetic Theory Model coefficients
106 
107  //- Use equilibrium approximation: generation == dissipation
108  Switch equilibrium_;
109 
110  //- Coefficient of restitution
112 
113  //- Maximum packing phase-fraction
114  dimensionedScalar alphaMax_;
115 
116  //- Min value for which the frictional stresses are zero
117  dimensionedScalar alphaMinFriction_;
118 
119  //- Residual phase fraction
120  dimensionedScalar residualAlpha_;
121 
122  //- Maximum turbulent viscosity
123  dimensionedScalar maxNut_;
124 
125 
126  // Kinetic Theory Model Fields
127 
128  //- The granular energy/temperature
129  volScalarField Theta_;
130 
131  //- The granular bulk viscosity
132  volScalarField lambda_;
133 
134  //- The granular radial distribution
135  volScalarField gs0_;
136 
137  //- The granular "thermal" conductivity
138  volScalarField kappa_;
139 
140  //- The frictional viscosity
141  volScalarField nuFric_;
142 
143 
144  // Private Member Functions
145 
146  void correctNut()
147  {}
148 
149  //- Return the continuous phase model
150  // which for two-phases is the "other" phase
151  // and for more than two phases must be specified
152  const phaseModel& continuousPhase() const;
153 
154 
155 public:
156 
157  //- Runtime type information
158  TypeName("kineticTheory");
159 
160 
161  // Constructors
162 
163  //- Construct from components
165  (
166  const volScalarField& alpha,
167  const volScalarField& rho,
168  const volVectorField& U,
170  const surfaceScalarField& phi,
171  const transportModel& transport,
172  const word& type = typeName
173  );
174 
175  //- Disallow default bitwise copy construction
176  kineticTheoryModel(const kineticTheoryModel&) = delete;
177 
178 
179  //- Destructor
180  virtual ~kineticTheoryModel();
181 
182 
183  // Member Functions
184 
185  //- Re-read model coefficients if they have changed
186  virtual bool read();
187 
188  //- Return the effective viscosity
189  virtual tmp<volScalarField> nuEff() const
190  {
191  return this->nut();
192  }
193 
194  //- Return the effective viscosity on patch
195  virtual tmp<scalarField> nuEff(const label patchi) const
196  {
197  return this->nut(patchi);
198  }
199 
200  //- Return the turbulence kinetic energy
201  virtual tmp<volScalarField> k() const;
202 
203  //- Return the turbulence kinetic energy dissipation rate
204  virtual tmp<volScalarField> epsilon() const;
205 
206  //- Return the stress tensor [m^2/s^2]
207  virtual tmp<volSymmTensorField> sigma() const;
208 
209  //- Return the phase-pressure'
210  // (derivative of phase-pressure w.r.t. phase-fraction)
211  virtual tmp<volScalarField> pPrime() const;
212 
213  //- Return the face-phase-pressure'
214  // (derivative of phase-pressure w.r.t. phase-fraction)
215  virtual tmp<surfaceScalarField> pPrimef() const;
216 
217  //- Return the effective stress tensor
218  virtual tmp<volSymmTensorField> devTau() const;
219 
220  //- Return the source term for the momentum equation
221  virtual tmp<fvVectorMatrix> divDevTau(volVectorField& U) const;
222 
223  //- Solve the kinetic theory equations and correct the viscosity
224  virtual void correct();
225 
226 
227  // Member Operators
228 
229  //- Disallow default bitwise assignment
230  void operator=(const kineticTheoryModel&) = delete;
231 };
232 
233 
234 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
235 
236 } // End namespace RASModels
237 } // End namespace Foam
238 
239 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
240 
241 #endif
242 
243 // ************************************************************************* //
virtual tmp< volSymmTensorField > devTau() const
Return the effective stress tensor.
const transportModel & transport() const
Access function to incompressible transport model.
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate.
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< volScalarField > nut() const
Return the turbulence viscosity.
virtual tmp< volSymmTensorField > sigma() const
Return the stress tensor [m^2/s^2].
virtual tmp< surfaceScalarField > pPrimef() const
Return the face-phase-pressure&#39;.
virtual tmp< fvVectorMatrix > divDevTau(volVectorField &U) const
Return the source term for the momentum equation.
TypeName("kineticTheory")
Runtime type information.
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none/any.
Definition: Switch.H:60
const volScalarField & rho() const
Return the density field.
Eddy viscosity turbulence model base class.
Definition: eddyViscosity.H:49
virtual tmp< volScalarField > nuEff() const
Return the effective viscosity.
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.
void operator=(const kineticTheoryModel &)=delete
Disallow default bitwise assignment.
virtual ~kineticTheoryModel()
Destructor.
Kinetic theory particle phase RAS model.
const volVectorField & U() const
Access function to velocity field.
const surfaceScalarField & alphaRhoPhi() const
Access function to phase flux field.
label patchi
virtual tmp< volScalarField > pPrime() const
Return the phase-pressure&#39;.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
virtual bool read()
Re-read model coefficients if they have changed.
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
A class for managing temporary objects.
Definition: PtrList.H:53
kineticTheoryModel(const volScalarField &alpha, const volScalarField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &type=typeName)
Construct from components.
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:53
Namespace for OpenFOAM.