Go to the documentation of this file.
37 namespace regionModels
39 namespace areaSurfaceFilmModels
53 { frictionMethodType::mquadraticProfile,
"quadraticProfile" },
54 { frictionMethodType::mlinearProfile,
"linearProfile" },
55 { frictionMethodType::mDarcyWeisbach,
"DarcyWeisbach" },
56 { frictionMethodType::mManningStrickler,
"ManningStrickler" }
66 { shearMethodType::msimple,
"simple" },
67 { shearMethodType::mwallFunction,
"wallFunction" }
73 filmTurbulenceModel::filmTurbulenceModel
75 const word& modelType,
82 method_(frictionMethodTypeNames_.get(
"friction", dict_)),
83 shearMethod_(shearMethodTypeNames_.get(
"shearStress", dict_)),
84 rhoName_(dict_.getOrDefault<
word>(
"rho",
"rho")),
87 if (rhoName_ ==
"rhoInf")
89 rhoRef_ = dict_.
get<scalar>(
"rhoInf");
116 auto&
Cw = tCw.ref();
153 const scalar Cf =
dict_.
get<scalar>(
"DarcyWeisbach");
170 Cw.primitiveFieldRef() =
178 <<
"Unimplemented method "
180 <<
"Please set 'frictionMethod' to one of "
202 switch (shearMethod_)
212 dict_.get<scalar>(
"Cf")
223 const volSymmTensorField::Boundary& devRhoReffb
224 = tdevRhoReff().boundaryField();
226 const label patchi = film_.patchID();
228 const surfaceVectorField::Boundary& Sfb =
229 film_.primaryMesh().Sf().boundaryField();
234 film_.regionMesh().faceAreaNormals().internalField();
237 fT -= nHat*(fT & nHat);
244 film_.primaryMesh().time().timeName(),
250 vectorField& aForce = taForce.ref().primitiveFieldRef();
253 const vectorField afT(film_.vsm().mapToSurface(fT));
256 film_.regionMesh().S();
258 aForce = afT/(film_.rho().primitiveField()*magSf);
260 tshearStress.
ref() += taForce();
262 if (film_.regionMesh().time().writeTime())
284 if (m.
foundObject<cmpTurbModel>(cmpTurbModel::propertiesName))
287 m.
lookupObject<cmpTurbModel>(cmpTurbModel::propertiesName);
289 return turb.devRhoReff();
291 else if (m.
foundObject<icoTurbModel>(icoTurbModel::propertiesName))
294 m.
lookupObject<icoTurbModel>(icoTurbModel::propertiesName);
307 const auto& laminarT =
324 <<
"No valid model for viscous stress calculation"
Defines the attributes of an object for which implicit objectRegistry management is supported,...
shearMethodType
Options for the shear stress models.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
A class for handling words, derived from Foam::string.
faMatrix< vector > faVectorMatrix
word rhoName_
Name of density field (optional)
const dimensionedScalar mu
Atomic mass unit.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
const fvMesh & primaryMesh() const
Return the reference to the primary mesh database.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
const dimensionSet dimVelocity
const dimensionSet dimDensity
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
IOdictionary transportProperties(IOobject("transportProperties", runTime.constant(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
static word timeName(const scalar t, const int precision=precision_)
autoPtr< surfaceVectorField > Uf
frictionMethodType
Options for the friction models.
Fundamental fluid thermodynamic properties.
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
tmp< volScalarField > rho() const
Return rho if specified otherwise rhoRef.
word UName() const
Name of the U field.
const Type & value() const
Return const reference to value.
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
EnumType get(const word &enumName) const
The enumeration corresponding to the given name.
const frictionMethodType method_
Method used.
bool foundObject(const word &name, const bool recursive=false) const
Is the named Type found?
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
const areaVectorField & Uf() const
Access const reference Uf.
scalar rhoRef_
Reference density needed for incompressible calculations.
static const Enum< frictionMethodType > frictionMethodTypeNames_
Names for friction models.
tmp< faMatrix< Type > > Sp(const areaScalarField &sp, const GeometricField< Type, faPatchField, areaMesh > &vf)
const liquidFilmBase & film_
Reference to liquidFilmBase.
const dictionary dict_
Model dictionary.
const dimensionedScalar h
Planck constant.
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
defineRunTimeSelectionTable(liquidFilmBase, dictionary)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual tmp< areaScalarField > Cw() const
Return the wall film surface friction.
const Type & lookupObject(const word &name, const bool recursive=false) const
const dimensionSet dimViscosity
const faMesh & regionMesh() const
Return the region mesh database.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
tmp< faVectorMatrix > primaryRegionFriction(areaVectorField &U) const
Return primary region friction.
Mesh data needed to do the Finite Volume discretisation.
const uniformDimensionedVectorField & g
FlatOutput::OutputAdaptor< Container, Delimiters > flatOutput(const Container &obj, Delimiters delim)
Global flatOutput() function with specified output delimiters.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Base-class for all transport models used by the incompressible turbulence models.
tmp< volSymmTensorField > devRhoReff() const
Return the effective viscous stress (laminar + turbulent)
virtual const areaScalarField & rho() const =0
Access const reference rho.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
static const Enum< shearMethodType > shearMethodTypeNames_
Names for shear stress models.
defineTypeNameAndDebug(kinematicThinFilm, 0)
const liquidFilmBase & film() const
Return film.
const dimensionedScalar & h0() const
Return h0.
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
static tmp< T > New(Args &&... args)
Construct tmp of T with forwarding arguments.
const areaScalarField & h() const
Access const reference h.
const Time & time() const
Return the top-level database.
List< word > sortedToc() const
The sorted list of enum names.
Templated abstract base class for single-phase incompressible turbulence models.
static const GeometricField< symmTensor, fvPatchField, volMesh > & null()
Return a null geometric field.
dimensionedScalar cbrt(const dimensionedScalar &ds)
compressible::turbulenceModel & turb
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
static const gravity & New(const Time &runTime)
Construct on Time.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
virtual const areaScalarField & mu() const =0
Access const reference mu.
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
word dictName() const
The local dictionary name (final part of scoped name)