Go to the documentation of this file.
39 Foam::reconstructedDistanceFunction::coupledFacesPatch()
const
45 for (
const polyPatch& pp :
patches)
47 if (isA<coupledPolyPatch>(pp))
49 nCoupled += pp.size();
55 for (
const polyPatch& pp :
patches)
57 if (isA<coupledPolyPatch>(pp))
59 label facei = pp.start();
63 nCoupledFaces[nCoupled++] = facei++;
83 const label neiRingLevel
88 if (mesh_.topoChanging())
91 if (nextToInterface_.size() != mesh_.nCells())
93 nextToInterface_.resize(mesh_.nCells());
95 coupledBoundaryPoints_ = coupledFacesPatch()().meshPoints();
101 boolList alreadyMarkedPoint(mesh_.nPoints(),
false);
102 nextToInterface_ =
false;
107 for (
int level=0;level<=neiRingLevel;level++)
112 forAll(coupledBoundaryPoints_, i)
114 const label
pi = coupledBoundaryPoints_[i];
117 const label celli = cPoints[
pi][j];
118 if (cellDistLevel_[celli] == level-1)
120 syncMap.insert(
pi,
true);
131 const label
pi = iter.key();
133 if (!alreadyMarkedPoint[
pi])
138 const label pCelli = cPoints[
pi][j];
139 if (cellDistLevel_[pCelli] == -1)
141 cellDistLevel_[pCelli] = level;
142 nextToInterface_[pCelli] =
true;
146 alreadyMarkedPoint[
pi] =
true;
151 forAll(cellDistLevel_, celli)
155 if (interfaceCells[celli])
157 cellDistLevel_[celli] = 0;
158 nextToInterface_[celli] =
true;
162 cellDistLevel_[celli] = -1;
167 if (cellDistLevel_[celli] == level-1)
171 const label pI = pCells[celli][i];
173 if (!alreadyMarkedPoint[pI])
177 const label pCelli = cPoints[pI][j];
178 if (cellDistLevel_[pCelli] == -1)
180 cellDistLevel_[pCelli] = level;
181 nextToInterface_[pCelli] =
true;
185 alreadyMarkedPoint[pI] =
true;
215 coupledBoundaryPoints_(coupledFacesPatch()().meshPoints()),
246 if (nextToInterface.size() != centre.size())
249 <<
"size of nextToInterface: " << nextToInterface.size()
250 <<
"size of centre:" << centre.size()
251 <<
"do not match. Did the mesh change?"
253 return reconDistFunc;
267 forAll(nextToInterface,celli)
269 if (nextToInterface[celli])
271 if (
mag(normal[celli]) != 0)
274 scalar dist = (centre[celli] - mesh_.C()[celli]) &
n;
275 reconDistFunc[celli] = dist;
279 scalar averageDist = 0;
280 scalar avgWeight = 0;
281 const point p = mesh_.C()[celli];
285 const label gblIdx = stencil[celli][i];
292 scalar distToSurf = distanceToIntSeg & (
n);
295 if (
mag(distanceToIntSeg) != 0)
297 distanceToIntSeg /=
mag(distanceToIntSeg);
298 weight =
sqr(
mag(distanceToIntSeg &
n));
304 averageDist += distToSurf * weight;
311 reconDistFunc[celli] = averageDist / avgWeight;
317 reconDistFunc[celli] = 0;
324 if (isA<calculatedFvPatchScalarField>(pRDF))
331 if (nextToInterface_[pCellI])
333 scalar averageDist = 0;
334 scalar avgWeight = 0;
337 forAll(stencil[pCellI], j)
339 const label gblIdx = stencil[pCellI][j];
345 distribute.
getValue(centre, mapCentres, gblIdx);
347 scalar distToSurf = distanceToIntSeg & (
n);
350 if (
mag(distanceToIntSeg) != 0)
352 distanceToIntSeg /=
mag(distanceToIntSeg);
353 weight =
sqr(
mag(distanceToIntSeg &
n));
359 averageDist += distToSurf * weight;
366 pRDF[i] = averageDist / avgWeight;
383 return reconDistFunc;
391 surfaceVectorField::Boundary& nHatb
395 const volScalarField::Boundary& abf =
alpha.boundaryField();
396 volScalarField::Boundary& RDFbf = this->boundaryFieldRef();
402 if (isA<alphaContactAngleTwoPhaseFvPatchScalarField>(abf[patchi]))
407 refCast<const alphaContactAngleTwoPhaseFvPatchScalarField>
420 1/acap.patch().deltaCoeffs()*
cos(theta)
421 + RDFbf[patchi].patchInternalField();
Type getValue(const GeometricField< Type, fvPatchField, volMesh > &phi, const Map< Type > &valuesFromOtherProc, const label gblIdx) const
static autoPtr< T > New(Args &&... args)
Construct autoPtr of T with forwarding arguments.
List< label > labelList
A List of labels.
virtual const pointField & points() const
Return raw points.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Map< Type > getDatafromOtherProc(const boolList &zone, const GeometricField< Type, fvPatchField, volMesh > &phi)
Returns stencil and provides a Map with globalNumbering.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
const volScalarField & constructRDF(const boolList &nextToInterface, const volVectorField ¢re, const volVectorField &normal, zoneDistribute &distribute, bool updateStencil=true)
static constexpr const zero Zero
Global zero (0)
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
static word timeName(const scalar t, const int precision=precision_)
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
A HashTable to objects of type <T> with a label key.
Unit conversion functions.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
void updateContactAngle(const volScalarField &alpha, const volVectorField &U, surfaceVectorField::Boundary &nHatb)
#define forAll(list, i)
Loop across all elements in list.
label nCells() const
Number of mesh cells.
A patch is a list of labels that address the faces in the global face list.
const labelListList & getStencil()
Stencil reference.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
void markCellsNearSurf(const boolList &interfaceCells, const label neiRingLevel)
Mesh data needed to do the Finite Volume discretisation.
const labelUList & faceCells() const
Return face-cell addressing.
constexpr scalar degToRad(const scalar deg) noexcept
Conversion from degrees to radians.
errorManipArg< error, int > exit(error &err, const int errNo=1)
void setUpCommforZone(const boolList &zone, bool updateStencil=true)
Update stencil with boolList the size has to match mesh nCells.
void correctBoundaryConditions()
Correct boundary field.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
virtual const faceList & faces() const
Return raw faces.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
constexpr scalar pi(M_PI)
forAllConstIters(mixture.phases(), phase)
Class for parallel communication in a narrow band. It either provides a Map with the neighbouring val...
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const polyBoundaryMesh & patches
const polyPatch & patch() const
Return the polyPatch.
const dimensionedScalar c
Speed of light in a vacuum.
const Time & time() const
Return the top-level database.
const fvPatch & patch() const
Return patch.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
reconstructedDistanceFunction(const fvMesh &mesh)
Construct from fvMesh.
dimensionedScalar cos(const dimensionedScalar &ds)