37 template<
class BasicTurbulenceModel>
50 template<
class BasicTurbulenceModel>
60 template<
class BasicTurbulenceModel>
68 min(
max(gammaIntEff_, scalar(0.1)), scalar(1))
73 template<
class BasicTurbulenceModel>
96 (1 -
sqr((gammaInt_() - 1.0/ce2_)/(1 - 1.0/ce2_)))
104 template<
class BasicTurbulenceModel>
121 const scalar ReThetat = ReThetat_[celli];
128 + 120.656e-4*ReThetat
129 - 868.230e-6*
sqr(ReThetat)
130 + 696.506e-9*
pow3(ReThetat)
131 - 174.105e-12*
pow4(ReThetat)
133 ReThetat - 593.11 - 0.482*(ReThetat - 1870);
140 template<
class BasicTurbulenceModel>
162 const scalar ReThetat = ReThetat_[celli];
168 - 119.270e-4*ReThetat
169 - 132.567e-6*
sqr(ReThetat);
171 else if (ReThetat < 596)
175 - 123.939e-2*ReThetat
176 + 194.548e-5*
sqr(ReThetat)
177 - 101.695e-8*
pow3(ReThetat);
179 else if (ReThetat < 1200)
181 Flength[celli] = 0.5 - 3
e-4*(ReThetat - 596);
185 Flength[celli] = 0.3188;
188 const scalar Fsublayer =
189 exp(-
sqr(
sqr(y[celli])*omega[celli]/(200*nu[celli])));
191 Flength[celli] = Flength[celli]*(1 - Fsublayer) + 40*Fsublayer;
198 template<
class BasicTurbulenceModel>
225 max(100*
sqrt((2.0/3.0)*k[celli])/Us[celli], scalar(0.027))
240 const scalar lambda0 =
lambda;
244 const scalar Flambda =
251 - 405.689*
pow3(lambda)
255 + 0.275*(1 -
exp(-35*lambda))
259 (1173.51 - 589.428*Tu + 0.2196/
sqr(Tu))
265 const scalar Flambda =
272 -405.689*
pow3(lambda)
276 + 0.275*(1 -
exp(-35*lambda))
280 331.50*
pow((Tu - 0.5658), -0.671)
281 *Flambda*nu[celli]/Us[celli];
284 lambda =
sqr(thetat)/nu[celli]*dUsds[celli];
285 lambda =
max(
min(lambda, 0.1), -0.1);
287 lambdaErr =
mag(lambda - lambda0);
289 maxIter =
max(maxIter, ++iter);
291 }
while (lambdaErr > lambdaErr_);
293 ReThetat0[celli] =
max(thetat*Us[celli]/nu[celli], scalar(20));
296 if (maxIter > maxLambdaIter_)
299 <<
"Number of lambda iterations exceeds maxLambdaIter(" 300 << maxLambdaIter_ <<
')'<<
endl;
307 template<
class BasicTurbulenceModel>
327 max(Fonset2 - Fonset3, scalar(0))
334 template<
class BasicTurbulenceModel>
343 const word& propertiesName,
414 this->coeffDict_.lookupOrDefault(
"lambdaErr", 1
e-6)
418 this->coeffDict_.lookupOrDefault(
"maxLambdaIter", 10)
427 this->runTime_.timeName(),
440 this->runTime_.timeName(),
453 this->runTime_.timeName(),
464 template<
class BasicTurbulenceModel>
469 ca1_.readIfPresent(this->coeffDict());
470 ca2_.readIfPresent(this->coeffDict());
471 ce1_.readIfPresent(this->coeffDict());
472 ce2_.readIfPresent(this->coeffDict());
473 sigmaThetat_.readIfPresent(this->coeffDict());
474 cThetat_.readIfPresent(this->coeffDict());
475 this->coeffDict().readIfPresent(
"lambdaErr", lambdaErr_);
476 this->coeffDict().readIfPresent(
"maxLambdaIter", maxLambdaIter_);
487 template<
class BasicTurbulenceModel>
516 alpha()*
rho()*(cThetat_/t)*(1 - Fthetat)
526 Pthetat*ReThetat0(Us, dUsds, nu) -
fvm::Sp(Pthetat, ReThetat_)
530 ReThetatEqn.
ref().relax();
545 *ca1_*Flength(nu)*S*
sqrt(gammaInt_()*Fonset(Rev, ReThetac, RT))
552 alpha()*
rho()*ca2_*Omega*Fturb*gammaInt_()
562 Pgamma -
fvm::Sp(ce1_*Pgamma, gammaInt_)
563 + Egamma -
fvm::Sp(ce2_*Egamma, gammaInt_)
567 gammaIntEqn.
ref().relax();
577 min(2*
max(Rev/(3.235*ReThetac) - 1, scalar(0))*Freattach, scalar(2))
581 gammaIntEff_ =
max(gammaInt_(), gammaSep);
585 template<
class BasicTurbulenceModel>
588 if (!this->turbulence_)
594 correctReThetatGammaInt();
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
BasicTurbulenceModel::transportModel transportModel
#define forAll(list, i)
Loop across all elements in list.
static dimensioned< Type > lookupOrAddToDict(const word &, dictionary &, const dimensionSet &dims=dimless, const Type &defaultValue=pTraits< Type >::zero)
Construct from dictionary, with default value.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
void clear() const
If object pointer points to valid object:
static tmp< DimensionedField< Type, GeoMesh > > New(const word &name, const Mesh &mesh, const dimensionSet &)
Return a temporary field constructed from name, mesh.
BasicTurbulenceModel::rhoField rhoField
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
tmp< volScalarField::Internal > ReThetat0(const volScalarField::Internal &Us, const volScalarField::Internal &dUsds, const volScalarField::Internal &nu) const
Return the transition onset momentum-thickness Reynolds number.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
kOmegaSSTLM(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName=turbulenceModel::propertiesName, const word &type=typeName)
Construct from components.
dimensionedTensor skew(const dimensionedTensor &dt)
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
T & ref() const
Return non-const reference or generate a fatal error.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
label k
Boltzmann constant.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
DimensionedField< scalar, volMesh > Internal
Type of the internal field from which this GeometricField is derived.
virtual tmp< volScalarField::Internal > epsilonByk(const volScalarField::Internal &F1, const volScalarField::Internal &F2) const
Modified form of the k-omega SST epsilon/k.
tmp< volScalarField::Internal > ReThetac() const
Empirical correlation for critical Reynolds number where the.
dimensionedScalar exp(const dimensionedScalar &ds)
A class for handling words, derived from string.
static word groupName(Name name, const word &group)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Specialisation for RAS of the generic kOmegaSSTBase base class. For more information, see Description of kOmegaSSTBase.H.
void constrain(fvMatrix< Type > &eqn)
Apply constraints to equation.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
virtual bool read()
Re-read model coefficients if they have changed.
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
tmp< volScalarField::Internal > Fonset(const volScalarField::Internal &Rev, const volScalarField::Internal &ReThetac, const volScalarField::Internal &RT) const
Transition onset location control function.
tmp< volScalarField::Internal > Flength(const volScalarField::Internal &nu) const
Empirical correlation that controls the length of the.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
virtual tmp< volScalarField::Internal > Pk(const volScalarField::Internal &G) const
Return k production rate.
tensor Ry(const scalar &omega)
Rotational transformation tensor about the y-axis by omega radians.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar pow3(const dimensionedScalar &ds)
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
#define WarningInFunction
Report a warning using Foam::Warning.
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
dimensionedScalar pow4(const dimensionedScalar &ds)
dimensionedScalar lambda(laminarTransport.lookup("lambda"))
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.
dimensioned< scalar > mag(const dimensioned< Type > &)
void correctReThetatGammaInt()
Solve the turbulence equations and correct the turbulence viscosity.
const doubleScalar e
Elementary charge.
virtual tmp< volScalarField::Internal > epsilonByk(const volScalarField::Internal &F1, const volScalarField::Internal &F2) const
Return epsilon/k which for standard RAS is betaStar*omega.
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
tmp< volScalarField::Internal > Fthetat(const volScalarField::Internal &Us, const volScalarField::Internal &Omega, const volScalarField::Internal &nu) const
Freestream blending-function.
A class for managing temporary objects.
static options & New(const fvMesh &mesh)
Construct fvOptions and register to datbase if not present.
BasicTurbulenceModel::alphaField alphaField
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
virtual tmp< volScalarField > F1(const volScalarField &CDkOmega) const
Modified form of the k-omega SST F1 function.
const dimensionSet dimVelocity
virtual tmp< volScalarField::Internal > Pk(const volScalarField::Internal &G) const
Modified form of the k-omega SST k production rate.