FreeFOAM The Cross-Platform CFD Toolkit
RASModel.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
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 Namespace
25  Foam::compressible::RASModels
26 
27 Description
28  Namespace for compressible RAS turbulence models.
29 
30 
31 Class
32  Foam::compressible::RASModel
33 
34 Description
35  Abstract base class for turbulence models for compressible and combusting
36  flows.
37 
38 SourceFiles
39  RASModel.C
40  newTurbulenceModel.C
41 
42 \*---------------------------------------------------------------------------*/
43 
44 #ifndef compressibleRASModel_H
45 #define compressibleRASModel_H
46 
48 #include <finiteVolume/volFields.H>
51 #include <finiteVolume/fvm.H>
52 #include <finiteVolume/fvc.H>
55 #include <OpenFOAM/IOdictionary.H>
56 #include <OpenFOAM/Switch.H>
57 #include <finiteVolume/bound.H>
58 #include <OpenFOAM/autoPtr.H>
60 
61 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
62 
63 namespace Foam
64 {
65 namespace compressible
66 {
67 
68 /*---------------------------------------------------------------------------*\
69  Class RASModel Declaration
70 \*---------------------------------------------------------------------------*/
71 
72 class RASModel
73 :
74  public turbulenceModel,
75  public IOdictionary
76 {
77 
78 protected:
79 
80  // Protected data
81 
82  //- Turbulence on/off flag
84 
85  //- Flag to print the model coeffs at run-time
87 
88  //- Model coefficients dictionary
90 
91  //- Value of y+ at the edge of the laminar sublayer
92  scalar yPlusLam_;
93 
94  //- Lower limit of k
96 
97  //- Lower limit of epsilon
99 
100  //- Small epsilon value used to avoid divide by zero
102 
103  //- Lower limit for omega
105 
106  //- Small omega value used to avoid divide by zero
108 
109  //- Near wall distance boundary field
111 
112 
113  // Protected member functions
114 
115  //- Print model coefficients
116  virtual void printCoeffs();
117 
118 
119 private:
120 
121  // Private Member Functions
122 
123  //- Disallow default bitwise copy construct
124  RASModel(const RASModel&);
125 
126  //- Disallow default bitwise assignment
127  void operator=(const RASModel&);
128 
129 
130 public:
131 
132  //- Runtime type information
133  TypeName("RASModel");
134 
135 
136  // Declare run-time constructor selection table
137 
139  (
140  autoPtr,
141  RASModel,
142  dictionary,
143  (
144  const volScalarField& rho,
145  const volVectorField& U,
146  const surfaceScalarField& phi,
147  const basicThermo& thermoPhysicalModel
148  ),
149  (rho, U, phi, thermoPhysicalModel)
150  );
151 
152 
153  // Constructors
154 
155  //- Construct from components
156  RASModel
157  (
158  const word& type,
159  const volScalarField& rho,
160  const volVectorField& U,
161  const surfaceScalarField& phi,
162  const basicThermo& thermoPhysicalModel
163  );
164 
165 
166  // Selectors
167 
168  //- Return a reference to the selected turbulence model
169  static autoPtr<RASModel> New
170  (
171  const volScalarField& rho,
172  const volVectorField& U,
173  const surfaceScalarField& phi,
174  const basicThermo& thermoPhysicalModel
175  );
176 
177 
178  //- Destructor
179  virtual ~RASModel()
180  {}
181 
182 
183  // Member Functions
184 
185  // Access
186 
187  //- Return the value of k0 which k is not allowed to be less than
188  const dimensionedScalar& k0() const
189  {
190  return k0_;
191  }
192 
193  //- Return the value of epsilon0 which epsilon is not allowed to be
194  // less than
195  const dimensionedScalar& epsilon0() const
196  {
197  return epsilon0_;
198  }
199 
200  //- Return the value of epsilonSmall which is added to epsilon when
201  // calculating nut
202  const dimensionedScalar& epsilonSmall() const
203  {
204  return epsilonSmall_;
205  }
206 
207  //- Return the value of omega0 which epsilon is not allowed to be
208  // less than
209  const dimensionedScalar& omega0() const
210  {
211  return omega0_;
212  }
213 
214  //- Return the value of omegaSmall which is added to epsilon when
215  // calculating nut
216  const dimensionedScalar& omegaSmall() const
217  {
218  return omegaSmall_;
219  }
220 
221  //- Allow k0 to be changed
223  {
224  return k0_;
225  }
226 
227  //- Allow epsilon0 to be changed
229  {
230  return epsilon0_;
231  }
232 
233  //- Allow epsilonSmall to be changed
235  {
236  return epsilonSmall_;
237  }
238 
239  //- Allow omega0 to be changed
241  {
242  return omega0_;
243  }
244 
245  //- Allow omegaSmall to be changed
247  {
248  return omegaSmall_;
249  }
250 
251  //- Return the near wall distances
252  const nearWallDist& y() const
253  {
254  return y_;
255  }
256 
257  //- Calculate y+ at the edge of the laminar sublayer
258  scalar yPlusLam(const scalar kappa, const scalar E) const;
259 
260  //- Const access to the coefficients dictionary
261  const dictionary& coeffDict() const
262  {
263  return coeffDict_;
264  }
265 
266 
267  //- Return the turbulence viscosity
268  virtual tmp<volScalarField> mut() const = 0;
269 
270  //- Return the effective viscosity
271  virtual tmp<volScalarField> muEff() const
272  {
273  return tmp<volScalarField>
274  (
275  new volScalarField("muEff", mut() + mu())
276  );
277  }
278 
279  //- Return the effective turbulent thermal diffusivity
280  virtual tmp<volScalarField> alphaEff() const = 0;
281 
282  //- Return the turbulence kinetic energy
283  virtual tmp<volScalarField> k() const = 0;
284 
285  //- Return the turbulence kinetic energy dissipation rate
286  virtual tmp<volScalarField> epsilon() const = 0;
287 
288  //- Return the Reynolds stress tensor
289  virtual tmp<volSymmTensorField> R() const = 0;
290 
291  //- Return the effective stress tensor including the laminar stress
292  virtual tmp<volSymmTensorField> devRhoReff() const = 0;
293 
294  //- Return the source term for the momentum equation
295  virtual tmp<fvVectorMatrix> divDevRhoReff(volVectorField& U) const = 0;
296 
297  //- Return yPlus for the given patch
298  virtual tmp<scalarField> yPlus
299  (
300  const label patchI,
301  const scalar Cmu
302  ) const;
303 
304  //- Solve the turbulence equations and correct the turbulence viscosity
305  virtual void correct() = 0;
306 
307  //- Read RASProperties dictionary
308  virtual bool read() = 0;
309 };
310 
311 
312 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
313 
314 } // End namespace compressible
315 } // End namespace Foam
316 
317 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
318 
319 #endif
320 
321 // ************************ vim: set sw=4 sts=4 et: ************************ //