From a78139ea1cb2577134e48df6b8e078415ce313ba Mon Sep 17 00:00:00 2001 From: mattijs Date: Tue, 25 Oct 2011 15:35:13 +0100 Subject: [PATCH 01/19] COMP: wmake/rules/linux*Gcc*: clean compilation on ubuntu 11.10 --- wmake/rules/linux64Gcc/c++ | 4 ++-- wmake/rules/linux64Gcc45/c++ | 4 ++-- wmake/rules/linux64Gcc46/c++ | 4 ++-- wmake/rules/linuxGcc/c++ | 4 ++-- wmake/rules/linuxGcc45/c++ | 4 ++-- wmake/rules/linuxGcc46/c++ | 4 ++-- 6 files changed, 12 insertions(+), 12 deletions(-) diff --git a/wmake/rules/linux64Gcc/c++ b/wmake/rules/linux64Gcc/c++ index 520e39c3b1..3dfb033360 100644 --- a/wmake/rules/linux64Gcc/c++ +++ b/wmake/rules/linux64Gcc/c++ @@ -17,5 +17,5 @@ cpptoo = $(Ctoo) LINK_LIBS = $(c++DBUG) -LINKLIBSO = $(CC) $(c++FLAGS) -shared -LINKEXE = $(CC) $(c++FLAGS) -Xlinker --add-needed +LINKLIBSO = $(CC) $(c++FLAGS) -shared -Xlinker --add-needed -Xlinker --no-as-needed +LINKEXE = $(CC) $(c++FLAGS) -Xlinker --add-needed -Xlinker --no-as-needed diff --git a/wmake/rules/linux64Gcc45/c++ b/wmake/rules/linux64Gcc45/c++ index 3ca193e296..98b25ed1fe 100644 --- a/wmake/rules/linux64Gcc45/c++ +++ b/wmake/rules/linux64Gcc45/c++ @@ -17,5 +17,5 @@ cpptoo = $(Ctoo) LINK_LIBS = $(c++DBUG) -LINKLIBSO = $(CC) $(c++FLAGS) -shared -LINKEXE = $(CC) $(c++FLAGS) -Xlinker --add-needed +LINKLIBSO = $(CC) $(c++FLAGS) -shared -Xlinker --add-needed -Xlinker --no-as-needed +LINKEXE = $(CC) $(c++FLAGS) -Xlinker --add-needed -Xlinker --no-as-needed diff --git a/wmake/rules/linux64Gcc46/c++ b/wmake/rules/linux64Gcc46/c++ index 3ca193e296..98b25ed1fe 100644 --- a/wmake/rules/linux64Gcc46/c++ +++ b/wmake/rules/linux64Gcc46/c++ @@ -17,5 +17,5 @@ cpptoo = $(Ctoo) LINK_LIBS = $(c++DBUG) -LINKLIBSO = $(CC) $(c++FLAGS) -shared -LINKEXE = $(CC) $(c++FLAGS) -Xlinker --add-needed +LINKLIBSO = $(CC) $(c++FLAGS) -shared -Xlinker --add-needed -Xlinker --no-as-needed +LINKEXE = $(CC) $(c++FLAGS) -Xlinker --add-needed -Xlinker --no-as-needed diff --git a/wmake/rules/linuxGcc/c++ b/wmake/rules/linuxGcc/c++ index e862181fc5..357f4106e1 100644 --- a/wmake/rules/linuxGcc/c++ +++ b/wmake/rules/linuxGcc/c++ @@ -17,5 +17,5 @@ cpptoo = $(Ctoo) LINK_LIBS = $(c++DBUG) -LINKLIBSO = $(CC) $(c++FLAGS) -shared -LINKEXE = $(CC) $(c++FLAGS) -Xlinker --add-needed +LINKLIBSO = $(CC) $(c++FLAGS) -shared -Xlinker --add-needed -Xlinker --no-as-needed +LINKEXE = $(CC) $(c++FLAGS) -Xlinker --add-needed -Xlinker --no-as-needed diff --git a/wmake/rules/linuxGcc45/c++ b/wmake/rules/linuxGcc45/c++ index e862181fc5..357f4106e1 100644 --- a/wmake/rules/linuxGcc45/c++ +++ b/wmake/rules/linuxGcc45/c++ @@ -17,5 +17,5 @@ cpptoo = $(Ctoo) LINK_LIBS = $(c++DBUG) -LINKLIBSO = $(CC) $(c++FLAGS) -shared -LINKEXE = $(CC) $(c++FLAGS) -Xlinker --add-needed +LINKLIBSO = $(CC) $(c++FLAGS) -shared -Xlinker --add-needed -Xlinker --no-as-needed +LINKEXE = $(CC) $(c++FLAGS) -Xlinker --add-needed -Xlinker --no-as-needed diff --git a/wmake/rules/linuxGcc46/c++ b/wmake/rules/linuxGcc46/c++ index e862181fc5..357f4106e1 100644 --- a/wmake/rules/linuxGcc46/c++ +++ b/wmake/rules/linuxGcc46/c++ @@ -17,5 +17,5 @@ cpptoo = $(Ctoo) LINK_LIBS = $(c++DBUG) -LINKLIBSO = $(CC) $(c++FLAGS) -shared -LINKEXE = $(CC) $(c++FLAGS) -Xlinker --add-needed +LINKLIBSO = $(CC) $(c++FLAGS) -shared -Xlinker --add-needed -Xlinker --no-as-needed +LINKEXE = $(CC) $(c++FLAGS) -Xlinker --add-needed -Xlinker --no-as-needed From 97ee0d2cf9b32e1a017e681e547ecc17d7a57cc1 Mon Sep 17 00:00:00 2001 From: Sergio Ferraris Date: Thu, 27 Oct 2011 15:37:30 +0100 Subject: [PATCH 02/19] ENH: Adding kklOmega incompressible turbuelence model --- .../incompressible/RAS/Make/files | 1 + .../incompressible/RAS/kkLOmega/kkLOmega.C | 779 ++++++++++++++++++ .../incompressible/RAS/kkLOmega/kkLOmega.H | 296 +++++++ 3 files changed, 1076 insertions(+) create mode 100644 src/turbulenceModels/incompressible/RAS/kkLOmega/kkLOmega.C create mode 100644 src/turbulenceModels/incompressible/RAS/kkLOmega/kkLOmega.H diff --git a/src/turbulenceModels/incompressible/RAS/Make/files b/src/turbulenceModels/incompressible/RAS/Make/files index f18b124e06..e28a685903 100644 --- a/src/turbulenceModels/incompressible/RAS/Make/files +++ b/src/turbulenceModels/incompressible/RAS/Make/files @@ -16,6 +16,7 @@ LienCubicKELowRe/LienCubicKELowRe.C NonlinearKEShih/NonlinearKEShih.C LienLeschzinerLowRe/LienLeschzinerLowRe.C LamBremhorstKE/LamBremhorstKE.C +kkLOmega/kkLOmega.C /* Wall functions */ wallFunctions = derivedFvPatchFields/wallFunctions diff --git a/src/turbulenceModels/incompressible/RAS/kkLOmega/kkLOmega.C b/src/turbulenceModels/incompressible/RAS/kkLOmega/kkLOmega.C new file mode 100644 index 0000000000..ad864bb3c7 --- /dev/null +++ b/src/turbulenceModels/incompressible/RAS/kkLOmega/kkLOmega.C @@ -0,0 +1,779 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +\*---------------------------------------------------------------------------*/ + +#include "kkLOmega.H" +#include "addToRunTimeSelectionTable.H" + +#include "backwardsCompatibilityWallFunctions.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace incompressible +{ +namespace RASModels +{ + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +defineTypeNameAndDebug(kkLOmega, 0); +addToRunTimeSelectionTable(RASModel, kkLOmega, dictionary); + +// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // + +tmp kkLOmega::fv(const volScalarField& Ret) const +{ + return(1.0 - exp(-sqrt(Ret)/Av_)); +} + + +tmp kkLOmega::fINT() const +{ + return + ( + min + ( + kl_/(Cint_*(kl_ + kt_ + kMin_)), + dimensionedScalar("1.0", dimless, 1.0) + ) + ); +} + + +tmp kkLOmega::fSS(const volScalarField& omega) const +{ + return(exp(-sqr(Css_*nu()*omega/(kt_ + kMin_)))); +} + + +tmp kkLOmega::Cmu(const volScalarField& S) const +{ + return(1.0/(A0_ + As_*(S/(omega_ + omegaMin_)))); +} + + +tmp kkLOmega::BetaTS(const volScalarField& Rew) const +{ + return(scalar(1.0) - exp(-sqr(max(Rew - CtsCrit_, 0.0))/Ats_)); +} + + +tmp kkLOmega::fTaul +( + const volScalarField& lambdaEff, + const volScalarField& ktL +) const +{ + return + ( + scalar(1.0) + - exp + ( + -CtauL_*ktL + / + ( + sqr + ( + lambdaEff*omega_ + + dimensionedScalar + ( + "ROTVSMALL", + dimLength*inv(dimTime), + ROOTVSMALL + ) + ) + ) + ) + ); +} + + +tmp kkLOmega::alphaT +( + const volScalarField& lambdaEff, + const volScalarField& fv, + const volScalarField& ktS +) const +{ + return(fv*CmuStd_*sqrt(ktS)*lambdaEff); +} + + +tmp kkLOmega::fOmega +( + const volScalarField& lambdaEff, + const volScalarField& lambdaT +) const +{ + return + ( + scalar(1.0) + - exp + ( + -0.41 + * pow4 + ( + lambdaEff + / ( + lambdaT + + dimensionedScalar + ( + "ROTVSMALL", + lambdaT.dimensions(), + ROOTVSMALL + ) + ) + ) + ) + ); +} + + +tmp kkLOmega::gammaBP(const volScalarField& omega) const +{ + return + ( + max + ( + kt_/nu() + / + ( + omega + + dimensionedScalar("ROTVSMALL", omega.dimensions(), ROOTVSMALL) + ) + - + CbpCrit_ + , + 0.0 + ) + ); +} + + +tmp kkLOmega::gammaNAT +( + const volScalarField& ReOmega, + const volScalarField& fNatCrit +) const +{ + return + ( + max + ( + ReOmega + - CnatCrit_ + / ( + fNatCrit + dimensionedScalar("ROTVSMALL", dimless, ROOTVSMALL) + ) + , + 0.0 + ) + ); +} + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +kkLOmega::kkLOmega +( + const volVectorField& U, + const surfaceScalarField& phi, + transportModel& transport, + const word& turbulenceModelName, + const word& modelName +) +: + RASModel(modelName, U, phi, transport, turbulenceModelName), + + A0_ + ( + dimensioned::lookupOrAddToDict + ( + "A0", + coeffDict_, + 4.04 + ) + ), + As_ + ( + dimensioned::lookupOrAddToDict + ( + "As", + coeffDict_, + 2.12 + ) + ), + Av_ + ( + dimensioned::lookupOrAddToDict + ( + "Av", + coeffDict_, + 6.75 + ) + ), + Abp_ + ( + dimensioned::lookupOrAddToDict + ( + "Abp", + coeffDict_, + 0.6 + ) + ), + Anat_ + ( + dimensioned::lookupOrAddToDict + ( + "Anat", + coeffDict_, + 200 + ) + ), + Ats_ + ( + dimensioned::lookupOrAddToDict + ( + "Ats", + coeffDict_, + 200 + ) + ), + CbpCrit_ + ( + dimensioned::lookupOrAddToDict + ( + "CbpCrit", + coeffDict_, + 1.2 + ) + ), + Cnc_ + ( + dimensioned::lookupOrAddToDict + ( + "Cnc", + coeffDict_, + 0.1 + ) + ), + CnatCrit_ + ( + dimensioned::lookupOrAddToDict + ( + "CnatCrit", + coeffDict_, + 1250 + ) + ), + Cint_ + ( + dimensioned::lookupOrAddToDict + ( + "Cint", + coeffDict_, + 0.75 + ) + ), + CtsCrit_ + ( + dimensioned::lookupOrAddToDict + ( + "CtsCrit", + coeffDict_, + 1000 + ) + ), + CrNat_ + ( + dimensioned::lookupOrAddToDict + ( + "CrNat", + coeffDict_, + 0.02 + ) + ), + C11_ + ( + dimensioned::lookupOrAddToDict + ( + "C11", + coeffDict_, + 3.4e-6 + ) + ), + C12_ + ( + dimensioned::lookupOrAddToDict + ( + "C12", + coeffDict_, + 1.0e-10 + ) + ), + CR_ + ( + dimensioned::lookupOrAddToDict + ( + "CR", + coeffDict_, + 0.12 + ) + ), + CalphaTheta_ + ( + dimensioned::lookupOrAddToDict + ( + "CalphaTheta", + coeffDict_, + 0.035 + ) + ), + Css_ + ( + dimensioned::lookupOrAddToDict + ( + "Css", + coeffDict_, + 1.5 + ) + ), + CtauL_ + ( + dimensioned::lookupOrAddToDict + ( + "CtauL", + coeffDict_, + 4360 + ) + ), + Cw1_ + ( + dimensioned::lookupOrAddToDict + ( + "Cw1", + coeffDict_, + 0.44 + ) + ), + Cw2_ + ( + dimensioned::lookupOrAddToDict + ( + "Cw2", + coeffDict_, + 0.92 + ) + ), + Cw3_ + ( + dimensioned::lookupOrAddToDict + ( + "Cw3", + coeffDict_, + 0.3 + ) + ), + CwR_ + ( + dimensioned::lookupOrAddToDict + ( + "CwR", + coeffDict_, + 1.5 + ) + ), + Clambda_ + ( + dimensioned::lookupOrAddToDict + ( + "Clambda", + coeffDict_, + 2.495 + ) + ), + CmuStd_ + ( + dimensioned::lookupOrAddToDict + ( + "CmuStd", + coeffDict_, + 0.09 + ) + ), + Prtheta_ + ( + dimensioned::lookupOrAddToDict + ( + "Prtheta", + coeffDict_, + 0.85 + ) + ), + Sigmak_ + ( + dimensioned::lookupOrAddToDict + ( + "Sigmak", + coeffDict_, + 1 + ) + ), + Sigmaw_ + ( + dimensioned::lookupOrAddToDict + ( + "Sigmaw", + coeffDict_, + 1.17 + ) + ), + kt_ + ( + IOobject + ( + "kt", + runTime_.timeName(), + mesh_, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + autoCreateK("kt", mesh_) + ), + omega_ + ( + IOobject + ( + "omega", + runTime_.timeName(), + mesh_, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + autoCreateOmega("omega", mesh_) + ), + kl_ + ( + IOobject + ( + "kl", + runTime_.timeName(), + mesh_, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + autoCreateK("kl", mesh_) + ), + nut_ + ( + IOobject + ( + "nut", + runTime_.timeName(), + mesh_, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + autoCreateNut("nut", mesh_) + ), + y_(mesh_) +{ + bound(kt_, kMin_); + bound(kl_, kMin_); + bound(omega_, omegaMin_); + + nut_ = kt_/(omega_ + omegaMin_); + nut_.correctBoundaryConditions(); + + printCoeffs(); +} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +tmp kkLOmega::R() const +{ + return tmp + ( + new volSymmTensorField + ( + IOobject + ( + "R", + runTime_.timeName(), + mesh_, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + ((2.0/3.0)*I)*(kt_) - nut_*twoSymm(fvc::grad(U_)), + kt_.boundaryField().types() + ) + ); +} + + +tmp kkLOmega::devReff() const +{ + return tmp + ( + new volSymmTensorField + ( + IOobject + ( + "devRhoReff", + runTime_.timeName(), + mesh_, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + -nuEff()*dev(twoSymm(fvc::grad(U_))) + ) + ); +} + + +tmp kkLOmega::divDevReff(volVectorField& U) const +{ + return + ( + - fvm::laplacian(nuEff(), U) + - fvc::div(nuEff()*dev(T(fvc::grad(U)))) + ); +} + + +bool kkLOmega::read() +{ + if (RASModel::read()) + { + A0_.readIfPresent(coeffDict()); + As_.readIfPresent(coeffDict()); + Av_.readIfPresent(coeffDict()); + Abp_.readIfPresent(coeffDict()); + Anat_.readIfPresent(coeffDict()); + Abp_.readIfPresent(coeffDict()); + Ats_.readIfPresent(coeffDict()); + CbpCrit_.readIfPresent(coeffDict()); + Cnc_.readIfPresent(coeffDict()); + CnatCrit_.readIfPresent(coeffDict()); + Cint_.readIfPresent(coeffDict()); + CtsCrit_.readIfPresent(coeffDict()); + CrNat_.readIfPresent(coeffDict()); + C11_.readIfPresent(coeffDict()); + C12_.readIfPresent(coeffDict()); + CR_.readIfPresent(coeffDict()); + CalphaTheta_.readIfPresent(coeffDict()); + Css_.readIfPresent(coeffDict()); + CtauL_.readIfPresent(coeffDict()); + Cw1_.readIfPresent(coeffDict()); + Cw2_.readIfPresent(coeffDict()); + Cw3_.readIfPresent(coeffDict()); + CwR_.readIfPresent(coeffDict()); + Clambda_.readIfPresent(coeffDict()); + CmuStd_.readIfPresent(coeffDict()); + Prtheta_.readIfPresent(coeffDict()); + Sigmak_.readIfPresent(coeffDict()); + Sigmaw_.readIfPresent(coeffDict()); + + return true; + } + else + { + return false; + } +} + + +void kkLOmega::correct() +{ + RASModel::correct(); + + if (!turbulence_) + { + return; + } + + if (mesh_.changing()) + { + y_.correct(); + y_.boundaryField() = max(y_.boundaryField(), VSMALL); + } + + + const volScalarField kT(kt_ + kl_); + + const volScalarField lambdaT(sqrt(kT)/(omega_ + omegaMin_)); + + const volScalarField lambdaEff(min(Clambda_*y_, lambdaT)); + + const volScalarField fw + ( + lambdaEff/(lambdaT + dimensionedScalar("SMALL", dimLength, ROOTVSMALL)) + ); + + const volTensorField gradU(fvc::grad(U_)); + + const volScalarField omega(sqrt(2.0)*mag(skew(gradU))); + + const volScalarField S2(2.0*magSqr(symm(gradU))); + + const volScalarField ktS(fSS(omega)*fw*kt_); + + const volScalarField nuts + ( + fv(sqr(fw)*kt_/nu()/(omega_ + omegaMin_)) + *fINT() + *Cmu(sqrt(S2))*sqrt(ktS)*lambdaEff + ); + const volScalarField Pkt(nuts*S2); + + const volScalarField ktL(kt_ - ktS); + const volScalarField ReOmega(sqr(y_)*omega/nu()); + const volScalarField nutl + ( + min + ( + C11_*fTaul(lambdaEff, ktL)*omega*sqr(lambdaEff) + * sqrt(ktL)*lambdaEff/nu() + + C12_*BetaTS(ReOmega)*ReOmega*sqr(y_)*omega + , + 0.5*(kl_ + ktL)/sqrt(S2) + ) + ); + + const volScalarField Pkl(nutl*S2); + + const volScalarField alphaTEff + ( + alphaT(lambdaEff, fv(sqr(fw)*kt_/nu()/(omega_ + omegaMin_)), ktS) + ); + + // By pass s0urce term divided by kl_ + + const dimensionedScalar fwMin("SMALL", dimless, ROOTVSMALL); + + const volScalarField Rbp + ( + CR_*(1.0 - exp(-gammaBP(omega)()/Abp_))*omega_ + / (fw + fwMin) + ); + + const volScalarField fNatCrit(1.0 - exp(-Cnc_*sqrt(kl_)*y_/nu())); + // Natural source term divided by kl_ + const volScalarField Rnat + ( + CrNat_*(1.0 - exp(-gammaNAT(ReOmega, fNatCrit)/Anat_))*omega + ); + + const volScalarField Dt(nu()*magSqr(fvc::grad(sqrt(kt_)))); + + // Turbulent kinetic energy equation + tmp ktEqn + ( + fvm::ddt(kt_) + + fvm::div(phi_, kt_) + - fvm::Sp(fvc::div(phi_), kt_) + - fvm::laplacian(DkEff(alphaTEff), kt_, "laplacian(alphaTEff,kt)") + == + Pkt + + (Rbp + Rnat)*kl_ + - Dt + - fvm::Sp(omega_, kt_) + ); + + ktEqn().relax(); + ktEqn().boundaryManipulate(kt_.boundaryField()); + + solve(ktEqn); + bound(kt_, kMin_); + + + const volScalarField Dl(nu()*magSqr(fvc::grad(sqrt(kl_)))); + + // Laminar kinetic energy equation + tmp klEqn + ( + fvm::ddt(kl_) + + fvm::div(phi_, kl_) + - fvm::Sp(fvc::div(phi_), kl_) + - fvm::laplacian(nu(), kl_, "laplacian(nu,kl)") + == + Pkl + - fvm::Sp(Rbp, kl_) + - fvm::Sp(Rnat, kl_) + - Dl + ); + + klEqn().relax(); + klEqn().boundaryManipulate(kl_.boundaryField()); + + solve(klEqn); + bound(kl_, kMin_); + + + omega_.boundaryField().updateCoeffs(); + // Turbulence specific dissipation rate equation + tmp omegaEqn + ( + fvm::ddt(omega_) + + fvm::div(phi_, omega_) + - fvm::Sp(fvc::div(phi_), omega_) + - fvm::laplacian + ( + DomegaEff(alphaTEff), + omega_, + "laplacian(alphaTEff,omega)" + ) + == + Cw1_*Pkt*omega_/(kt_ + kMin_) + + fvm::SuSp + ( + (CwR_/(fw + fwMin) - 1.0)*kl_*(Rbp + Rnat)/(kt_ + kMin_) + , omega_ + ) + - fvm::Sp(Cw2_*omega_, omega_) + + Cw3_*fOmega(lambdaEff, lambdaT)*alphaTEff*sqr(fw)*sqrt(kt_)/pow3(y_) + ); + + + omegaEqn().relax(); + omegaEqn().boundaryManipulate(omega_.boundaryField()); + + solve(omegaEqn); + bound(omega_, omegaMin_); + + // Re-calculate viscosity + nut_ = nuts + nutl; + nut_.correctBoundaryConditions(); +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace RASModels +} // End namespace incompressible +} // End namespace Foam + +// ************************************************************************* // diff --git a/src/turbulenceModels/incompressible/RAS/kkLOmega/kkLOmega.H b/src/turbulenceModels/incompressible/RAS/kkLOmega/kkLOmega.H new file mode 100644 index 0000000000..b7d451b315 --- /dev/null +++ b/src/turbulenceModels/incompressible/RAS/kkLOmega/kkLOmega.H @@ -0,0 +1,296 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +Class + Foam::incompressible::RASModels::kkLOmega + +Description + Low Reynolds-number k-kl-omega turbulence model for + incompressible flows. + + Turbulence model described in: + \verbatim + D. Keith Walters, Davor Cokljat + "A Three-Equation Eddy-Viscosity Model for Reynold-Averaged + Navier-Stokes Simulations of Transitional Flow" + \endverbatim + + The default model coefficients correspond to the following: + \verbatim + kkLOmegaCoeffs + { + A0 4.04 + As 2.12 + Av 6.75 + Abp 0.6 + Anat 200 + Ats 200 + CbpCrit 1.2 + Cnc 0.1 + CnatCrit 1250 + Cint 0.75 + CtsCrit 1000 + CrNat 0.02 + C11 3.4e-6 + C12 1.0e-10 + CR 0.12 + CalphaTheta 0.035 + Css 1.5 + CtauL 4360 + Cw1 0.44 + Cw2 0.92 + Cw3 0.3 + CwR 1.5 + Clambda 2.495 + CmuStd 0.09 + Prtheta 0.85 + Sigmak 1 + Sigmaw 1.17 + } + \endverbatim + +SourceFiles + kkLOmega.C + +\*---------------------------------------------------------------------------*/ + +#ifndef kkLOmega_H +#define kkLOmega_H + +#include "RASModel.H" +#include "wallDist.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace incompressible +{ +namespace RASModels +{ + +/*---------------------------------------------------------------------------*\ + Class kkLOmega Declaration +\*---------------------------------------------------------------------------*/ + +class kkLOmega +: + public RASModel +{ + // Private memmber functions + + tmp fv(const volScalarField& Ret) const; + + tmp fINT() const; + + tmp fSS(const volScalarField& omega) const; + + tmp Cmu(const volScalarField& S) const; + + tmp BetaTS(const volScalarField& Rew) const; + + tmp fTaul + ( + const volScalarField& lambdaEff, + const volScalarField& ktL + ) const; + + tmp alphaT + ( + const volScalarField& lambdaEff, + const volScalarField& fv, + const volScalarField& ktS + ) const; + + tmp fOmega + ( + const volScalarField& lambdaEff, + const volScalarField& lambdaT + ) const; + + tmp gammaBP(const volScalarField& omega) const; + + tmp gammaNAT + ( + const volScalarField& ReOmega, + const volScalarField& fNatCrit + ) const; + +protected: + + // Protected data + + // Model coefficients + + dimensionedScalar A0_; + dimensionedScalar As_; + dimensionedScalar Av_; + dimensionedScalar Abp_; + dimensionedScalar Anat_; + dimensionedScalar Ats_; + dimensionedScalar CbpCrit_; + dimensionedScalar Cnc_; + dimensionedScalar CnatCrit_; + dimensionedScalar Cint_; + dimensionedScalar CtsCrit_; + dimensionedScalar CrNat_; + dimensionedScalar C11_; + dimensionedScalar C12_; + dimensionedScalar CR_; + dimensionedScalar CalphaTheta_; + dimensionedScalar Css_; + dimensionedScalar CtauL_; + dimensionedScalar Cw1_; + dimensionedScalar Cw2_; + dimensionedScalar Cw3_; + dimensionedScalar CwR_; + dimensionedScalar Clambda_; + dimensionedScalar CmuStd_; + dimensionedScalar Prtheta_; + dimensionedScalar Sigmak_; + dimensionedScalar Sigmaw_; + + + // Fields + + volScalarField kt_; + volScalarField omega_; + volScalarField kl_; + volScalarField nut_; + wallDist y_; + + +public: + + //- Runtime type information + TypeName("kkLOmega"); + + // Constructors + + //- Construct from components + kkLOmega + ( + const volVectorField& U, + const surfaceScalarField& phi, + transportModel& transport, + const word& turbulenceModelName = turbulenceModel::typeName, + const word& modelName = typeName + ); + + + //- Destructor + virtual ~kkLOmega() + {} + + + // Member Functions + + //- Return the turbulence viscosity + virtual tmp nut() const + { + return nut_; + } + + //- Return the effective diffusivity for k + tmp DkEff(const volScalarField& alphaT) const + { + return tmp + ( + new volScalarField("DkEff", alphaT/Sigmak_ + nu()) + ); + } + + //- Return the effective diffusivity for omega + tmp DomegaEff(const volScalarField& alphaT) const + { + return tmp + ( + new volScalarField("DomegaEff", alphaT/Sigmaw_ + nu()) + ); + } + + //- Return the laminar kinetic energy + virtual tmp kl() const + { + return kl_; + } + + //- Return the turbulence kinetic energy + virtual tmp k() const + { + return kt_; + } + + //- Return the turbulence specific dissipation rate + virtual tmp omega() const + { + return omega_; + } + + //- Return the turbulence kinetic energy dissipation rate + virtual tmp epsilon() const + { + return tmp + ( + new volScalarField + ( + IOobject + ( + "epsilon", + mesh_.time().timeName(), + mesh_ + ), + kt_*omega_, + omega_.boundaryField().types() + ) + ); + } + + //- Return the Reynolds stress tensor + virtual tmp R() const; + + //- Return the effective stress tensor including the laminar stress + virtual tmp devReff() const; + + //- Return the source term for the momentum equation + virtual tmp divDevReff(volVectorField& U) const; + + //- Solve the turbulence equations and correct the turbulence viscosity + virtual void correct(); + + //- Read RASProperties dictionary + virtual bool read(); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace RASModels +} // End namespace incompressible +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // From 242684dc3d396df725b2480193f71e04ae500129 Mon Sep 17 00:00:00 2001 From: Henry Date: Thu, 27 Oct 2011 17:19:26 +0100 Subject: [PATCH 03/19] Thermo: limit the temperature during energy->temperature conversion --- .../specie/thermo/eConst/eConstThermo.H | 3 ++ .../specie/thermo/eConst/eConstThermoI.H | 10 +++++ .../specie/thermo/hConst/hConstThermo.H | 3 ++ .../specie/thermo/hConst/hConstThermoI.H | 10 +++++ .../thermo/hPolynomial/hPolynomialThermo.H | 3 ++ .../thermo/hPolynomial/hPolynomialThermoI.H | 10 +++++ .../specie/thermo/janaf/janafThermo.H | 6 +-- .../specie/thermo/janaf/janafThermoI.H | 40 ++++++++++--------- .../specie/thermo/specieThermo/specieThermo.H | 5 ++- .../thermo/specieThermo/specieThermoI.H | 33 ++++++++++++--- 10 files changed, 95 insertions(+), 28 deletions(-) diff --git a/src/thermophysicalModels/specie/thermo/eConst/eConstThermo.H b/src/thermophysicalModels/specie/thermo/eConst/eConstThermo.H index 9154126ada..bd8bdb4573 100644 --- a/src/thermophysicalModels/specie/thermo/eConst/eConstThermo.H +++ b/src/thermophysicalModels/specie/thermo/eConst/eConstThermo.H @@ -135,6 +135,9 @@ public: // Member Functions + //- Limit the temperature to be in the range Tlow_ to Thigh_ + inline scalar limit(const scalar T) const; + // Fundamental properties //- Heat capacity at constant pressure [J/(kmol K)] diff --git a/src/thermophysicalModels/specie/thermo/eConst/eConstThermoI.H b/src/thermophysicalModels/specie/thermo/eConst/eConstThermoI.H index f031aa63f6..31bd392b5d 100644 --- a/src/thermophysicalModels/specie/thermo/eConst/eConstThermoI.H +++ b/src/thermophysicalModels/specie/thermo/eConst/eConstThermoI.H @@ -89,6 +89,16 @@ Foam::eConstThermo::New(const dictionary& dict) // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // +template +inline Foam::scalar Foam::eConstThermo::limit +( + const scalar T +) const +{ + return T; +} + + template inline Foam::scalar Foam::eConstThermo::cp ( diff --git a/src/thermophysicalModels/specie/thermo/hConst/hConstThermo.H b/src/thermophysicalModels/specie/thermo/hConst/hConstThermo.H index 014de5eef2..9a8352a064 100644 --- a/src/thermophysicalModels/specie/thermo/hConst/hConstThermo.H +++ b/src/thermophysicalModels/specie/thermo/hConst/hConstThermo.H @@ -133,6 +133,9 @@ public: // Member Functions + //- Limit the temperature to be in the range Tlow_ to Thigh_ + inline scalar limit(const scalar T) const; + // Fundamental properties //- Heat capacity at constant pressure [J/(kmol K)] diff --git a/src/thermophysicalModels/specie/thermo/hConst/hConstThermoI.H b/src/thermophysicalModels/specie/thermo/hConst/hConstThermoI.H index 76ec41b620..0a4a35358c 100644 --- a/src/thermophysicalModels/specie/thermo/hConst/hConstThermoI.H +++ b/src/thermophysicalModels/specie/thermo/hConst/hConstThermoI.H @@ -89,6 +89,16 @@ Foam::hConstThermo::New(const dictionary& dict) // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // +template +inline Foam::scalar Foam::hConstThermo::limit +( + const scalar T +) const +{ + return T; +} + + template inline Foam::scalar Foam::hConstThermo::cp ( diff --git a/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermo.H b/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermo.H index 7a4dcb3867..31954a4a56 100644 --- a/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermo.H +++ b/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermo.H @@ -151,6 +151,9 @@ public: // Member Functions + //- Limit the temperature to be in the range Tlow_ to Thigh_ + inline scalar limit(const scalar T) const; + // Fundamental properties //- Heat capacity at constant pressure [J/(kmol K)] diff --git a/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermoI.H b/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermoI.H index 620019e6d1..a164392c61 100644 --- a/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermoI.H +++ b/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermoI.H @@ -82,6 +82,16 @@ inline Foam::hPolynomialThermo::hPolynomialThermo // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // +template +inline Foam::scalar Foam::hPolynomialThermo::limit +( + const scalar T +) const +{ + return T; +} + + template inline Foam::scalar Foam::hPolynomialThermo::cp ( diff --git a/src/thermophysicalModels/specie/thermo/janaf/janafThermo.H b/src/thermophysicalModels/specie/thermo/janaf/janafThermo.H index 9c414afaac..b552fad094 100644 --- a/src/thermophysicalModels/specie/thermo/janaf/janafThermo.H +++ b/src/thermophysicalModels/specie/thermo/janaf/janafThermo.H @@ -119,9 +119,6 @@ private: //- Check that input data is valid void checkInputData() const; - //- Check given temperature is within the range of the fitted coeffs - inline void checkT(const scalar T) const; - //- Return the coefficients corresponding to the given temperature inline const coeffArray& coeffs(const scalar T) const; @@ -153,6 +150,9 @@ public: // Member Functions + //- Limit the temperature to be in the range Tlow_ to Thigh_ + inline scalar limit(const scalar T) const; + // Fundamental properties //- Heat capacity at constant pressure [J/(kmol K)] diff --git a/src/thermophysicalModels/specie/thermo/janaf/janafThermoI.H b/src/thermophysicalModels/specie/thermo/janaf/janafThermoI.H index 05f0c0c0cc..5740816c23 100644 --- a/src/thermophysicalModels/specie/thermo/janaf/janafThermoI.H +++ b/src/thermophysicalModels/specie/thermo/janaf/janafThermoI.H @@ -52,22 +52,6 @@ inline Foam::janafThermo::janafThermo } -template -inline void Foam::janafThermo::checkT(const scalar T) const -{ - if (T < Tlow_ || T > Thigh_) - { - FatalErrorIn - ( - "janafThermo::checkT(const scalar T) const" - ) << "attempt to use janafThermo" - " out of temperature range " - << Tlow_ << " -> " << Thigh_ << "; T = " << T - << abort(FatalError); - } -} - - template inline const typename Foam::janafThermo::coeffArray& Foam::janafThermo::coeffs @@ -75,8 +59,6 @@ Foam::janafThermo::coeffs const scalar T ) const { - checkT(T); - if (T < Tcommon_) { return lowCpCoeffs_; @@ -112,6 +94,28 @@ inline Foam::janafThermo::janafThermo // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // +template +inline Foam::scalar Foam::janafThermo::limit(const scalar T) const +{ + if (T < Tlow_ || T > Thigh_) + { + WarningIn + ( + "janafThermo::limit(const scalar T) const" + ) << "attempt to use janafThermo" + " out of temperature range " + << Tlow_ << " -> " << Thigh_ << "; T = " << T + << endl; + + return min(max(T, Tlow_), Thigh_); + } + else + { + return T; + } +} + + template inline Foam::scalar Foam::janafThermo::cp ( diff --git a/src/thermophysicalModels/specie/thermo/specieThermo/specieThermo.H b/src/thermophysicalModels/specie/thermo/specieThermo/specieThermo.H index c6f0f7a316..93268b9e6b 100644 --- a/src/thermophysicalModels/specie/thermo/specieThermo/specieThermo.H +++ b/src/thermophysicalModels/specie/thermo/specieThermo/specieThermo.H @@ -103,14 +103,15 @@ class specieThermo // Private Member Functions - //- return the temperature corresponding to the value of the + //- Return the temperature corresponding to the value of the // thermodynamic property f, given the function f = F(T) and dF(T)/dT inline scalar T ( scalar f, scalar T0, scalar (specieThermo::*F)(const scalar) const, - scalar (specieThermo::*dFdT)(const scalar) const + scalar (specieThermo::*dFdT)(const scalar) const, + scalar (specieThermo::*limit)(const scalar) const ) const; diff --git a/src/thermophysicalModels/specie/thermo/specieThermo/specieThermoI.H b/src/thermophysicalModels/specie/thermo/specieThermo/specieThermoI.H index f8632376c2..b2d2edd5c1 100644 --- a/src/thermophysicalModels/specie/thermo/specieThermo/specieThermoI.H +++ b/src/thermophysicalModels/specie/thermo/specieThermo/specieThermoI.H @@ -43,7 +43,8 @@ inline Foam::scalar Foam::specieThermo::T scalar f, scalar T0, scalar (specieThermo::*F)(const scalar) const, - scalar (specieThermo::*dFdT)(const scalar) const + scalar (specieThermo::*dFdT)(const scalar) const, + scalar (specieThermo::*limit)(const scalar) const ) const { scalar Test = T0; @@ -54,7 +55,8 @@ inline Foam::scalar Foam::specieThermo::T do { Test = Tnew; - Tnew = Test - ((this->*F)(Test) - f)/(this->*dFdT)(Test); + Tnew = + (this->*limit)(Test - ((this->*F)(Test) - f)/(this->*dFdT)(Test)); if (iter++ > maxIter_) { @@ -276,7 +278,14 @@ inline Foam::scalar Foam::specieThermo::TH const scalar T0 ) const { - return T(h, T0, &specieThermo::H, &specieThermo::Cp); + return T + ( + h, + T0, + &specieThermo::H, + &specieThermo::Cp, + &specieThermo::limit + ); } @@ -287,7 +296,14 @@ inline Foam::scalar Foam::specieThermo::THs const scalar T0 ) const { - return T(hs, T0, &specieThermo::Hs, &specieThermo::Cp); + return T + ( + hs, + T0, + &specieThermo::Hs, + &specieThermo::Cp, + &specieThermo::limit + ); } @@ -298,7 +314,14 @@ inline Foam::scalar Foam::specieThermo::TE const scalar T0 ) const { - return T(e, T0, &specieThermo::E, &specieThermo::Cv); + return T + ( + e, + T0, + &specieThermo::E, + &specieThermo::Cv, + &specieThermo::limit + ); } From a0459cb2532e3c37af512ce771cdd152327171d1 Mon Sep 17 00:00:00 2001 From: Henry Date: Thu, 27 Oct 2011 17:19:57 +0100 Subject: [PATCH 04/19] Minor reformatting --- tutorials/incompressible/MRFSimpleFoam/mixerVessel2D/Allrun | 1 - 1 file changed, 1 deletion(-) diff --git a/tutorials/incompressible/MRFSimpleFoam/mixerVessel2D/Allrun b/tutorials/incompressible/MRFSimpleFoam/mixerVessel2D/Allrun index 9210b03c39..caad132d7b 100755 --- a/tutorials/incompressible/MRFSimpleFoam/mixerVessel2D/Allrun +++ b/tutorials/incompressible/MRFSimpleFoam/mixerVessel2D/Allrun @@ -1,7 +1,6 @@ #!/bin/sh cd ${0%/*} || exit 1 # run from this directory - # Source tutorial run functions . $WM_PROJECT_DIR/bin/tools/RunFunctions From 9a7bdf6ba97fd1d4d38bd8e64c38f7590961c40f Mon Sep 17 00:00:00 2001 From: Henry Date: Thu, 27 Oct 2011 17:22:02 +0100 Subject: [PATCH 05/19] Corrected line length --- src/thermophysicalModels/specie/thermo/janaf/janafThermoI.H | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/src/thermophysicalModels/specie/thermo/janaf/janafThermoI.H b/src/thermophysicalModels/specie/thermo/janaf/janafThermoI.H index 5740816c23..de54035ea4 100644 --- a/src/thermophysicalModels/specie/thermo/janaf/janafThermoI.H +++ b/src/thermophysicalModels/specie/thermo/janaf/janafThermoI.H @@ -95,7 +95,10 @@ inline Foam::janafThermo::janafThermo // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // template -inline Foam::scalar Foam::janafThermo::limit(const scalar T) const +inline Foam::scalar Foam::janafThermo::limit +( + const scalar T +) const { if (T < Tlow_ || T > Thigh_) { From a0352b0ae0c585447e70edffaac123e78a33c2d9 Mon Sep 17 00:00:00 2001 From: mattijs Date: Fri, 28 Oct 2011 11:31:23 +0100 Subject: [PATCH 06/19] ENH: polyTopoChange: efficient sorting --- .../polyTopoChange/polyTopoChange.C | 26 ++++++++++++++----- 1 file changed, 20 insertions(+), 6 deletions(-) diff --git a/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.C b/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.C index 622748be1e..63f533e090 100644 --- a/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.C +++ b/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.C @@ -680,13 +680,17 @@ void Foam::polyTopoChange::getFaceOrder // First unassigned face label newFaceI = 0; + labelList nbr; + labelList order; + forAll(cellMap_, cellI) { label startOfCell = cellFaceOffsets[cellI]; label nFaces = cellFaceOffsets[cellI+1] - startOfCell; // Neighbouring cells - SortableList