1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
5 \\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
7 -------------------------------------------------------------------------------
9 This file is part of OpenFOAM.
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 the
13 Free Software Foundation; either version 2 of the License, or (at your
14 option) any later version.
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
21 You should have received a copy of the GNU General Public License
22 along with OpenFOAM; if not, write to the Free Software Foundation,
23 Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
25 \*---------------------------------------------------------------------------*/
27 #include "lowReOneEqEddy.H"
28 #include "addToRunTimeSelectionTable.H"
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34 namespace compressible
39 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
41 defineTypeNameAndDebug(lowReOneEqEddy, 0);
42 addToRunTimeSelectionTable(LESModel, lowReOneEqEddy, dictionary);
44 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
46 void lowReOneEqEddy::updateSubGridScaleFields()
48 // High Re eddy viscosity
49 muSgs_ = ck_*rho()*sqrt(k_)*delta();
51 // low Re no corrected eddy viscosity
52 muSgs_ -= (mu()/beta_)*(scalar(1) - exp(-beta_*muSgs_/mu()));
53 muSgs_.correctBoundaryConditions();
55 alphaSgs_ = muSgs_/Prt_;
56 alphaSgs_.correctBoundaryConditions();
60 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
62 lowReOneEqEddy::lowReOneEqEddy
64 const volScalarField& rho,
65 const volVectorField& U,
66 const surfaceScalarField& phi,
67 const basicThermo& thermoPhysicalModel
70 LESModel(typeName, rho, U, phi, thermoPhysicalModel),
71 GenEddyVisc(rho, U, phi, thermoPhysicalModel),
75 dimensioned<scalar>::lookupOrAddToDict
84 dimensioned<scalar>::lookupOrAddToDict
92 updateSubGridScaleFields();
98 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
100 void lowReOneEqEddy::correct(const tmp<volTensorField>& tgradU)
102 const volTensorField& gradU = tgradU();
104 GenEddyVisc::correct(gradU);
106 volScalarField divU = fvc::div(phi()/fvc::interpolate(rho()));
107 volScalarField G = 2*muSgs_*(gradU && dev(symm(gradU)));
112 + fvm::div(phi(), k_)
113 - fvm::laplacian(DkEff(), k_)
116 - fvm::SuSp(2.0/3.0*rho()*divU, k_)
117 - fvm::Sp(ce_*rho()*sqrt(k_)/delta(), k_)
122 updateSubGridScaleFields();
126 bool lowReOneEqEddy::read()
128 if (GenEddyVisc::read())
130 ck_.readIfPresent(coeffDict());
131 beta_.readIfPresent(coeffDict());
142 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
144 } // End namespace LESModels
145 } // End namespace compressible
146 } // End namespace Foam
148 // ************************************************************************* //