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 compressible
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  muSgs_ = ck_(D)*rho()*sqrt(k_)*delta();
48  muSgs_.correctBoundaryConditions();
49 
50  alphaSgs_ = muSgs_/Prt_;
51  alphaSgs_.correctBoundaryConditions();
52 }
53 
54 
55 dimensionedScalar dynOneEqEddy::ck_(const volSymmTensorField& D) const
56 {
57  volScalarField KK = 0.5*(filter_(magSqr(U())) - magSqr(filter_(U())));
58 
59  volSymmTensorField LL = dev(filter_(sqr(U())) - (sqr(filter_(U()))));
60 
62  delta()*(filter_(sqrt(k_)*D) - 2*sqrt(KK + filter_(k_))*filter_(D));
63 
64  return average(LL && MM)/average(magSqr(MM));
65 }
66 
67 
68 dimensionedScalar dynOneEqEddy::ce_(const volSymmTensorField& D) const
69 {
70  volScalarField KK = 0.5*(filter_(magSqr(U())) - magSqr(filter_(U())));
71 
72  volScalarField mm =
73  pow(KK + filter_(k_), 1.5)/(2*delta()) - filter_(pow(k_, 1.5))/delta();
74 
75  volScalarField ee =
76  2*delta()*ck_(D)*
77  (
78  filter_(sqrt(k_)*magSqr(D))
79  - 2*sqrt(KK + filter_(k_))*magSqr(filter_(D))
80  );
81 
82  return average(ee*mm)/average(mm*mm);
83 }
84 
85 
86 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
87 
88 dynOneEqEddy::dynOneEqEddy
89 (
90  const volScalarField& rho,
91  const volVectorField& U,
92  const surfaceScalarField& phi,
93  const basicThermo& thermoPhysicalModel
94 )
95 :
96  LESModel(typeName, rho, U, phi, thermoPhysicalModel),
97  GenEddyVisc(rho, U, phi, thermoPhysicalModel),
98 
99  filterPtr_(LESfilter::New(U.mesh(), coeffDict())),
100  filter_(filterPtr_())
101 {
102  updateSubGridScaleFields(dev(symm(fvc::grad(U))));
103 
104  printCoeffs();
105 }
106 
107 
108 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
109 
111 {
112  const volTensorField& gradU = tgradU();
113 
114  GenEddyVisc::correct(gradU);
115 
116  volSymmTensorField D = dev(symm(gradU));
118  volScalarField G = 2*muSgs_*(gradU && D);
119 
120  solve
121  (
122  fvm::ddt(rho(), k_)
123  + fvm::div(phi(), k_)
124  - fvm::laplacian(DkEff(), k_)
125  ==
126  G
127  - fvm::SuSp(2.0/3.0*rho()*divU, k_)
128  - fvm::Sp(ce_(D)*rho()*sqrt(k_)/delta(), k_)
129  );
130 
131  bound(k_, dimensionedScalar("0", k_.dimensions(), 1.0e-10));
132 
133  updateSubGridScaleFields(D);
134 }
135 
136 
137 bool dynOneEqEddy::read()
138 {
139  if (GenEddyVisc::read())
140  {
141  filter_.read(coeffDict());
142 
143  return true;
144  }
145  else
146  {
147  return false;
148  }
149 }
150 
151 
152 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
153 
154 } // End namespace LESModels
155 } // End namespace compressible
156 } // End namespace Foam
157 
158 // ************************ vim: set sw=4 sts=4 et: ************************ //