Forward compatibility: flex
[foam-extend-3.2.git] / src / turbulenceModels / incompressible / LES / locDynOneEqEddy / locDynOneEqEddy.C
blobd6deafc6ea109427051d51f04a6835875a43c79e
1 /*---------------------------------------------------------------------------*\
2   =========                 |
3   \\      /  F ield         | foam-extend: Open Source CFD
4    \\    /   O peration     | Version:     3.2
5     \\  /    A nd           | Web:         http://www.foam-extend.org
6      \\/     M anipulation  | For copyright notice see file Copyright
7 -------------------------------------------------------------------------------
8 License
9     This file is part of foam-extend.
11     foam-extend is free software: you can redistribute it and/or modify it
12     under the terms of the GNU General Public License as published by the
13     Free Software Foundation, either version 3 of the License, or (at your
14     option) any later version.
16     foam-extend is distributed in the hope that it will be useful, but
17     WITHOUT ANY WARRANTY; without even the implied warranty of
18     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
19     General Public License for more details.
21     You should have received a copy of the GNU General Public License
22     along with foam-extend.  If not, see <http://www.gnu.org/licenses/>.
24 \*---------------------------------------------------------------------------*/
26 #include "locDynOneEqEddy.H"
27 #include "addToRunTimeSelectionTable.H"
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 namespace Foam
33 namespace incompressible
35 namespace LESModels
38 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
40 defineTypeNameAndDebug(locDynOneEqEddy, 0);
41 addToRunTimeSelectionTable(LESModel, locDynOneEqEddy, dictionary);
43 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
45 void locDynOneEqEddy::updateSubGridScaleFields
47     const volSymmTensorField& D,
48     const volScalarField& KK
51     nuSgs_ = ck(D, KK)*sqrt(k_)*delta();
52     nuSgs_.correctBoundaryConditions();
56 volScalarField locDynOneEqEddy::ck
58     const volSymmTensorField& D,
59     const volScalarField& KK
60 ) const
62     volSymmTensorField LL =
63         simpleFilter_(dev(filter_(sqr(U())) - (sqr(filter_(U())))));
65     volSymmTensorField MM = simpleFilter_(-2.0*delta()*pow(KK, 0.5)*filter_(D));
67     volScalarField ck =
68         simpleFilter_(0.5*(LL && MM))
69        /(
70             simpleFilter_(magSqr(MM))
71           + dimensionedScalar("small", sqr(MM.dimensions()), VSMALL)
72         );
74     return 0.5*(mag(ck) + ck);
78 volScalarField locDynOneEqEddy::ce
80     const volSymmTensorField& D,
81     const volScalarField& KK
82 ) const
84     volScalarField ce =
85         simpleFilter_(nuEff()*(filter_(magSqr(D)) - magSqr(filter_(D))))
86        /simpleFilter_(pow(KK, 1.5)/(2.0*delta()));
88     return 0.5*(mag(ce) + ce);
92 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
94 locDynOneEqEddy::locDynOneEqEddy
96     const volVectorField& U,
97     const surfaceScalarField& phi,
98     transportModel& transport
101     LESModel(typeName, U, phi, transport),
102     GenEddyVisc(U, phi, transport),
104     k_
105     (
106         IOobject
107         (
108             "k",
109             runTime_.timeName(),
110             U_.db(),
111             IOobject::MUST_READ,
112             IOobject::AUTO_WRITE
113         ),
114         mesh_
115     ),
117     simpleFilter_(U.mesh()),
118     filterPtr_(LESfilter::New(U.mesh(), coeffDict())),
119     filter_(filterPtr_())
121     volScalarField KK = 0.5*(filter_(magSqr(U)) - magSqr(filter_(U)));
122     updateSubGridScaleFields(symm(fvc::grad(U)), KK);
124     printCoeffs();
128 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
130 void locDynOneEqEddy::correct(const tmp<volTensorField>& gradU)
132     LESModel::correct(gradU);
134     volSymmTensorField D = symm(gradU);
136     volScalarField KK = 0.5*(filter_(magSqr(U())) - magSqr(filter_(U())));
137     KK.max(dimensionedScalar("small", KK.dimensions(), SMALL));
139     volScalarField P = 2.0*nuSgs_*magSqr(D);
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     );
151     kEqn.relax();
152     kEqn.solve();
154     bound(k_, k0());
156     updateSubGridScaleFields(D, KK);
160 bool locDynOneEqEddy::read()
162     if (GenEddyVisc::read())
163     {
164         filter_.read(coeffDict());
166         return true;
167     }
168     else
169     {
170         return false;
171     }
175 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
177 } // End namespace LESModels
178 } // End namespace incompressible
179 } // End namespace Foam
181 // ************************************************************************* //