37 namespace regionModels
39 namespace surfaceFilmModels
49 void contactAngleForce::initialise()
53 if (zeroForcePatches.size())
58 Info<<
" Assigning zero contact force within " << dLim
59 <<
" of patches:" <<
endl;
65 label patchi = iter.key();
94 contactAngleForce::contactAngleForce
100 force(typeName, owner, dict),
102 rndGen_(
label(0), -1),
115 typeName +
":contactForceMask",
145 typeName +
":contactForce",
170 const label cellO = own[facei];
171 const label cellN = nbr[facei];
174 if ((alpha[cellO] > 0.5) && (alpha[cellN] < 0.5))
178 else if ((alpha[cellO] < 0.5) && (alpha[cellN] > 0.5))
183 if (celli != -1 && mask_[celli] > 0.5)
187 gradAlpha[celli]/(
mag(gradAlpha[celli]) + ROOTVSMALL);
188 scalar theta =
cos(
degToRad(distribution_->sample()));
189 force[celli] += Ccf_*n*sigma[celli]*(1.0 - theta)/invDx;
204 if (maskf[facei] > 0.5)
206 label cellO = faceCells[facei];
208 if ((alpha[cellO] > 0.5) && (alphaf[facei] < 0.5))
212 /(
mag(gradAlpha[cellO]) + ROOTVSMALL);
213 scalar theta =
cos(
degToRad(distribution_->sample()));
215 Ccf_*n*sigma[cellO]*(1.0 - theta)/invDx[facei];
232 tfvm.
ref() += tForce;
const fvPatch & patch() const
Return patch.
virtual const volScalarField & magSf() const
Return the face area magnitudes / [m2].
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
const Time & time() const
Return the reference to the time database.
Base class for film (stress-based) force models.
A list of keyword definitions, which are a keyword followed by any number of values (e...
surfaceFilmModel & owner_
Reference to the owner surface film model.
Unit conversion functions.
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
Ostream & endl(Ostream &os)
Add newline and flush stream.
bool isCoupledPatch(const label regionPatchi) const
Return true if patchi on the local region is a coupled.
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
bool writeTime() const
Return true if this is a write time.
const dictionary coeffDict_
Coefficients dictionary.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Macros for easy insertion into run-time selection tables.
scalar degToRad(const scalar deg)
Conversion from degrees to radians.
const fvMesh & regionMesh() const
Return the region mesh database.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys.
dimensionedScalar pos(const dimensionedScalar &ds)
const dimensionSet dimVolume(pow3(dimLength))
static autoPtr< distributionModel > New(const dictionary &dict, cachedRandom &rndGen)
Selector.
addToRunTimeSelectionTable(surfaceFilmModel, kinematicSingleLayer, mesh)
dimensionedScalar cos(const dimensionedScalar &ds)
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const labelUList & neighbour() const
Internal face neighbour.
Calculate the gradient of the given field.
virtual const volScalarField & alpha() const =0
Return the film coverage, 1 = covered, 0 = uncovered / [].
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
virtual const volScalarField & sigma() const =0
Return the film surface tension [N/m].
const dimensionSet dimForce
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const surfaceScalarField & deltaCoeffs() const
Return reference to cell-centre difference coefficients.
const surfaceFilmModel & owner() const
Return const access to the owner surface film model.
fvMatrix< vector > fvVectorMatrix
Base class for surface film models.
dimensioned< scalar > mag(const dimensioned< Type > &)
List< wordRe > wordReList
A List of wordRe (word or regular expression)
const labelUList & owner() const
Internal face owner.
A class for managing temporary objects.
T & ref() const
Return non-const reference or generate a fatal error.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
defineTypeNameAndDebug(kinematicSingleLayer, 0)
const dimensionSet dimArea(sqr(dimLength))
const Time & time() const
Return the top-level database.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.