FreeFOAM The Cross-Platform CFD Toolkit
kEpsilon.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 "kEpsilon.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(kEpsilon, 0);
43 addToRunTimeSelectionTable(RASModel, kEpsilon, 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.09
63  )
64  ),
65  C1_
66  (
68  (
69  "C1",
70  coeffDict_,
71  1.44
72  )
73  ),
74  C2_
75  (
77  (
78  "C2",
79  coeffDict_,
80  1.92
81  )
82  ),
83  sigmaEps_
84  (
86  (
87  "sigmaEps",
88  coeffDict_,
89  1.3
90  )
91  ),
92 
93  k_
94  (
95  IOobject
96  (
97  "k",
98  runTime_.timeName(),
99  mesh_,
102  ),
103  autoCreateK("k", mesh_)
104  ),
105  epsilon_
106  (
107  IOobject
108  (
109  "epsilon",
110  runTime_.timeName(),
111  mesh_,
114  ),
115  autoCreateEpsilon("epsilon", mesh_)
116  ),
117  nut_
118  (
119  IOobject
120  (
121  "nut",
122  runTime_.timeName(),
123  mesh_,
126  ),
127  autoCreateNut("nut", mesh_)
128  )
129 {
130  nut_ = Cmu_*sqr(k_)/(epsilon_ + epsilonSmall_);
131  nut_.correctBoundaryConditions();
132 
133  printCoeffs();
134 }
135 
136 
137 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
138 
140 {
142  (
144  (
145  IOobject
146  (
147  "R",
148  runTime_.timeName(),
149  mesh_,
152  ),
153  ((2.0/3.0)*I)*k_ - nut_*twoSymm(fvc::grad(U_)),
154  k_.boundaryField().types()
155  )
156  );
157 }
158 
159 
161 {
163  (
165  (
166  IOobject
167  (
168  "devRhoReff",
169  runTime_.timeName(),
170  mesh_,
173  ),
175  )
176  );
177 }
178 
179 
181 {
182  return
183  (
184  - fvm::laplacian(nuEff(), U)
185  - fvc::div(nuEff()*dev(fvc::grad(U)().T()))
186  );
187 }
188 
189 
191 {
192  if (RASModel::read())
193  {
194  Cmu_.readIfPresent(coeffDict());
195  C1_.readIfPresent(coeffDict());
196  C2_.readIfPresent(coeffDict());
197  sigmaEps_.readIfPresent(coeffDict());
198 
199  return true;
200  }
201  else
202  {
203  return false;
204  }
205 }
206 
207 
209 {
211 
212  if (!turbulence_)
213  {
214  return;
215  }
216 
217  volScalarField G("RASModel::G", nut_*2*magSqr(symm(fvc::grad(U_))));
218 
219  // Update espsilon and G at the wall
220  epsilon_.boundaryField().updateCoeffs();
221 
222  // Dissipation equation
223  tmp<fvScalarMatrix> epsEqn
224  (
225  fvm::ddt(epsilon_)
226  + fvm::div(phi_, epsilon_)
227  - fvm::Sp(fvc::div(phi_), epsilon_)
228  - fvm::laplacian(DepsilonEff(), epsilon_)
229  ==
230  C1_*G*epsilon_/k_
231  - fvm::Sp(C2_*epsilon_/k_, epsilon_)
232  );
233 
234  epsEqn().relax();
235 
236  epsEqn().boundaryManipulate(epsilon_.boundaryField());
237 
238  solve(epsEqn);
239  bound(epsilon_, epsilon0_);
240 
241 
242  // Turbulent kinetic energy equation
244  (
245  fvm::ddt(k_)
246  + fvm::div(phi_, k_)
247  - fvm::Sp(fvc::div(phi_), k_)
248  - fvm::laplacian(DkEff(), k_)
249  ==
250  G
251  - fvm::Sp(epsilon_/k_, k_)
252  );
253 
254  kEqn().relax();
255  solve(kEqn);
256  bound(k_, k0_);
257 
258 
259  // Re-calculate viscosity
260  nut_ = Cmu_*sqr(k_)/epsilon_;
262 }
263 
264 
265 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
266 
267 } // End namespace RASModels
268 } // End namespace incompressible
269 } // End namespace Foam
270 
271 // ************************ vim: set sw=4 sts=4 et: ************************ //