FreeFOAM The Cross-Platform CFD Toolkit
kOmegaSSTSAS.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) 2008-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 "kOmegaSSTSAS.H"
28 #include <finiteVolume/wallDist.H>
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace incompressible
35 {
36 namespace LESModels
37 {
38 
39 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
40 
41 defineTypeNameAndDebug(kOmegaSSTSAS, 0);
42 addToRunTimeSelectionTable(LESModel, kOmegaSSTSAS, dictionary);
43 
44 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
45 
46 void kOmegaSSTSAS::updateSubGridScaleFields(const volScalarField& S2)
47 {
48  nuSgs_ == a1_*k_/max(a1_*omega_, F2()*sqrt(S2));
50 }
51 
52 
54 {
55  volScalarField CDkOmegaPlus = max
56  (
57  CDkOmega,
58  dimensionedScalar("1.0e-10", dimless/sqr(dimTime), 1.0e-10)
59  );
60 
61  volScalarField arg1 = min
62  (
63  min
64  (
65  max
66  (
67  (scalar(1)/betaStar_)*sqrt(k_)/(omega_*y_),
68  scalar(500)*nu()/(sqr(y_)*omega_)
69  ),
70  (4*alphaOmega2_)*k_/(CDkOmegaPlus*sqr(y_))
71  ),
72  scalar(10)
73  );
74 
75  return tanh(pow4(arg1));
76 }
77 
78 
80 {
81  volScalarField arg2 = min
82  (
83  max
84  (
85  (scalar(2)/betaStar_)*sqrt(k_)/(omega_*y_),
86  scalar(500)*nu()/(sqr(y_)*omega_)
87  ),
88  scalar(100)
89  );
90 
91  return tanh(sqr(arg2));
92 }
93 
94 
96 (
97  const volScalarField& S2
98 ) const
99 {
100  return max
101  (
102  kappa_*sqrt(S2)
103  /(
104  mag(fvc::laplacian(U()))
106  (
107  "ROOTVSMALL",
108  dimensionSet(0, -1 , -1, 0, 0, 0, 0),
109  ROOTVSMALL
110  )
111  ),
112  Cs_*delta()
113  );
114 }
115 
116 
117 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
118 
119 kOmegaSSTSAS::kOmegaSSTSAS
120 (
121  const volVectorField& U,
122  const surfaceScalarField& phi,
123  transportModel& transport,
124  const word& modelName
125 )
126 :
127  LESModel(modelName, U, phi, transport),
128 
129  alphaK1_
130  (
132  (
133  "alphaK1",
134  coeffDict_,
135  0.85034
136  )
137  ),
138  alphaK2_
139  (
141  (
142  "alphaK2",
143  coeffDict_,
144  1.0
145  )
146  ),
147  alphaOmega1_
148  (
150  (
151  "alphaOmega1",
152  coeffDict_,
153  0.5
154  )
155  ),
156  alphaOmega2_
157  (
159  (
160  "alphaOmega2",
161  coeffDict_,
162  0.85616
163  )
164  ),
165  gamma1_
166  (
168  (
169  "gamma1",
170  coeffDict_,
171  0.5532
172  )
173  ),
174  gamma2_
175  (
177  (
178  "gamma2",
179  coeffDict_,
180  0.4403
181  )
182  ),
183  beta1_
184  (
186  (
187  "beta1",
188  coeffDict_,
189  0.075
190  )
191  ),
192  beta2_
193  (
195  (
196  "beta2",
197  coeffDict_,
198  0.0828
199  )
200  ),
201  betaStar_
202  (
204  (
205  "betaStar",
206  coeffDict_,
207  0.09
208  )
209  ),
210  a1_
211  (
213  (
214  "a1",
215  coeffDict_,
216  0.31
217  )
218  ),
219  c1_
220  (
222  (
223  "c1",
224  coeffDict_,
225  10.0
226  )
227  ),
228  Cs_
229  (
231  (
232  "Cs",
233  coeffDict_,
234  0.262
235  )
236  ),
237  alphaPhi_
238  (
240  (
241  "alphaPhi",
242  coeffDict_,
243  0.666667
244  )
245  ),
246  zetaTilda2_
247  (
249  (
250  "zetaTilda2",
251  coeffDict_,
252  1.755
253  )
254  ),
255  FSAS_
256  (
258  (
259  "FSAS",
260  coeffDict_,
261  1.25
262  )
263  ),
264 
265  omega0_("omega0", dimless/dimTime, SMALL),
266  omegaSmall_("omegaSmall", dimless/dimTime, SMALL),
267  y_(mesh_),
268  Cmu_
269  (
271  (
272  "Cmu",
273  coeffDict_,
274  0.09
275  )
276  ),
277  kappa_
278  (
280  (
281  "kappa",
282  *this,
283  0.41
284  )
285  ),
286 
287  k_
288  (
289  IOobject
290  (
291  "k",
292  runTime_.timeName(),
293  mesh_,
296  ),
297  mesh_
298  ),
299 
300  omega_
301  (
302  IOobject
303  (
304  "omega",
305  runTime_.timeName(),
306  mesh_,
309  ),
310  mesh_
311  ),
312 
313  nuSgs_
314  (
315  IOobject
316  (
317  "nuSgs",
318  runTime_.timeName(),
319  mesh_,
322  ),
323  mesh_
324  )
325 {
326  updateSubGridScaleFields(magSqr(2.0*symm(fvc::grad(U))));
327 
328  printCoeffs();
329 }
330 
331 
332 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
333 
335 {
336  LESModel::correct(gradU);
337 
338  if (mesh_.changing())
339  {
340  y_.correct();
341  }
342 
343  volScalarField S2 = magSqr(2.0*symm(gradU()));
344  gradU.clear();
345 
346  volVectorField gradK = fvc::grad(k_);
347  volVectorField gradOmega = fvc::grad(omega_);
348  volScalarField L = sqrt(k_)/(pow(Cmu_, 0.25)*(omega_ + omegaSmall_));
349  volScalarField CDkOmega =
350  (2.0*alphaOmega2_)*(gradK & gradOmega)/(omega_ + omegaSmall_);
351  volScalarField F1 = this->F1(CDkOmega);
352  volScalarField G = nuSgs_*0.5*S2;
353 
354  // Turbulent kinetic energy equation
355  {
356  fvScalarMatrix kEqn
357  (
358  fvm::ddt(k_)
359  + fvm::div(phi(), k_)
360  - fvm::Sp(fvc::div(phi()), k_)
361  - fvm::laplacian(DkEff(F1), k_)
362  ==
363  min(G, c1_*betaStar_*k_*omega_)
364  - fvm::Sp(betaStar_*omega_, k_)
365  );
366 
367  kEqn.relax();
368  kEqn.solve();
369  }
370  bound(k_, k0());
371 
372  volScalarField grad_omega_k = max
373  (
374  magSqr(gradOmega)/
375  sqr(omega_ + omegaSmall_),
376  magSqr(gradK)/
377  sqr(k_ + k0())
378  );
379 
380  // Turbulent frequency equation
381  {
382  fvScalarMatrix omegaEqn
383  (
384  fvm::ddt(omega_)
385  + fvm::div(phi(), omega_)
386  - fvm::Sp(fvc::div(phi()), omega_)
387  - fvm::laplacian(DomegaEff(F1), omega_)
388  ==
389  gamma(F1)*0.5*S2
390  - fvm::Sp(beta(F1)*omega_, omega_)
391  - fvm::SuSp // cross diffusion term
392  (
393  (F1 - scalar(1))*CDkOmega/omega_,
394  omega_
395  )
396  + FSAS_
397  *max
398  (
399  dimensionedScalar("zero",dimensionSet(0, 0 , -2, 0, 0),0. ),
400  zetaTilda2_*kappa_*S2*sqr(L/Lvk2(S2))
401  - 2.0/alphaPhi_*k_*grad_omega_k
402  )
403  );
404 
405  omegaEqn.relax();
406  omegaEqn.solve();
407  }
408  bound(omega_, omega0_);
409 
410  updateSubGridScaleFields(S2);
411 }
412 
413 
415 {
416  return 2.0*nuEff()*magSqr(symm(fvc::grad(U())));
417 }
418 
419 
421 {
422  return ((2.0/3.0)*I)*k() - nuSgs()*twoSymm(fvc::grad(U()));
423 }
424 
425 
427 {
428  return -nuEff()*dev(twoSymm(fvc::grad(U())));
429 }
430 
431 
433 {
434  return
435  (
436  - fvm::laplacian(nuEff(), U) - fvc::div(nuEff()*dev(fvc::grad(U)().T()))
437  );
438 }
439 
440 
442 {
443  if (LESModel::read())
444  {
460 
461  return true;
462  }
463  else
464  {
465  return false;
466  }
467 }
468 
469 
470 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
471 
472 } // End namespace LESModels
473 } // End namespace incompressible
474 } // End namespace Foam
475 
476 // ************************ vim: set sw=4 sts=4 et: ************************ //