FreeFOAM The Cross-Platform CFD Toolkit
LamBremhorstKE.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 "LamBremhorstKE.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(LamBremhorstKE, 0);
43 addToRunTimeSelectionTable(RASModel, LamBremhorstKE, 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  "alphaEps",
88  coeffDict_,
89  1.3
90  )
91  ),
92 
93  k_
94  (
95  IOobject
96  (
97  "k",
98  runTime_.timeName(),
99  mesh_,
102  ),
103  mesh_
104  ),
105 
106  epsilon_
107  (
108  IOobject
109  (
110  "epsilon",
111  runTime_.timeName(),
112  mesh_,
115  ),
116  mesh_
117  ),
118 
119  y_(mesh_),
120 
121  Rt_(sqr(k_)/(nu()*epsilon_)),
122 
123  fMu_
124  (
125  sqr(scalar(1) - exp(-0.0165*(sqrt(k_)*y_/nu())))
126  *(scalar(1) + 20.5/(Rt_ + SMALL))
127  ),
128 
129  nut_
130  (
131  IOobject
132  (
133  "nut",
134  runTime_.timeName(),
135  mesh_,
138  ),
139  autoCreateLowReNut("nut", mesh_)
140  )
141 {
142  nut_ = Cmu_*fMu_*sqr(k_)/(epsilon_ + epsilonSmall_);
143  nut_.correctBoundaryConditions();
144 
145  printCoeffs();
146 }
147 
148 
149 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
150 
152 {
154  (
156  (
157  IOobject
158  (
159  "R",
160  runTime_.timeName(),
161  mesh_,
164  ),
165  ((2.0/3.0)*I)*k_ - nut_*twoSymm(fvc::grad(U_)),
166  k_.boundaryField().types()
167  )
168  );
169 }
170 
171 
173 {
175  (
177  (
178  IOobject
179  (
180  "devRhoReff",
181  runTime_.timeName(),
182  mesh_,
185  ),
187  )
188  );
189 }
190 
191 
193 {
194  return
195  (
196  - fvm::laplacian(nuEff(), U)
197  - fvc::div(nuEff()*dev(fvc::grad(U)().T()))
198  );
199 }
200 
201 
203 {
204  if (RASModel::read())
205  {
206  Cmu_.readIfPresent(coeffDict());
207  C1_.readIfPresent(coeffDict());
208  C2_.readIfPresent(coeffDict());
209  sigmaEps_.readIfPresent(coeffDict());
210 
211  return true;
212  }
213  else
214  {
215  return false;
216  }
217 }
218 
219 
221 {
223 
224  if (!turbulence_)
225  {
226  return;
227  }
228 
229  if (mesh_.changing())
230  {
231  y_.correct();
232  }
233 
234  volScalarField G("RASModel::G", nut_*2*magSqr(symm(fvc::grad(U_))));
235 
236 
237  // Calculate parameters and coefficients for low-Reynolds number model
238 
239  Rt_ = sqr(k_)/(nu()*epsilon_);
240  volScalarField Ry = sqrt(k_)*y_/nu();
241 
242  fMu_ = sqr(scalar(1) - exp(-0.0165*Ry))*(scalar(1) + 20.5/(Rt_ + SMALL));
243 
244  volScalarField f1 = scalar(1) + pow(0.05/(fMu_ + SMALL), 3);
245  volScalarField f2 = scalar(1) - exp(-sqr(Rt_));
246 
247 
248  // Dissipation equation
249 
250  tmp<fvScalarMatrix> epsEqn
251  (
252  fvm::ddt(epsilon_)
253  + fvm::div(phi_, epsilon_)
254  - fvm::laplacian(DepsilonEff(), epsilon_)
255  ==
256  C1_*f1*G*epsilon_/k_
257  - fvm::Sp(C2_*f2*epsilon_/k_, epsilon_)
258  );
259 
260  epsEqn().relax();
261  solve(epsEqn);
262  bound(epsilon_, epsilon0_);
263 
264 
265  // Turbulent kinetic energy equation
266 
268  (
269  fvm::ddt(k_)
270  + fvm::div(phi_, k_)
271  - fvm::laplacian(DkEff(), k_)
272  ==
273  G - fvm::Sp(epsilon_/k_, k_)
274  );
275 
276  kEqn().relax();
277  solve(kEqn);
278  bound(k_, k0_);
279 
280 
281  // Re-calculate viscosity
282  nut_ == Cmu_*fMu_*sqr(k_)/epsilon_;
283 }
284 
285 
286 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
287 
288 } // End namespace RASModels
289 } // End namespace incompressible
290 } // End namespace Foam
291 
292 // ************************ vim: set sw=4 sts=4 et: ************************ //