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

RAS/kOmegaSST: Added F3 coefficient for rough-walls from

        Hellsten, A.
        "Some Improvements in Menter’s k-omega-SST turbulence model"
        29th AIAA Fluid Dynamics Conference,
        AIAA-98-2554,
        June 1998.
parent a1af190a
No related merge requests found
......@@ -85,6 +85,18 @@ tmp<volScalarField> kOmegaSST::F2() const
}
tmp<volScalarField> kOmegaSST::F3() const
{
tmp<volScalarField> arg3 = min
(
150*nu()/(omega_*sqr(y_)),
scalar(10)
);
return 1 - tanh(pow4(arg3));
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
kOmegaSST::kOmegaSST
......@@ -188,6 +200,15 @@ kOmegaSST::kOmegaSST
0.31
)
),
b1_
(
dimensioned<scalar>::lookupOrAddToDict
(
"b1",
coeffDict_,
1.0
)
),
c1_
(
dimensioned<scalar>::lookupOrAddToDict
......@@ -246,7 +267,7 @@ kOmegaSST::kOmegaSST
/ max
(
a1_*omega_,
F2()*sqrt(2.0)*mag(symm(fvc::grad(U_)))
b1_*F2()*F3()*sqrt(2.0)*mag(symm(fvc::grad(U_)))
)
);
nut_.correctBoundaryConditions();
......@@ -416,7 +437,7 @@ void kOmegaSST::correct()
// Re-calculate viscosity
nut_ = a1_*k_/max(a1_*omega_, F2()*sqrt(S2));
nut_ = a1_*k_/max(a1_*omega_, b1_*F2()*F3()*sqrt(S2));
nut_.correctBoundaryConditions();
}
......
......@@ -24,19 +24,25 @@ License
Class
Foam::incompressible::RASModels::kOmegaSST
Group
grpIcoRASTurbulence
Description
Implementation of the k-omega-SST turbulence model for incompressible
flows.
Turbulence model described in:
\verbatim
Menter, F., Esch, T.
"Elements of Industrial Heat Transfer Prediction"
Menter, F., Esch, T.,
"Elements of Industrial Heat Transfer Prediction",
16th Brazilian Congress of Mechanical Engineering (COBEM),
Nov. 2001
Nov. 2001.
\endverbatim
with the addition of the F3 term for rough walls from
\verbatim
Hellsten, A.
"Some Improvements in Menter’s k-omega-SST turbulence model"
29th AIAA Fluid Dynamics Conference,
AIAA-98-2554,
June 1998.
\endverbatim
Note that this implementation is written in terms of alpha diffusion
......@@ -69,6 +75,7 @@ Description
gamma1 0.5532;
gamma2 0.4403;
a1 0.31;
b1 1.0;
c1 10.0;
}
\endverbatim
......@@ -122,6 +129,7 @@ protected:
dimensionedScalar betaStar_;
dimensionedScalar a1_;
dimensionedScalar b1_;
dimensionedScalar c1_;
//- Wall distance field
......@@ -139,6 +147,7 @@ protected:
tmp<volScalarField> F1(const volScalarField& CDkOmega) const;
tmp<volScalarField> F2() const;
tmp<volScalarField> F3() const;
tmp<volScalarField> blend
(
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment