Skip to content
Snippets Groups Projects
Commit a9b8f174 authored by Henry's avatar Henry
Browse files

non-linear turbulence models: corrected generation term and tested

parent d05b8e6d
Branches
Tags
No related merge requests found
...@@ -306,15 +306,13 @@ void LienCubicKE::correct() ...@@ -306,15 +306,13 @@ void LienCubicKE::correct()
tmp<volTensorField> tgradU = fvc::grad(U_); tmp<volTensorField> tgradU = fvc::grad(U_);
const volTensorField& gradU = tgradU(); const volTensorField& gradU = tgradU();
// generation term
tmp<volScalarField> S2 = symm(gradU) && gradU;
volScalarField G volScalarField G
( (
GName(), GName(),
Cmu_*sqr(k_)/epsilon_*S2 - (nonlinearStress_ && gradU) (nut_*twoSymm(gradU) - nonlinearStress_) && gradU
); );
// Update epsilon and G at the wall // Update epsilon and G at the wall
epsilon_.boundaryField().updateCoeffs(); epsilon_.boundaryField().updateCoeffs();
...@@ -330,15 +328,12 @@ void LienCubicKE::correct() ...@@ -330,15 +328,12 @@ void LienCubicKE::correct()
); );
epsEqn().relax(); epsEqn().relax();
epsEqn().boundaryManipulate(epsilon_.boundaryField()); epsEqn().boundaryManipulate(epsilon_.boundaryField());
solve(epsEqn); solve(epsEqn);
bound(epsilon_, epsilonMin_); bound(epsilon_, epsilonMin_);
// Turbulent kinetic energy equation // Turbulent kinetic energy equation
tmp<fvScalarMatrix> kEqn tmp<fvScalarMatrix> kEqn
( (
fvm::ddt(k_) fvm::ddt(k_)
......
...@@ -371,9 +371,6 @@ void LienCubicKELowRe::correct() ...@@ -371,9 +371,6 @@ void LienCubicKELowRe::correct()
tmp<volTensorField> tgradU = fvc::grad(U_); tmp<volTensorField> tgradU = fvc::grad(U_);
const volTensorField& gradU = tgradU(); const volTensorField& gradU = tgradU();
// generation term
tmp<volScalarField> S2 = symm(gradU) && gradU;
yStar_ = sqrt(k_)*y_/nu() + SMALL; yStar_ = sqrt(k_)*y_/nu() + SMALL;
tmp<volScalarField> Rt = sqr(k_)/(nu()*epsilon_); tmp<volScalarField> Rt = sqr(k_)/(nu()*epsilon_);
...@@ -385,7 +382,7 @@ void LienCubicKELowRe::correct() ...@@ -385,7 +382,7 @@ void LienCubicKELowRe::correct()
volScalarField G volScalarField G
( (
GName(), GName(),
Cmu_*fMu()*sqr(k_)/epsilon_*S2 - (nonlinearStress_ && gradU) (nut_*twoSymm(gradU) - nonlinearStress_) && gradU
); );
// Dissipation equation // Dissipation equation
...@@ -404,16 +401,11 @@ void LienCubicKELowRe::correct() ...@@ -404,16 +401,11 @@ void LienCubicKELowRe::correct()
); );
epsEqn().relax(); epsEqn().relax();
#include "LienCubicKELowReSetWallDissipation.H"
#include "wallDissipationI.H"
solve(epsEqn); solve(epsEqn);
bound(epsilon_, epsilonMin_); bound(epsilon_, epsilonMin_);
// Turbulent kinetic energy equation // Turbulent kinetic energy equation
tmp<fvScalarMatrix> kEqn tmp<fvScalarMatrix> kEqn
( (
fvm::ddt(k_) fvm::ddt(k_)
......
...@@ -243,21 +243,15 @@ void LienLeschzinerLowRe::correct() ...@@ -243,21 +243,15 @@ void LienLeschzinerLowRe::correct()
return; return;
} }
scalar Cmu75 = pow(Cmu_.value(), 0.75); tmp<volTensorField> tgradU = fvc::grad(U_);
volScalarField G(GName(), nut_*(tgradU() && twoSymm(tgradU())));
const volTensorField gradU(fvc::grad(U_)); tgradU.clear();
// generation term
tmp<volScalarField> S2 = symm(gradU) && gradU;
scalar Cmu75 = pow(Cmu_.value(), 0.75);
yStar_ = sqrt(k_)*y_/nu() + SMALL; yStar_ = sqrt(k_)*y_/nu() + SMALL;
tmp<volScalarField> Rt = sqr(k_)/(nu()*epsilon_); tmp<volScalarField> Rt = sqr(k_)/(nu()*epsilon_);
const volScalarField f2(scalar(1) - 0.3*exp(-sqr(Rt))); const volScalarField f2(scalar(1) - 0.3*exp(-sqr(Rt)));
volScalarField G(GName(), Cmu_*fMu()*sqr(k_)/epsilon_*S2);
// Dissipation equation // Dissipation equation
tmp<fvScalarMatrix> epsEqn tmp<fvScalarMatrix> epsEqn
( (
...@@ -266,24 +260,21 @@ void LienLeschzinerLowRe::correct() ...@@ -266,24 +260,21 @@ void LienLeschzinerLowRe::correct()
- fvm::laplacian(DepsilonEff(), epsilon_) - fvm::laplacian(DepsilonEff(), epsilon_)
== ==
C1_*G*epsilon_/k_ C1_*G*epsilon_/k_
// E-term // E-term
+ C2_*f2*Cmu75*sqrt(k_) + C2_*f2*Cmu75*sqrt(k_)
/(kappa_*y_*(scalar(1) - exp(-Aepsilon_*yStar_))) /(kappa_*y_*(scalar(1) - exp(-Aepsilon_*yStar_)))
*exp(-Amu_*sqr(yStar_))*epsilon_ *exp(-Amu_*sqr(yStar_))*epsilon_
- fvm::Sp(C2_*f2*epsilon_/k_, epsilon_) - fvm::Sp(C2_*f2*epsilon_/k_, epsilon_)
); );
epsEqn().relax(); epsEqn().relax();
#include "LienLeschzinerLowReSetWallDissipation.H"
#include "wallDissipationI.H"
solve(epsEqn); solve(epsEqn);
bound(epsilon_, epsilonMin_); bound(epsilon_, epsilonMin_);
// Turbulent kinetic energy equation // Turbulent kinetic energy equation
tmp<fvScalarMatrix> kEqn tmp<fvScalarMatrix> kEqn
( (
fvm::ddt(k_) fvm::ddt(k_)
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment