tutorials/incompressible/pimpleFoam/elipsekkLOmega: improved variable naming and updated test-case

This commit is contained in:
Henry 2015-01-27 08:43:20 +00:00
parent f4b870a3d3
commit 8ecf9e447a
8 changed files with 104 additions and 111 deletions

View file

@ -26,8 +26,6 @@ License
#include "kkLOmega.H"
#include "addToRunTimeSelectionTable.H"
#include "backwardsCompatibilityWallFunctions.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
@ -63,9 +61,9 @@ tmp<volScalarField> kkLOmega::fINT() const
}
tmp<volScalarField> kkLOmega::fSS(const volScalarField& omega) const
tmp<volScalarField> kkLOmega::fSS(const volScalarField& Omega) const
{
return(exp(-sqr(Css_*nu()*omega/(kt_ + kMin_))));
return(exp(-sqr(Css_*nu()*Omega/(kt_ + kMin_))));
}
@ -85,7 +83,7 @@ tmp<volScalarField> kkLOmega::fTaul
(
const volScalarField& lambdaEff,
const volScalarField& ktL,
const volScalarField& omega
const volScalarField& Omega
) const
{
return
@ -98,7 +96,7 @@ tmp<volScalarField> kkLOmega::fTaul
(
sqr
(
lambdaEff*omega
lambdaEff*Omega
+ dimensionedScalar
(
"ROOTVSMALL",
@ -153,7 +151,7 @@ tmp<volScalarField> kkLOmega::fOmega
}
tmp<volScalarField> kkLOmega::phiBP(const volScalarField& omega) const
tmp<volScalarField> kkLOmega::phiBP(const volScalarField& Omega) const
{
return
(
@ -163,11 +161,11 @@ tmp<volScalarField> kkLOmega::phiBP(const volScalarField& omega) const
(
kt_/nu()
/ (
omega
Omega
+ dimensionedScalar
(
"ROTVSMALL",
omega.dimensions(),
Omega.dimensions(),
ROOTVSMALL
)
)
@ -201,6 +199,12 @@ tmp<volScalarField> kkLOmega::phiNAT
}
tmp<volScalarField> kkLOmega::D(const volScalarField& k) const
{
return nu()*magSqr(fvc::grad(sqrt(k)));
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
kkLOmega::kkLOmega
@ -464,22 +468,10 @@ kkLOmega::kkLOmega
"kt",
runTime_.timeName(),
mesh_,
IOobject::NO_READ,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
autoCreateK("kt", mesh_)
),
omega_
(
IOobject
(
"omega",
runTime_.timeName(),
mesh_,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
autoCreateOmega("omega", mesh_)
mesh_
),
kl_
(
@ -488,10 +480,32 @@ kkLOmega::kkLOmega
"kl",
runTime_.timeName(),
mesh_,
IOobject::NO_READ,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
autoCreateK("kl", mesh_)
mesh_
),
omega_
(
IOobject
(
"omega",
runTime_.timeName(),
mesh_,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh_
),
epsilon_
(
IOobject
(
"epsilon",
runTime_.timeName(),
mesh_
),
kt_*omega_ + D(kl_) + D(kt_)
),
nut_
(
@ -500,18 +514,19 @@ kkLOmega::kkLOmega
"nut",
runTime_.timeName(),
mesh_,
IOobject::NO_READ,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
autoCreateNut("nut", mesh_)
mesh_
),
y_(mesh_)
{
bound(kt_, kMin_);
bound(kl_, kMin_);
bound(omega_, omegaMin_);
bound(epsilon_, epsilonMin_);
nut_ = kt_/(omega_ + omegaMin_);
// Evaluating nut_ is complex so start from the field read from file
nut_.correctBoundaryConditions();
printCoeffs();
@ -534,7 +549,7 @@ tmp<volSymmTensorField> kkLOmega::R() const
IOobject::NO_READ,
IOobject::NO_WRITE
),
((2.0/3.0)*I)*(kt_ + kl_) - nut_*twoSymm(fvc::grad(U_)),
((2.0/3.0)*I)*k() - nut_*twoSymm(fvc::grad(U_)),
kt_.boundaryField().types()
)
);
@ -658,36 +673,34 @@ void kkLOmega::correct()
)
);
const volTensorField gradU(fvc::grad(U_));
tmp<volTensorField> tgradU(fvc::grad(U_));
const volTensorField& gradU = tgradU();
const volScalarField omega(sqrt(2.0)*mag(skew(gradU)));
const volScalarField Omega(sqrt(2.0)*mag(skew(gradU)));
const volScalarField S2(2.0*magSqr(symm(gradU)));
const volScalarField S2(2.0*magSqr(dev(symm(gradU))));
const volScalarField ktS(fSS(omega)*fw*kt_);
const volScalarField ktS(fSS(Omega)*fw*kt_);
const volScalarField nuts
(
fv
(
sqr(fw)*kt_/nu()/(omega_ + omegaMin_)
)
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 ReOmega(sqr(y_)*Omega/nu());
const volScalarField nutl
(
min
(
C11_*fTaul(lambdaEff, ktL, omega)*omega*sqr(lambdaEff)
C11_*fTaul(lambdaEff, ktL, Omega)*Omega*sqr(lambdaEff)
*sqrt(ktL)*lambdaEff/nu()
+ C12_*BetaTS(ReOmega)*ReOmega*sqr(y_)*omega
+ C12_*BetaTS(ReOmega)*ReOmega*sqr(y_)*Omega
,
0.5*(kl_ + ktL)/sqrt(S2)
0.5*(kl_ + ktL)/(sqrt(S2) + omegaMin_)
)
);
@ -704,7 +717,7 @@ void kkLOmega::correct()
const volScalarField Rbp
(
CR_*(1.0 - exp(-phiBP(omega)()/Abp_))*omega_
CR_*(1.0 - exp(-phiBP(Omega)()/Abp_))*omega_
/(fw + fwMin)
);
@ -713,7 +726,7 @@ void kkLOmega::correct()
// Natural source term divided by kl_
const volScalarField Rnat
(
CrNat_*(1.0 - exp(-phiNAT(ReOmega, fNatCrit)/Anat_))*omega
CrNat_*(1.0 - exp(-phiNAT(ReOmega, fNatCrit)/Anat_))*Omega
);
@ -727,9 +740,9 @@ void kkLOmega::correct()
- fvm::laplacian(DomegaEff(alphaTEff), omega_)
==
Cw1_*Pkt*omega_/(kt_ + kMin_)
+ fvm::SuSp
- fvm::SuSp
(
(CwR_/(fw + fwMin) - 1.0)*kl_*(Rbp + Rnat)/(kt_ + kMin_)
(1.0 - CwR_/(fw + fwMin))*kl_*(Rbp + Rnat)/(kt_ + kMin_)
, omega_
)
- fvm::Sp(Cw2_*sqr(fw)*omega_, omega_)
@ -745,7 +758,7 @@ void kkLOmega::correct()
bound(omega_, omegaMin_);
const volScalarField Dl(nu()*magSqr(fvc::grad(sqrt(kl_))));
const volScalarField Dl(D(kl_));
// Laminar kinetic energy equation
tmp<fvScalarMatrix> klEqn
@ -755,9 +768,7 @@ void kkLOmega::correct()
- fvm::laplacian(nu(), kl_)
==
Pkl
- fvm::Sp(Rbp, kl_)
- fvm::Sp(Rnat, kl_)
- fvm::Sp(Dl/max(kl_, kMin_), kl_)
- fvm::Sp(Rbp + Rnat + Dl/(kl_ + kMin_), kl_)
);
klEqn().relax();
@ -767,7 +778,7 @@ void kkLOmega::correct()
bound(kl_, kMin_);
const volScalarField Dt(nu()*magSqr(fvc::grad(sqrt(kt_))));
const volScalarField Dt(D(kt_));
// Turbulent kinetic energy equation
tmp<fvScalarMatrix> ktEqn
@ -778,8 +789,7 @@ void kkLOmega::correct()
==
Pkt
+ (Rbp + Rnat)*kl_
- fvm::Sp(Dt/max(kt_, kMin_), kt_)
- fvm::Sp(omega_, kt_)
- fvm::Sp(omega_ + Dt/(kt_+ kMin_), kt_)
);
ktEqn().relax();
@ -789,6 +799,11 @@ void kkLOmega::correct()
bound(kt_, kMin_);
// Update total fluctuation kinetic energy dissipation rate
epsilon_ = kt_*omega_ + Dl + Dt;
bound(epsilon_, epsilonMin_);
// Re-calculate turbulent viscosity
nut_ = nuts + nutl;
nut_.correctBoundaryConditions();

View file

@ -152,6 +152,8 @@ class kkLOmega
const volScalarField& fNatCrit
) const;
tmp<volScalarField> D(const volScalarField& k) const;
protected:
@ -191,8 +193,9 @@ protected:
// Fields
volScalarField kt_;
volScalarField omega_;
volScalarField kl_;
volScalarField omega_;
volScalarField epsilon_;
volScalarField nut_;
//- Wall distance
@ -260,6 +263,12 @@ public:
return kt_;
}
//- Return the turbulence specific dissipation rate
virtual tmp<volScalarField> omega() const
{
return omega_;
}
//- Return the total fluctuation kinetic energy
virtual tmp<volScalarField> k() const
{
@ -279,29 +288,10 @@ public:
);
}
//- Return the turbulence specific dissipation rate
virtual tmp<volScalarField> omega() const
{
return omega_;
}
//- Return the turbulence kinetic energy dissipation rate
//- Return the total fluctuation kinetic energy dissipation rate
virtual tmp<volScalarField> epsilon() const
{
return tmp<volScalarField>
(
new volScalarField
(
IOobject
(
"epsilon",
mesh_.time().timeName(),
mesh_
),
kt_*omega_,
omega_.boundaryField().types()
)
);
return epsilon_;
}
//- Return the Reynolds stress tensor

View file

@ -15,7 +15,7 @@ FoamFile
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
dimensions [ 0 2 -2 0 0 0 0 ];
dimensions [0 2 -2 0 0 0 0];
internalField uniform 0;
@ -23,23 +23,23 @@ boundaryField
{
inlet
{
type fixedValue;
value uniform 0;
type fixedValue;
value uniform 0;
}
outlet
{
type inletOutlet;
inletValue uniform 0;
value uniform 0;
type inletOutlet;
inletValue uniform 0;
value uniform 0;
}
up
{
type symmetry;
type symmetry;
}
hole
{
type fixedValue;
value uniform 0;
type fixedValue;
value uniform 0;
}
frontAndBack
{

View file

@ -15,7 +15,7 @@ FoamFile
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
dimensions [ 0 2 -2 0 0 0 0 ];
dimensions [0 2 -2 0 0 0 0];
internalField uniform 0;

View file

@ -15,7 +15,7 @@ FoamFile
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
dimensions [ 0 0 -1 0 0 0 0 ];
dimensions [0 0 -1 0 0 0 0];
internalField uniform 1e-5;
@ -23,23 +23,22 @@ boundaryField
{
inlet
{
type fixedValue;
value $internalField;
type fixedValue;
value $internalField;
}
outlet
{
type inletOutlet;
inletValue $internalField;
value $internalField;
type inletOutlet;
inletValue $internalField;
value $internalField;
}
up
{
type symmetry;
type symmetry;
}
hole
{
type fixedValue;
value $internalField;
type zeroGradient;
}
frontAndBack
{

View file

@ -23,32 +23,27 @@ ddtSchemes
gradSchemes
{
default Gauss linear;
grad(p) Gauss linear;
grad(U) Gauss linear;
grad(U) cellLimited Gauss linear 1;
}
divSchemes
{
default none;
div(phi,U) Gauss limitedLinearV 1;
div(phi,U) Gauss linearUpwindV grad(U);
div(phi,kl) Gauss limitedLinear 1;
div(phi,kt) Gauss limitedLinear 1;
div(phi,omega) Gauss limitedLinear 1;
div(phi,R) Gauss limitedLinear 1;
div(R) Gauss linear;
div(phi,nuTilda) Gauss limitedLinear 1;
div((nuEff*dev(T(grad(U))))) Gauss linear;
}
laplacianSchemes
{
default Gauss linear corrected;
default Gauss linear corrected;
}
interpolationSchemes
{
default linear;
interpolate(U) linear;
}
snGradSchemes

View file

@ -21,7 +21,7 @@ solvers
{
solver GAMG;
tolerance 1e-06;
relTol 0.1;
relTol 0.01;
smoother DICGaussSeidel;
cacheAgglomeration true;
nCellsInCoarsestLevel 10;
@ -56,22 +56,16 @@ PIMPLE
{
nOuterCorrectors 1;
nCorrectors 2;
nNonOrthogonalCorrectors 0;
nNonOrthogonalCorrectors 1;
pRefCell 0;
pRefValue 0;
}
relaxationFactors
{
fields
{
}
equations
{
"U.*" 1;
"kl.*" 1;
"kt.*" 1;
"omega.*" 1;
".*" 1;
}
}

View file

@ -1,7 +1,7 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 2.3.0 |
| \\ / O peration | Version: 2.3.x |
| \\ / A nd | Web: www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/