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 "LienCubicKELowRe.H"
28 #include <finiteVolume/wallFvPatch.H>
29 #include <OpenFOAM/addToRunTimeSelectionTable.H>
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35 namespace incompressible
40 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
42 defineTypeNameAndDebug(LienCubicKELowRe, 0);
43 addToRunTimeSelectionTable(RASModel, LienCubicKELowRe, dictionary);
45 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
47 LienCubicKELowRe::LienCubicKELowRe
49 const volVectorField& U,
50 const surfaceScalarField& phi,
51 transportModel& lamTransportModel
54 RASModel(typeName, U, phi, lamTransportModel),
58 dimensioned<scalar>::lookupOrAddToDict
67 dimensioned<scalar>::lookupOrAddToDict
76 dimensioned<scalar>::lookupOrAddToDict
85 dimensioned<scalar>::lookupOrAddToDict
94 dimensioned<scalar>::lookupOrAddToDict
103 dimensioned<scalar>::lookupOrAddToDict
112 dimensioned<scalar>::lookupOrAddToDict
121 dimensioned<scalar>::lookupOrAddToDict
130 dimensioned<scalar>::lookupOrAddToDict
139 dimensioned<scalar>::lookupOrAddToDict
148 dimensioned<scalar>::lookupOrAddToDict
157 dimensioned<scalar>::lookupOrAddToDict
166 dimensioned<scalar>::lookupOrAddToDict
175 dimensioned<scalar>::lookupOrAddToDict
184 dimensioned<scalar>::lookupOrAddToDict
220 gradU_(fvc::grad(U)),
221 eta_(k_/epsilon_*sqrt(2.0*magSqr(0.5*(gradU_ + gradU_.T())))),
222 ksi_(k_/epsilon_*sqrt(2.0*magSqr(0.5*(gradU_ - gradU_.T())))),
223 Cmu_(2.0/(3.0*(A1_ + eta_ + alphaKsi_*ksi_))),
224 fEta_(A2_ + pow(eta_, 3.0)),
228 -2.0*pow(Cmu_, 3.0)*pow(k_, 4.0)/pow(epsilon_, 3.0)
229 *(magSqr(gradU_ + gradU_.T()) - magSqr(gradU_ - gradU_.T()))
232 yStar_(sqrt(k_)*y_/nu() + SMALL),
238 scalar(1) - exp(-Am_*yStar_))
239 /(scalar(1) - exp(-Aepsilon_*yStar_) + SMALL
241 *sqr(k_)/(epsilon_ + epsilonSmall_)
242 // cubic term C5, implicit part
246 dimensionedScalar("0", C5viscosity_.dimensions(), 0.0)
249 // turbulent viscosity, with implicit part of C5
257 pow(k_, 3.0)/sqr(epsilon_)
262 + (gradU_ & gradU_)().T()
264 + Ctau2_/fEta_*(gradU_ & gradU_.T())
265 + Ctau3_/fEta_*(gradU_.T() & gradU_)
268 - 20.0*pow(k_, 4.0)/pow(epsilon_, 3.0)
271 ((gradU_ & gradU_) & gradU_.T())
272 + ((gradU_ & gradU_.T()) & gradU_.T())
273 - ((gradU_.T() & gradU_) & gradU_)
274 - ((gradU_.T() & gradU_.T()) & gradU_)
276 // cubic term C5, explicit part
280 dimensionedScalar("0", C5viscosity_.dimensions(), 0.0)
289 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
291 tmp<volSymmTensorField> LienCubicKELowRe::R() const
293 return tmp<volSymmTensorField>
295 new volSymmTensorField
305 ((2.0/3.0)*I)*k_ - nut_*twoSymm(gradU_) + nonlinearStress_,
306 k_.boundaryField().types()
312 tmp<volSymmTensorField> LienCubicKELowRe::devReff() const
314 return tmp<volSymmTensorField>
316 new volSymmTensorField
326 -nuEff()*dev(twoSymm(fvc::grad(U_))) + nonlinearStress_
332 tmp<fvVectorMatrix> LienCubicKELowRe::divDevReff(volVectorField& U) const
336 fvc::div(nonlinearStress_)
337 - fvm::laplacian(nuEff(), U)
338 - fvc::div(nuEff()*dev(fvc::grad(U)().T()))
343 bool LienCubicKELowRe::read()
345 if (RASModel::read())
347 C1_.readIfPresent(coeffDict());
348 C2_.readIfPresent(coeffDict());
349 sigmak_.readIfPresent(coeffDict());
350 sigmaEps_.readIfPresent(coeffDict());
351 A1_.readIfPresent(coeffDict());
352 A2_.readIfPresent(coeffDict());
353 Ctau1_.readIfPresent(coeffDict());
354 Ctau2_.readIfPresent(coeffDict());
355 Ctau3_.readIfPresent(coeffDict());
356 alphaKsi_.readIfPresent(coeffDict());
357 CmuWall_.readIfPresent(coeffDict());
358 kappa_.readIfPresent(coeffDict());
359 Am_.readIfPresent(coeffDict());
360 Aepsilon_.readIfPresent(coeffDict());
361 Amu_.readIfPresent(coeffDict());
372 void LienCubicKELowRe::correct()
381 if (mesh_.changing())
386 gradU_ = fvc::grad(U_);
389 volScalarField S2 = symm(gradU_) && gradU_;
391 yStar_ = sqrt(k_)*y_/nu() + SMALL;
392 volScalarField Rt = sqr(k_)/(nu()*epsilon_);
395 (scalar(1) - exp(-Am_*yStar_))
396 /(scalar(1) - exp(-Aepsilon_*yStar_) + SMALL);
398 volScalarField f2 = scalar(1) - 0.3*exp(-sqr(Rt));
401 Cmu_*fMu*sqr(k_)/epsilon_*S2 - (nonlinearStress_ && gradU_);
403 // Dissipation equation
404 tmp<fvScalarMatrix> epsEqn
407 + fvm::div(phi_, epsilon_)
408 - fvm::laplacian(DepsilonEff(), epsilon_)
412 + C2_*f2*pow(Cmu_, 0.75)*pow(k_, scalar(0.5))
413 /(kappa_*y_*(scalar(1) - exp(-Aepsilon_*yStar_)))
414 *exp(-Amu_*sqr(yStar_))*epsilon_
415 - fvm::Sp(C2_*f2*epsilon_/k_, epsilon_)
420 # include "LienCubicKELowReSetWallDissipation.H"
421 # include <incompressibleRASModels/wallDissipationI.H>
424 bound(epsilon_, epsilon0_);
427 // Turbulent kinetic energy equation
429 tmp<fvScalarMatrix> kEqn
433 - fvm::laplacian(DkEff(), k_)
436 - fvm::Sp(epsilon_/k_, k_)
444 // Re-calculate viscosity
446 eta_ = k_/epsilon_*sqrt(2.0*magSqr(0.5*(gradU_ + gradU_.T())));
447 ksi_ = k_/epsilon_*sqrt(2.0*magSqr(0.5*(gradU_ - gradU_.T())));
448 Cmu_ = 2.0/(3.0*(A1_ + eta_ + alphaKsi_*ksi_));
449 fEta_ = A2_ + pow(eta_, 3.0);
452 - 2.0*pow(Cmu_, 3.0)*pow(k_, 4.0)/pow(epsilon_, 3.0)
453 *(magSqr(gradU_ + gradU_.T()) - magSqr(gradU_ - gradU_.T()));
456 Cmu_*fMu*sqr(k_)/epsilon_
461 dimensionedScalar("0", C5viscosity_.dimensions(), 0.0)
464 nonlinearStress_ = symm
467 pow(k_, 3.0)/sqr(epsilon_)
472 + (gradU_ & gradU_)().T()
474 + Ctau2_/fEta_*(gradU_ & gradU_.T())
475 + Ctau3_/fEta_*(gradU_.T() & gradU_)
478 - 20.0*pow(k_, 4.0)/pow(epsilon_, 3.0)
481 ((gradU_ & gradU_) & gradU_.T())
482 + ((gradU_ & gradU_.T()) & gradU_.T())
483 - ((gradU_.T() & gradU_) & gradU_)
484 - ((gradU_.T() & gradU_.T()) & gradU_)
486 // cubic term C5, explicit part
490 dimensionedScalar("0", C5viscosity_.dimensions(), 0.0)
496 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
498 } // End namespace RASModels
499 } // End namespace incompressible
500 } // End namespace Foam
502 // ************************ vim: set sw=4 sts=4 et: ************************ //