ENH: patchCloud: return pTraits<Type>::max for unfound points
[OpenFOAM-1.7.x.git] / src / turbulenceModels / incompressible / LES / LRRDiffStress / LRRDiffStress.C
bloba913cb2a60b5b7a05126ba825b90786ffbf1b0f7
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.
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.
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.
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/>.
24 \*---------------------------------------------------------------------------*/
26 #include "LRRDiffStress.H"
27 #include "addToRunTimeSelectionTable.H"
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 namespace Foam
33 namespace incompressible
35 namespace LESModels
38 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
40 defineTypeNameAndDebug(LRRDiffStress, 0);
41 addToRunTimeSelectionTable(LESModel, LRRDiffStress, dictionary);
43 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
45 void LRRDiffStress::updateSubGridScaleFields(const volScalarField& K)
47     nuSgs_ = ck_*sqrt(K)*delta();
48     nuSgs_.correctBoundaryConditions();
52 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
54 LRRDiffStress::LRRDiffStress
56     const volVectorField& U,
57     const surfaceScalarField& phi,
58     transportModel& transport
61     LESModel(typeName, U, phi, transport),
62     GenSGSStress(U, phi, transport),
64     ck_
65     (
66         dimensioned<scalar>::lookupOrAddToDict
67         (
68             "ck",
69             coeffDict_,
70             0.09
71         )
72     ),
73     c1_
74     (
75         dimensioned<scalar>::lookupOrAddToDict
76         (
77             "c1",
78             coeffDict_,
79             1.8
80         )
81     ),
82     c2_
83     (
84         dimensioned<scalar>::lookupOrAddToDict
85         (
86             "c2",
87             coeffDict_,
88             0.6
89         )
90     )
92     updateSubGridScaleFields(0.5*tr(B_));
94     printCoeffs();
98 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
100 void LRRDiffStress::correct(const tmp<volTensorField>& tgradU)
102     const volTensorField& gradU = tgradU();
104     GenSGSStress::correct(gradU);
106     volSymmTensorField D = symm(gradU);
108     volSymmTensorField P = -twoSymm(B_ & gradU);
110     volScalarField K = 0.5*tr(B_);
111     volScalarField Epsilon = 2*nuEff()*magSqr(D);
113     fvSymmTensorMatrix BEqn
114     (
115         fvm::ddt(B_)
116       + fvm::div(phi(), B_)
117       - fvm::laplacian(DBEff(), B_)
118       + fvm::Sp(c1_*Epsilon/K, B_)
119      ==
120         P
121       - (0.667*(1.0 - c1_)*I)*Epsilon
122       - c2_*(P - 0.333*I*tr(P))
123       - (0.667 - 2*c1_)*I*pow(K, 1.5)/delta()
124     );
126     BEqn.relax();
127     BEqn.solve();
129     // Bounding the component kinetic energies
131     forAll(B_, celli)
132     {
133         B_[celli].component(symmTensor::XX) =
134             max(B_[celli].component(symmTensor::XX), k0().value());
135         B_[celli].component(symmTensor::YY) =
136             max(B_[celli].component(symmTensor::YY), k0().value());
137         B_[celli].component(symmTensor::ZZ) =
138             max(B_[celli].component(symmTensor::ZZ), k0().value());
139     }
141     K = 0.5*tr(B_);
142     bound(K, k0());
144     updateSubGridScaleFields(K);
148 bool LRRDiffStress::read()
150     if (GenSGSStress::read())
151     {
152         ck_.readIfPresent(coeffDict());
153         c1_.readIfPresent(coeffDict());
154         c2_.readIfPresent(coeffDict());
156         return true;
157     }
158     else
159     {
160         return false;
161     }
165 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
167 } // End namespace LESModels
168 } // End namespace incompressible
169 } // End namespace Foam
171 // ************************************************************************* //