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 "vanDriestDelta.H"
28 #include <compressibleLESModels/LESModel.H>
29 #include <finiteVolume/wallFvPatch.H>
30 #include <finiteVolume/wallDistData.H>
31 #include <finiteVolume/wallPointYPlus.H>
32 #include <OpenFOAM/addToRunTimeSelectionTable.H>
34 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
38 namespace compressible
43 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
45 defineTypeNameAndDebug(vanDriestDelta, 0);
46 addToRunTimeSelectionTable(LESdelta, vanDriestDelta, dictionary);
48 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
50 void vanDriestDelta::calcDelta()
52 const LESModel& lesModel = mesh_.lookupObject<LESModel>("LESProperties");
54 const volVectorField& U = lesModel.U();
55 const volScalarField& rho = lesModel.rho();
56 const volScalarField& mu = lesModel.mu();
57 tmp<volScalarField> muSgs = lesModel.muSgs();
64 mesh_.time().constant(),
68 dimensionedScalar("ystar", dimLength, GREAT)
71 const fvPatchList& patches = mesh_.boundary();
72 forAll(patches, patchi)
74 if (isType<wallFvPatch>(patches[patchi]))
76 const fvPatchVectorField& Uw = U.boundaryField()[patchi];
77 const scalarField& rhow = rho.boundaryField()[patchi];
78 const scalarField& muw = mu.boundaryField()[patchi];
79 const scalarField& muSgsw = muSgs().boundaryField()[patchi];
81 ystar.boundaryField()[patchi] =
82 muw/(rhow*sqrt((muw + muSgsw)*mag(Uw.snGrad())/rhow + VSMALL));
86 wallPointYPlus::yPlusCutOff = 500;
87 wallDistData<wallPointYPlus> y(mesh_, ystar);
91 static_cast<const volScalarField&>(geometricDelta_()),
92 (kappa_/Cdelta_)*((scalar(1) + SMALL) - exp(-y/ystar/Aplus_))*y
97 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
99 vanDriestDelta::vanDriestDelta
106 LESdelta(name, mesh),
109 LESdelta::New("geometricDelta", mesh, dd.subDict(type() + "Coeffs"))
111 kappa_(dd.lookupOrDefault<scalar>("kappa", 0.41)),
114 dd.subDict(type() + "Coeffs").lookupOrDefault<scalar>("Aplus", 26.0)
118 dd.subDict(type() + "Coeffs").lookupOrDefault<scalar>("Cdelta", 0.158)
122 dd.subDict(type() + "Coeffs").lookupOrDefault<label>("calcInterval", 1)
125 delta_ = geometricDelta_();
129 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
131 void vanDriestDelta::read(const dictionary& d)
133 const dictionary& dd(d.subDict(type() + "Coeffs"));
135 geometricDelta_().read(dd);
136 d.readIfPresent<scalar>("kappa", kappa_);
137 dd.readIfPresent<scalar>("Aplus", Aplus_);
138 dd.readIfPresent<scalar>("Cdelta", Cdelta_);
139 dd.readIfPresent<label>("calcInterval", calcInterval_);
144 void vanDriestDelta::correct()
146 if (mesh().time().timeIndex() % calcInterval_ == 0)
148 geometricDelta_().correct();
154 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
156 } // End namespace LESModels
157 } // End namespace compressible
158 } // End namespace Foam
160 // ************************ vim: set sw=4 sts=4 et: ************************ //