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