117#include "surfaceInterpolate.H"
135 public surfaceInterpolationScheme<Type>,
136 public blendedSchemeBase<Type>
141 tmp<surfaceInterpolationScheme<Type>> tScheme1_;
176 void operator=(
const DEShybrid&) =
delete;
197 CH3_*Omega*
max(S, Omega)
212 nuEff/(
pow(0.09, 1.5)*
K),
225 CDES_*
delta/
max(lTurb*
g, SMALL*L0_) - 0.5
242 if (blendedSchemeBaseName::debug)
259 auto& factor = *factorPtr;
261 factor =
max(sigmaMax_*
tanh(
pow(
A, CH1_)), sigmaMin_);
265 vf.
name() +
"BlendingFactor",
279 vf.
name() +
"BlendingFactor",
309 CDES_(readScalar(is)),
312 sigmaMin_(readScalar(is)),
313 sigmaMax_(readScalar(is)),
314 OmegaLim_(readScalar(is)),
319 if (U0_.
value() <= 0)
322 <<
"U0 coefficient must be > 0. "
325 if (L0_.
value() <= 0)
328 <<
"L0 coefficient must be > 0. "
334 <<
"sigmaMin coefficient must be >= 0. "
340 <<
"sigmaMax coefficient must be >= 0. "
346 <<
"sigmaMin coefficient must be <= 1. "
352 <<
"sigmaMax coefficient must be <= 1. "
375 CDES_(readScalar(is)),
378 sigmaMin_(readScalar(is)),
379 sigmaMax_(readScalar(is)),
380 OmegaLim_(readScalar(is)),
385 if (U0_.
value() <= 0)
388 <<
"U0 coefficient must be > 0. "
391 if (L0_.
value() <= 0)
394 <<
"L0 coefficient must be > 0. "
400 <<
"sigmaMin coefficient must be >= 0. "
406 <<
"sigmaMax coefficient must be >= 0. "
412 <<
"sigmaMin coefficient must be <= 1. "
418 <<
"sigmaMax coefficient must be <= 1. "
439 lookupObject<const volScalarField>(deltaName_);
448 const icoModel& model =
451 return calcBlendingFactor(vf, model.nuEff(), model.U(),
delta);
455 const cmpModel& model =
458 return calcBlendingFactor
460 vf, model.muEff()/model.rho(), model.U(),
delta
465 <<
"Scheme requires a turbulence model to be present. "
466 <<
"Unable to retrieve turbulence model from the mesh "
474 tmp<surfaceScalarField>
weights
476 const GeometricField<Type, fvPatchField, volMesh>& vf
482 (scalar(1) - bf)*tScheme1_().weights(vf)
483 + bf*tScheme2_().weights(vf);
498 (scalar(1) - bf)*tScheme1_().interpolate(vf)
499 + bf*tScheme2_().interpolate(vf);
506 return tScheme1_().corrected() || tScheme2_().corrected();
static const Foam::dimensionedScalar A("", Foam::dimPressure, 611.21)
static const Foam::dimensionedScalar B("", Foam::dimless, 18.678)
CGAL::Exact_predicates_exact_constructions_kernel K
const uniformDimensionedVectorField & g
Hybrid convection scheme of Travin et al. for hybrid RAS/LES calculations.
virtual bool corrected() const
Return true if this scheme uses an explicit correction.
tmp< surfaceScalarField > weights(const GeometricField< Type, fvPatchField, volMesh > &vf) const
Return the interpolation weighting factors.
virtual tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > correction(const GeometricField< Type, fvPatchField, volMesh > &vf) const
Return the explicit correction to the face-interpolate.
DEShybrid(const fvMesh &mesh, const surfaceScalarField &faceFlux, Istream &is)
Construct from mesh, faceFlux and Istream.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &vf) const
Return the face-interpolate of the given cell field.
virtual tmp< surfaceScalarField > blendingFactor(const GeometricField< Type, fvPatchField, volMesh > &vf) const
Return the face-based blending factor.
TypeName("DEShybrid")
Runtime type information.
Generic GeometricField class.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const word & name() const noexcept
Return the object name.
static word scopedName(const std::string &scope, const word &name)
Create scope:name or scope_name string.
Templated abstract base class for single-phase incompressible turbulence models.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
static word timeName(const scalar t, const int precision=precision_)
Base class for blended schemes to provide access to the blending factor surface field.
const Type & value() const
Return const reference to value.
Mesh data needed to do the Finite Volume discretisation.
const Time & time() const
Return the top-level database.
bool foundObject(const word &name, const bool recursive=false) const
Is the named Type found?
const Type & lookupObject(const word &name, const bool recursive=false) const
Type * getObjectPtr(const word &name, const bool recursive=false) const
Abstract base class for surface interpolation schemes.
const fvMesh & mesh() const
Return mesh reference.
static tmp< surfaceInterpolationScheme< Type > > New(const fvMesh &mesh, Istream &schemeData)
Return new tmp interpolation scheme.
A class for managing temporary objects.
A class for handling words, derived from Foam::string.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Calculate the gradient of the given field.
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
dimensionedScalar tanh(const dimensionedScalar &ds)
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
dimensionedScalar pow4(const dimensionedScalar &ds)
static constexpr const zero Zero
Global zero (0)
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensionedTensor skew(const dimensionedTensor &dt)
#define TypeName(TypeNameString)
Declare a ClassName() with extra virtual type info.