FreeFOAM The Cross-Platform CFD Toolkit
dynOneEqEddy.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 "dynOneEqEddy.H"
28 
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace incompressible
34 {
35 namespace LESModels
36 {
37 
38 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
39 
40 defineTypeNameAndDebug(dynOneEqEddy, 0);
41 addToRunTimeSelectionTable(LESModel, dynOneEqEddy, dictionary);
42 
43 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
44 
45 void dynOneEqEddy::updateSubGridScaleFields(const volSymmTensorField& D)
46 {
47  nuSgs_ = ck(D)*sqrt(k_)*delta();
49 }
50 
51 
52 dimensionedScalar dynOneEqEddy::ck(const volSymmTensorField& D) const
53 {
54  volScalarField KK = 0.5*(filter_(magSqr(U())) - magSqr(filter_(U())));
55 
56  volSymmTensorField LL = dev(filter_(sqr(U())) - sqr(filter_(U())));
57 
59  delta()*(filter_(sqrt(k_)*D) - 2*sqrt(KK + filter_(k_))*filter_(D));
60 
61  dimensionedScalar MMMM = average(magSqr(MM));
62 
63  if (MMMM.value() > VSMALL)
64  {
65  return average(LL && MM)/MMMM;
66  }
67  else
68  {
69  return 0.0;
70  }
71 }
72 
73 
74 dimensionedScalar dynOneEqEddy::ce(const volSymmTensorField& D) const
75 {
76  volScalarField KK = 0.5*(filter_(magSqr(U())) - magSqr(filter_(U())));
77 
78  volScalarField mm =
79  pow(KK + filter_(k_), 1.5)/(2*delta()) - filter_(pow(k_, 1.5))/delta();
80 
81  volScalarField ee =
82  2*delta()*ck(D)
83  *(
84  filter_(sqrt(k_)*magSqr(D))
85  - 2*sqrt(KK + filter_(k_))*magSqr(filter_(D))
86  );
87 
88  dimensionedScalar mmmm = average(magSqr(mm));
89 
90  if (mmmm.value() > VSMALL)
91  {
92  return average(ee*mm)/mmmm;
93  }
94  else
95  {
96  return 0.0;
97  }
98 }
99 
100 
101 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
102 
103 dynOneEqEddy::dynOneEqEddy
104 (
105  const volVectorField& U,
106  const surfaceScalarField& phi,
107  transportModel& transport
108 )
109 :
110  LESModel(typeName, U, phi, transport),
111  GenEddyVisc(U, phi, transport),
112 
113  k_
114  (
115  IOobject
116  (
117  "k",
118  runTime_.timeName(),
119  mesh_,
122  ),
123  mesh_
124  ),
125 
126  filterPtr_(LESfilter::New(U.mesh(), coeffDict())),
127  filter_(filterPtr_())
128 {
129  updateSubGridScaleFields(symm(fvc::grad(U)));
130 
131  printCoeffs();
132 }
133 
134 
135 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
136 
138 {
139  GenEddyVisc::correct(gradU);
140 
141  volSymmTensorField D = symm(gradU);
142 
143  volScalarField P = 2.0*nuSgs_*magSqr(D);
144 
145  fvScalarMatrix kEqn
146  (
147  fvm::ddt(k_)
148  + fvm::div(phi(), k_)
149  - fvm::laplacian(DkEff(), k_)
150  ==
151  P
152  - fvm::Sp(ce(D)*sqrt(k_)/delta(), k_)
153  );
154 
155  kEqn.relax();
156  kEqn.solve();
157 
158  bound(k_, k0());
159 
160  updateSubGridScaleFields(D);
161 }
162 
163 
165 {
166  if (GenEddyVisc::read())
167  {
168  filter_.read(coeffDict());
169 
170  return true;
171  }
172  else
173  {
174  return false;
175  }
176 }
177 
178 
179 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
180 
181 } // End namespace LESModels
182 } // End namespace incompressible
183 } // End namespace Foam
184 
185 // ************************ vim: set sw=4 sts=4 et: ************************ //