Go to the documentation of this file.
41 namespace regionModels
43 namespace surfaceFilmModels
60 thixotropicViscosity::thixotropicViscosity
71 c_(
"c",
pow(
dimTime, d_.value() - scalar(1)), coeffDict_),
73 muInf_(
"muInf", mu0_.dimensions(), coeffDict_),
74 K_(1 -
sqrt(muInf_/mu0_)),
93 mu_.correctBoundaryConditions();
136 deltaRho.dimensions(),
151 a_*
pow((1 - lambda_), b_)
161 )/(deltaRho + deltaRho0),
172 mu_ = muInf_/(
sqr(1 - K_*lambda_) + ROOTVSMALL);
173 mu_.correctBoundaryConditions();
int debug
Static debugging option.
const dimensionSet dimPressure
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
const dimensionedScalar mu
Atomic mass unit.
static constexpr const zero Zero
Global zero (0)
virtual ~thixotropicViscosity()
Destructor.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
virtual const volVectorField & U() const
Return the film velocity [m/s].
static word timeName(const scalar t, const int precision=precision_)
Calculate the divergence of the given field.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Kinematic form of single-cell layer surface film model.
volScalarField & rhoSp()
Mass [kg/m2/s].
addToRunTimeSelectionTable(surfaceFilmRegionModel, kinematicSingleLayer, mesh)
const dimensionedScalar & deltaSmall() const
Return small delta.
Calculate the matrix for the divergence of the given field and flux.
Base class for surface film viscosity models.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
Base class for surface film models.
virtual void correct(const volScalarField &p, const volScalarField &T)
Correct.
void relax(const scalar alpha)
Relax matrix (for steady-state solution).
const fvMesh & regionMesh() const
Return the region mesh database.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
virtual const surfaceScalarField & phi() const
Return the film flux [kg.m/s].
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Macros for easy insertion into run-time selection tables.
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.
Calculate the matrix for implicit and explicit sources.
defineTypeNameAndDebug(kinematicSingleLayer, 0)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
virtual const volVectorField & Uw() const
Return the film wall velocity [m/s].
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const volScalarField & delta() const
Return const access to the film thickness [m].
const volScalarField & alpha() const
Return the film coverage, 1 = covered, 0 = uncovered [].
A special matrix type and solver, designed for finite volume solutions of scalar equations....
const Time & time() const
Return the top-level database.
virtual const volScalarField & deltaRho() const
Return the film thickness*density (helper field) [kg/m3].
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Calculate the matrix for the first temporal derivative.