FreeFOAM The Cross-Platform CFD Toolkit
RNGkEpsilon.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 "RNGkEpsilon.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(RNGkEpsilon, 0);
43 addToRunTimeSelectionTable(RASModel, RNGkEpsilon, dictionary);
44 
45 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
46 
48 (
49  const volVectorField& U,
50  const surfaceScalarField& phi,
51  transportModel& lamTransportModel
52 )
53 :
54  RASModel(typeName, U, phi, lamTransportModel),
55 
56  Cmu_
57  (
59  (
60  "Cmu",
61  coeffDict_,
62  0.0845
63  )
64  ),
65  C1_
66  (
68  (
69  "C1",
70  coeffDict_,
71  1.42
72  )
73  ),
74  C2_
75  (
77  (
78  "C2",
79  coeffDict_,
80  1.68
81  )
82  ),
83  sigmak_
84  (
86  (
87  "sigmak",
88  coeffDict_,
89  0.71942
90  )
91  ),
92  sigmaEps_
93  (
95  (
96  "sigmaEps",
97  coeffDict_,
98  0.71942
99  )
100  ),
101  eta0_
102  (
104  (
105  "eta0",
106  coeffDict_,
107  4.38
108  )
109  ),
110  beta_
111  (
113  (
114  "beta",
115  coeffDict_,
116  0.012
117  )
118  ),
119 
120  k_
121  (
122  IOobject
123  (
124  "k",
125  runTime_.timeName(),
126  mesh_,
129  ),
130  autoCreateK("k", mesh_)
131  ),
132  epsilon_
133  (
134  IOobject
135  (
136  "epsilon",
137  runTime_.timeName(),
138  mesh_,
141  ),
142  autoCreateEpsilon("epsilon", mesh_)
143  ),
144  nut_
145  (
146  IOobject
147  (
148  "nut",
149  runTime_.timeName(),
150  mesh_,
153  ),
154  autoCreateNut("nut", mesh_)
155  )
156 {
157  nut_ = Cmu_*sqr(k_)/(epsilon_ + epsilonSmall_);
158  nut_.correctBoundaryConditions();
159 
160  printCoeffs();
161 }
162 
163 
164 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
165 
166 
168 {
170  (
172  (
173  IOobject
174  (
175  "R",
176  runTime_.timeName(),
177  mesh_,
180  ),
181  ((2.0/3.0)*I)*k_ - nut_*twoSymm(fvc::grad(U_)),
182  k_.boundaryField().types()
183  )
184  );
185 }
186 
187 
189 {
191  (
193  (
194  IOobject
195  (
196  "devRhoReff",
197  runTime_.timeName(),
198  mesh_,
201  ),
203  )
204  );
205 }
206 
207 
209 {
210  return
211  (
212  - fvm::laplacian(nuEff(), U)
213  - fvc::div(nuEff()*dev(fvc::grad(U)().T()))
214  );
215 }
216 
217 
219 {
220  if (RASModel::read())
221  {
222  Cmu_.readIfPresent(coeffDict());
223  C1_.readIfPresent(coeffDict());
224  C2_.readIfPresent(coeffDict());
225  sigmak_.readIfPresent(coeffDict());
226  sigmaEps_.readIfPresent(coeffDict());
227  eta0_.readIfPresent(coeffDict());
228  beta_.readIfPresent(coeffDict());
229 
230  return true;
231  }
232  else
233  {
234  return false;
235  }
236 }
237 
238 
240 {
242 
243  if (!turbulence_)
244  {
245  return;
246  }
247 
249 
250  volScalarField G("RASModel::G", nut_*S2);
251 
252  volScalarField eta = sqrt(S2)*k_/epsilon_;
253  volScalarField R =
254  ((eta*(scalar(1) - eta/eta0_))/(scalar(1) + beta_*eta*sqr(eta)));
255 
256  // Update espsilon and G at the wall
257  epsilon_.boundaryField().updateCoeffs();
258 
259  // Dissipation equation
260  tmp<fvScalarMatrix> epsEqn
261  (
262  fvm::ddt(epsilon_)
263  + fvm::div(phi_, epsilon_)
264  - fvm::laplacian(DepsilonEff(), epsilon_)
265  ==
266  (C1_ - R)*G*epsilon_/k_
267  //- fvm::SuSp(R*G/k_, epsilon_)
268  - fvm::Sp(C2_*epsilon_/k_, epsilon_)
269  );
270 
271  epsEqn().relax();
272 
273  epsEqn().boundaryManipulate(epsilon_.boundaryField());
274 
275  solve(epsEqn);
276  bound(epsilon_, epsilon0_);
277 
278 
279  // Turbulent kinetic energy equation
280 
282  (
283  fvm::ddt(k_)
284  + fvm::div(phi_, k_)
285  - fvm::laplacian(DkEff(), k_)
286  ==
287  G - fvm::Sp(epsilon_/k_, k_)
288  );
289 
290  kEqn().relax();
291  solve(kEqn);
292  bound(k_, k0_);
293 
294 
295  // Re-calculate viscosity
296  nut_ = Cmu_*sqr(k_)/epsilon_;
298 }
299 
300 
301 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
302 
303 } // End namespace RASModels
304 } // End namespace incompressible
305 } // End namespace Foam
306 
307 // ************************ vim: set sw=4 sts=4 et: ************************ //