FreeFOAM The Cross-Platform CFD Toolkit
spectEddyVisc.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 "spectEddyVisc.H"
28 
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace incompressible
34 {
35 namespace LESModels
36 {
37 
38 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
39 
40 defineTypeNameAndDebug(spectEddyVisc, 0);
41 addToRunTimeSelectionTable(LESModel, spectEddyVisc, dictionary);
42 
43 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
44 
45 void spectEddyVisc::updateSubGridScaleFields(const volTensorField& gradU)
46 {
47  volScalarField Re = sqr(delta())*mag(symm(gradU))/nu();
48  for (label i=0; i<5; i++)
49  {
50  nuSgs_ =
51  nu()
52  /(
53  scalar(1)
54  - exp(-cB_*pow(nu()/(nuSgs_ + nu()), 1.0/3.0)*pow(Re, -2.0/3.0))
55  );
56  }
57 
59 }
60 
61 
62 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
63 
64 spectEddyVisc::spectEddyVisc
65 (
66  const volVectorField& U,
67  const surfaceScalarField& phi,
68  transportModel& transport
69 )
70 :
71  LESModel(typeName, U, phi, transport),
72  GenEddyVisc(U, phi, transport),
73 
74 
75  cB_
76  (
78  (
79  "cB",
80  coeffDict_,
81  8.22
82  )
83  ),
84  cK1_
85  (
87  (
88  "cK1",
89  coeffDict_,
90  0.83
91  )
92  ),
93  cK2_
94  (
96  (
97  "cK2",
98  coeffDict_,
99  1.03
100  )
101  ),
102  cK3_
103  (
105  (
106  "cK3",
107  coeffDict_,
108  4.75
109  )
110  ),
111  cK4_
112  (
114  (
115  "cK4",
116  coeffDict_,
117  2.55
118  )
119  )
120 {
121  printCoeffs();
122 
123  updateSubGridScaleFields(fvc::grad(U));
124 }
125 
126 
127 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
128 
130 {
131  volScalarField eps = 2*nuEff()*magSqr(symm(fvc::grad(U())));
132 
133  return
134  cK1_*pow(delta()*eps, 2.0/3.0)
135  *exp(-cK2_*pow(delta(), -4.0/3.0)*nu()/pow(eps, 1.0/3.0))
136  - cK3_*sqrt(eps*nu())
137  *erfc(cK4_*pow(delta(), -2.0/3.0)*sqrt(nu())*pow(eps, -1.0/6.0));
138 }
139 
140 
142 {
143  GenEddyVisc::correct(gradU);
144  updateSubGridScaleFields(gradU());
145 }
146 
147 
149 {
150  if (GenEddyVisc::read())
151  {
152  cB_.readIfPresent(coeffDict());
153  cK1_.readIfPresent(coeffDict());
154  cK2_.readIfPresent(coeffDict());
155  cK3_.readIfPresent(coeffDict());
156  cK4_.readIfPresent(coeffDict());
157 
158  return true;
159  }
160  else
161  {
162  return false;
163  }
164 }
165 
166 
167 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
168 
169 } // End namespace LESModels
170 } // End namespace incompressible
171 } // End namespace Foam
172 
173 // ************************ vim: set sw=4 sts=4 et: ************************ //