Go to the documentation of this file.
34 #include "surfaceInterpolate.H"
44 namespace regionModels
46 namespace surfaceFilmModels
86 const scalar rMin = 1
e-6;
89 forAll(definedPatchRadii_, i)
91 label patchi = definedPatchRadii_[i].first();
92 scalar definedInvR1 = 1.0/
max(rMin, definedPatchRadii_[i].second());
97 const scalar rMax = 1e6;
100 if (
mag(invR1[i]) < 1/rMax)
129 label cellO = own[facei];
130 label cellN = nbr[facei];
132 if (
phi[facei] > phiMax[cellO])
134 phiMax[cellO] =
phi[facei];
135 cosAngle[cellO] = -gHat_ & nf[facei];
137 if (-
phi[facei] > phiMax[cellN])
139 phiMax[cellN] = -
phi[facei];
140 cosAngle[cellN] = -gHat_ & -nf[facei];
153 if (phip[i] > phiMax[celli])
155 phiMax[celli] = phip[i];
156 cosAngle[celli] = -gHat_ & nf[i];
204 mesh.time().timeName(),
210 zeroGradientFvPatchScalarField::typeName
213 volCosAngle.correctBoundaryConditions();
217 return max(
min(cosAngle, scalar(1)), scalar(-1));
223 curvatureSeparation::curvatureSeparation
231 deltaByR1Min_(coeffDict_.lookupOrDefault<scalar>(
"deltaByR1Min", 0)),
232 definedPatchRadii_(),
236 if (magG_ < ROOTVSMALL)
239 <<
"Acceleration due to gravity must be non-zero"
243 gHat_ = film.
g().
value()/magG_;
254 labelList patchIDs = findIndices(allPatchNames, prIn[i].first());
257 const label patchi = patchIDs[j];
259 if (!uniquePatchIDs.found(patchi))
261 const scalar radius = prIn[i].second();
264 uniquePatchIDs.insert(patchi);
269 definedPatchRadii_.
transfer(prData);
289 refCast<const kinematicSingleLayer>(this->film());
303 const scalar Fthreshold = 1
e-10;
308 if ((invR1[i] > 0) && (
delta[i]*invR1[i] > deltaByR1Min_))
310 scalar R1 = 1.0/(invR1[i] + ROOTVSMALL);
311 scalar R2 = R1 +
delta[i];
314 scalar Fi = -
delta[i]*
rho[i]*magSqrU[i]*72.0/60.0*invR1[i];
318 - 0.5*
rho[i]*magG_*invR1[i]*(
sqr(R1) -
sqr(R2))*cosAngle[i];
321 scalar Fs =
sigma[i]/R2;
323 Fnet[i] = Fi +
Fb + Fs;
325 if (Fnet[i] + Fthreshold < 0)
333 massToInject = separated*availableMass;
334 diameterToInject = separated*
delta;
335 availableMass -= separated*availableMass;
337 addToInjectedMass(
sum(separated*availableMass));
346 mesh.time().timeName(),
352 zeroGradientFvPatchScalarField::typeName
355 volFnet.correctBoundaryConditions();
int debug
Static debugging option.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
virtual ~curvatureSeparation()
Destructor.
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
static constexpr const zero Zero
Global zero.
const dimensionSet dimVelocity
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
virtual const volVectorField & U() const
Return the film velocity [m/s].
tmp< volScalarField > calcInvR1(const volVectorField &U) const
Calculate local (inverse) radius of curvature.
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Calculate the divergence of the given field.
Kinematic form of single-cell layer surface film model.
void append(const T &val)
Append an element at the end of the list.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
const Type & value() const
Return const reference to value.
const dimensionSet dimForce
addToRunTimeSelectionTable(surfaceFilmRegionModel, kinematicSingleLayer, mesh)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
wordList names() const
Return a list of patch names.
const dimensionedVector & g() const
Return the acceleration due to gravity.
#define forAll(list, i)
Loop across all elements in list.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
An indexed form of CGAL::Triangulation_face_base_2<K> used to keep track of the vertices in the trian...
tmp< vectorField > nf() const
Return face normals.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Base class for surface film models.
const fvMesh & regionMesh() const
Return the region mesh database.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
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,...
Base class for film injection models, handling mass transfer from the film.
Macros for easy insertion into run-time selection tables.
Mesh data needed to do the Finite Volume discretisation.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
const volScalarField & sigma() const
Return const access to the surface tension [kg/s2].
defineTypeNameAndDebug(kinematicSingleLayer, 0)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const fvPatch & patch() const
Return patch.
Calculate the gradient of the given field.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
virtual const volVectorField & nHat() const
Return the patch normal vectors.
const volScalarField & delta() const
Return const access to the film thickness [m].
void transfer(List< T > &lst)
Transfer contents of the argument List into this.
const dimensionedScalar e
Elementary charge.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
virtual const labelUList & faceCells() const
Return faceCells.
#define forAllReverse(list, i)
Reverse loop across all elements in list.
A List with indirect addressing.
Operations on lists of strings.
dimensionedScalar sigma("sigma", dimMass/sqr(dimTime), transportProperties)
A 2-tuple for storing two objects of dissimilar types. The container is similar in purpose to std::pa...
virtual const volScalarField & rho() const
Return the film density [kg/m3].
Smooth ATC in cells next to a set of patches supplied by type.
tmp< scalarField > calcCosAngle(const surfaceScalarField &phi) const
Calculate the cosine of the angle between gravity vector and.