FreeFOAM The Cross-Platform CFD Toolkit
LaunderSharmaKE.C
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 \*---------------------------------------------------------------------------*/
25 
26 #include "LaunderSharmaKE.H"
28 
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace incompressible
36 {
37 namespace RASModels
38 {
39 
40 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
41 
42 defineTypeNameAndDebug(LaunderSharmaKE, 0);
43 addToRunTimeSelectionTable(RASModel, LaunderSharmaKE, dictionary);
44 
45 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
46 
47 tmp<volScalarField> LaunderSharmaKE::fMu() const
48 {
49  return exp(-3.4/sqr(scalar(1) + sqr(k_)/(nu()*epsilonTilda_)/50.0));
50 }
51 
52 
53 tmp<volScalarField> LaunderSharmaKE::f2() const
54 {
55  return
56  scalar(1)
57  - 0.3*exp(-min(sqr(sqr(k_)/(nu()*epsilonTilda_)), scalar(50.0)));
58 }
59 
60 
61 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
62 
64 (
65  const volVectorField& U,
66  const surfaceScalarField& phi,
67  transportModel& lamTransportModel
68 )
69 :
70  RASModel(typeName, U, phi, lamTransportModel),
71 
72  Cmu_
73  (
75  (
76  "Cmu",
77  coeffDict_,
78  0.09
79  )
80  ),
81  C1_
82  (
84  (
85  "C1",
86  coeffDict_,
87  1.44
88  )
89  ),
90  C2_
91  (
93  (
94  "C2",
95  coeffDict_,
96  1.92
97  )
98  ),
99  sigmaEps_
100  (
102  (
103  "sigmaEps",
104  coeffDict_,
105  1.3
106  )
107  ),
108 
109  k_
110  (
111  IOobject
112  (
113  "k",
114  runTime_.timeName(),
115  mesh_,
118  ),
119  mesh_
120  ),
121 
122  epsilonTilda_
123  (
124  IOobject
125  (
126  "epsilon",
127  runTime_.timeName(),
128  mesh_,
131  ),
132  mesh_
133  ),
134 
135  nut_
136  (
137  IOobject
138  (
139  "nut",
140  runTime_.timeName(),
141  mesh_,
144  ),
145  autoCreateLowReNut("nut", mesh_)
146  )
147 {
148  nut_ = Cmu_*fMu()*sqr(k_)/(epsilonTilda_ + epsilonSmall_);
149  nut_.correctBoundaryConditions();
150 
151  printCoeffs();
152 }
153 
154 
155 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
156 
158 {
160  (
162  (
163  IOobject
164  (
165  "R",
166  runTime_.timeName(),
167  mesh_,
170  ),
171  ((2.0/3.0)*I)*k_ - nut_*twoSymm(fvc::grad(U_)),
172  k_.boundaryField().types()
173  )
174  );
175 }
176 
177 
179 {
181  (
183  (
184  IOobject
185  (
186  "devRhoReff",
187  runTime_.timeName(),
188  mesh_,
191  ),
193  )
194  );
195 }
196 
197 
199 {
200  return
201  (
202  - fvm::laplacian(nuEff(), U)
203  - fvc::div(nuEff()*dev(fvc::grad(U)().T()))
204  );
205 }
206 
207 
209 {
210  if (RASModel::read())
211  {
212  Cmu_.readIfPresent(coeffDict());
213  C1_.readIfPresent(coeffDict());
214  C2_.readIfPresent(coeffDict());
215  sigmaEps_.readIfPresent(coeffDict());
216 
217  return true;
218  }
219  else
220  {
221  return false;
222  }
223 }
224 
225 
227 {
229 
230  if (!turbulence_)
231  {
232  return;
233  }
234 
236 
237  volScalarField G("RASModel::G", nut_*S2);
238 
239  volScalarField E = 2.0*nu()*nut_*fvc::magSqrGradGrad(U_);
240  volScalarField D = 2.0*nu()*magSqr(fvc::grad(sqrt(k_)));
241 
242 
243  // Dissipation rate equation
244 
245  tmp<fvScalarMatrix> epsEqn
246  (
247  fvm::ddt(epsilonTilda_)
248  + fvm::div(phi_, epsilonTilda_)
249  - fvm::laplacian(DepsilonEff(), epsilonTilda_)
250  ==
251  C1_*G*epsilonTilda_/k_
252  - fvm::Sp(C2_*f2()*epsilonTilda_/k_, epsilonTilda_)
253  + E
254  );
255 
256  epsEqn().relax();
257  solve(epsEqn);
258  bound(epsilonTilda_, epsilon0_);
259 
260 
261  // Turbulent kinetic energy equation
262 
264  (
265  fvm::ddt(k_)
266  + fvm::div(phi_, k_)
267  - fvm::laplacian(DkEff(), k_)
268  ==
269  G - fvm::Sp((epsilonTilda_ + D)/k_, k_)
270  );
271 
272  kEqn().relax();
273  solve(kEqn);
274  bound(k_, k0_);
275 
276 
277  // Re-calculate viscosity
278  nut_ == Cmu_*fMu()*sqr(k_)/epsilonTilda_;
279 }
280 
281 
282 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
283 
284 } // End namespace RASModels
285 } // End namespace incompressible
286 } // End namespace Foam
287 
288 // ************************ vim: set sw=4 sts=4 et: ************************ //