41namespace areaSurfaceFilmModels
72 const scalar rMin = 1
e-6;
81 const scalar rMax = 1e6;
84 if ((
mag(invR1[i]) < 1/rMax))
114 const label cellO = own[edgei];
115 const label cellN = nbr[edgei];
117 if (
phi[edgei] > phiMax[cellO])
119 phiMax[cellO] =
phi[edgei];
120 cosAngle[cellO] = -
gHat_ & UHat[cellN];
122 if (-
phi[edgei] > phiMax[cellN])
124 phiMax[cellN] = -
phi[edgei];
125 cosAngle[cellN] = -
gHat_ & UHat[cellO];
129 cosAngle *=
pos(invR1);
140 film().primaryMesh(),
151 return max(
min(cosAngle, scalar(1)), scalar(-1));
164 gradNHat_(
fac::grad(film.regionMesh().faceAreaNormals())),
165 deltaByR1Min_(coeffDict_.getOrDefault<scalar>(
"deltaByR1Min", 0)),
168 coeffDict_.getOrDefault<scalar>(
"definedPatchRadii", 0)
170 magG_(
mag(film.
g().value())),
174 coeffDict_.getOrDefault<scalar>(
"fThreshold", 1
e-8)
178 coeffDict_.getOrDefault<scalar>(
"minInvR1", 5)
181 if (
magG_ < ROOTVSMALL)
184 <<
"Acceleration due to gravity must be non-zero"
222 const scalar R1 = 1.0/(invR1[i] + ROOTVSMALL);
223 const scalar R2 = R1 +
delta[i];
226 const scalar Fi = -
delta[i]*
rho[i]*magSqrU[i]*72.0/60.0*invR1[i];
233 const scalar Fs = sigma[i]/R2;
235 Fnet[i] = Fi +
Fb + Fs;
245 massToInject = separated*availableMass;
246 diameterToInject = separated*
delta;
247 availableMass -= separated*availableMass;
259 film().primaryMesh(),
274 film().primaryMesh(),
281 areaSeparated.
write();
289 film().primaryMesh(),
296 areaMassToInject.
write();
304 film().primaryMesh(),
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
const uniformDimensionedVectorField & g
An indexed form of CGAL::Triangulation_face_base_2<K> used to keep track of the vertices in the trian...
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
void correctBoundaryConditions()
Correct boundary field.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
bool writeTime() const noexcept
True if this is a write time.
void size(const label n)
Older name for setAddressableSize.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
const Type & value() const
Return const reference to value.
Finite area mesh (used for 2-D non-Euclidian finite area method) defined using a patch of faces on a ...
const Time & time() const
Return the top-level database.
label nFaces() const noexcept
Number of mesh faces.
virtual bool write(const bool valid=true) const
Write using setting from DB.
Curvature film separation model.
areaTensorField gradNHat_
Gradient of surface normals.
scalar deltaByR1Min_
Minimum gravity driven film thickness (non-dimensionalised delta/R1)
scalar magG_
Magnitude of gravity vector.
tmp< areaScalarField > calcInvR1(const areaVectorField &U, const scalarField &calcCosAngle) const
Calculate local (inverse) radius of curvature.
tmp< scalarField > calcCosAngle(const edgeScalarField &phi) const
scalar minInvR1_
Minimum inv R1 for separation.
vector gHat_
Direction of gravity vector.
scalar definedPatchRadii_
scalar fThreshold_
Threshold force for separation.
const liquidFilmBase & film() const
Return const access to the film surface film model.
Base class for film injection models, handling mass transfer from the film.
void addToInjectedMass(const scalar dMass)
Add to injected mass.
virtual const areaScalarField & sigma() const =0
Access const reference sigma.
const areaVectorField & Uf() const
Access const reference Uf.
const areaScalarField & h() const
Access const reference h.
virtual const areaScalarField & rho() const =0
Access const reference rho.
const edgeScalarField & phi2s() const
Access continuity flux.
const uniformDimensionedVectorField & g() const
Gravity.
const faMesh & regionMesh() const
Return the region mesh database.
A class for managing temporary objects.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
autoPtr< surfaceVectorField > Uf
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedScalar pos(const dimensionedScalar &ds)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
const dimensionSet dimVelocity
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
GeometricField< scalar, faPatchField, areaMesh > areaScalarField
const dimensionSet dimForce
static constexpr const zero Zero
Global zero (0)
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Calculate the second temporal derivative.
#define forAll(list, i)
Loop across all elements in list.
Operations on lists of strings.