FreeFOAM The Cross-Platform CFD Toolkit
kOmega.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 "kOmega.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(kOmega, 0);
43 addToRunTimeSelectionTable(RASModel, kOmega, 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  "betaStar",
61  coeffDict_,
62  0.09
63  )
64  ),
65  beta_
66  (
68  (
69  "beta",
70  coeffDict_,
71  0.072
72  )
73  ),
74  alpha_
75  (
77  (
78  "alpha",
79  coeffDict_,
80  0.52
81  )
82  ),
83  alphaK_
84  (
86  (
87  "alphaK",
88  coeffDict_,
89  0.5
90  )
91  ),
92  alphaOmega_
93  (
95  (
96  "alphaOmega",
97  coeffDict_,
98  0.5
99  )
100  ),
101 
102  k_
103  (
104  IOobject
105  (
106  "k",
107  runTime_.timeName(),
108  mesh_,
111  ),
112  autoCreateK("k", mesh_)
113  ),
114  omega_
115  (
116  IOobject
117  (
118  "omega",
119  runTime_.timeName(),
120  mesh_,
123  ),
124  autoCreateOmega("omega", mesh_)
125  ),
126  nut_
127  (
128  IOobject
129  (
130  "nut",
131  runTime_.timeName(),
132  mesh_,
135  ),
136  autoCreateNut("nut", mesh_)
137  )
138 {
139  nut_ = k_/(omega_ + omegaSmall_);
140  nut_.correctBoundaryConditions();
141 
142  printCoeffs();
143 }
144 
145 
146 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
147 
149 {
151  (
153  (
154  IOobject
155  (
156  "R",
157  runTime_.timeName(),
158  mesh_,
161  ),
162  ((2.0/3.0)*I)*k_ - nut_*twoSymm(fvc::grad(U_)),
163  k_.boundaryField().types()
164  )
165  );
166 }
167 
168 
170 {
172  (
174  (
175  IOobject
176  (
177  "devRhoReff",
178  runTime_.timeName(),
179  mesh_,
182  ),
184  )
185  );
186 }
187 
188 
190 {
191  return
192  (
193  - fvm::laplacian(nuEff(), U)
194  - fvc::div(nuEff()*dev(fvc::grad(U)().T()))
195  );
196 }
197 
198 
200 {
201  if (RASModel::read())
202  {
203  Cmu_.readIfPresent(coeffDict());
204  beta_.readIfPresent(coeffDict());
205  alphaK_.readIfPresent(coeffDict());
206  alphaOmega_.readIfPresent(coeffDict());
207 
208  return true;
209  }
210  else
211  {
212  return false;
213  }
214 }
215 
216 
218 {
220 
221  if (!turbulence_)
222  {
223  return;
224  }
225 
226  volScalarField G("RASModel::G", nut_*2*magSqr(symm(fvc::grad(U_))));
227 
228  // Update omega and G at the wall
229  omega_.boundaryField().updateCoeffs();
230 
231  // Turbulence specific dissipation rate equation
232  tmp<fvScalarMatrix> omegaEqn
233  (
234  fvm::ddt(omega_)
235  + fvm::div(phi_, omega_)
236  - fvm::Sp(fvc::div(phi_), omega_)
237  - fvm::laplacian(DomegaEff(), omega_)
238  ==
239  alpha_*G*omega_/k_
240  - fvm::Sp(beta_*omega_, omega_)
241  );
242 
243  omegaEqn().relax();
244 
245  omegaEqn().boundaryManipulate(omega_.boundaryField());
246 
247  solve(omegaEqn);
248  bound(omega_, omega0_);
249 
250 
251  // Turbulent kinetic energy equation
253  (
254  fvm::ddt(k_)
255  + fvm::div(phi_, k_)
256  - fvm::Sp(fvc::div(phi_), k_)
257  - fvm::laplacian(DkEff(), k_)
258  ==
259  G
260  - fvm::Sp(Cmu_*omega_, k_)
261  );
262 
263  kEqn().relax();
264  solve(kEqn);
265  bound(k_, k0_);
266 
267 
268  // Re-calculate viscosity
269  nut_ = k_/(omega_ + omegaSmall_);
271 }
272 
273 
274 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
275 
276 } // End namespace RASModels
277 } // End namespace incompressible
278 } // End namespace Foam
279 
280 // ************************ vim: set sw=4 sts=4 et: ************************ //