76void Foam::interfaceTrackingFvMesh::initializeData()
80 const word fsPatchName(
motion().get<word>(
"fsPatchName"));
87 <<
"Patch name " << fsPatchName <<
" not found."
91 fsPatchIndex_ =
patch.index();
98 for (
const word& patchName : pointNormalsCorrectionPatches_)
105 <<
"Patch name '" << patchName
106 <<
"' for point normals correction does not exist"
115 if (!normalMotionDir_)
125 "nonReflectingFreeSurfacePatches",
126 nonReflectingFreeSurfacePatches_
131void Foam::interfaceTrackingFvMesh::makeUs()
const
134 <<
"making free-surface velocity field" <<
nl;
139 <<
"free-surface velocity field already exists"
146 zeroGradientFaPatchVectorField::typeName
157 patchFieldTypes[patchI] =
158 wedgeFaPatchVectorField::typeName;
162 label ngbPolyPatchID =
165 if (ngbPolyPatchID != -1)
173 patchFieldTypes[patchI] =
174 slipFaPatchVectorField::typeName;
180 for (
const word& patchName : fixedFreeSurfacePatches_)
182 const label fixedPatchID =
185 if (fixedPatchID == -1)
188 <<
"Wrong faPatch name '" << patchName
189 <<
"' in the fixedFreeSurfacePatches list"
190 <<
" defined in the dynamicMeshDict dictionary"
194 label ngbPolyPatchID =
197 if (ngbPolyPatchID != -1)
205 patchFieldTypes[fixedPatchID] =
206 fixedValueFaPatchVectorField::typeName;
226 for (
const word& patchName : fixedFreeSurfacePatches_)
230 if (fixedPatchID == -1)
233 <<
"Wrong faPatch name '" << patchName
234 <<
"' in the fixedFreeSurfacePatches list"
235 <<
" defined in the dynamicMeshDict dictionary"
239 label ngbPolyPatchID =
242 if (ngbPolyPatchID != -1)
250 UsPtr_->boundaryFieldRef()[fixedPatchID] ==
Zero;
257void Foam::interfaceTrackingFvMesh::makeFsNetPhi()
const
260 <<
"making free-surface net flux" <<
nl;
265 <<
"free-surface net flux already exists"
285void Foam::interfaceTrackingFvMesh::makeControlPoints()
288 <<
"making control points" <<
nl;
290 if (controlPointsPtr_)
293 <<
"control points already exists"
307 Info<<
"Reading control points" <<
endl;
323 Info<<
"Creating new control points" <<
endl;
335 aMesh().areaCentres().internalField()
338 initializeControlPointsPosition();
343void Foam::interfaceTrackingFvMesh::makeMotionPointsMask()
346 <<
"making motion points mask" <<
nl;
348 if (motionPointsMaskPtr_)
351 <<
"motion points mask already exists"
374 forAll(patchPoints, pointI)
376 motionPointsMask()[patchPoints[pointI]] = -1;
382 for (
const word& patchName : fixedFreeSurfacePatches_)
386 if (fixedPatchID == -1)
389 <<
"Wrong faPatch name in the fixedFreeSurfacePatches list"
390 <<
" defined in the dynamicMeshDict dictionary"
397 forAll(patchPoints, pointI)
399 motionPointsMask()[patchPoints[pointI]] = 0;
405void Foam::interfaceTrackingFvMesh::makeDirections()
408 <<
"make displacement directions for points and control points" <<
nl;
410 if (pointsDisplacementDirPtr_ || facesDisplacementDirPtr_)
413 <<
"points, control points displacement directions already exist"
417 pointsDisplacementDirPtr_ =
424 facesDisplacementDirPtr_ =
431 if (!normalMotionDir())
433 if (
mag(motionDir_) < SMALL)
436 <<
"Zero motion direction"
440 facesDisplacementDir() = motionDir_;
441 pointsDisplacementDir() = motionDir_;
444 updateDisplacementDirections();
448void Foam::interfaceTrackingFvMesh::makePhis()
451 <<
"making free-surface flux" <<
nl;
456 <<
"free-surface flux already exists"
475void Foam::interfaceTrackingFvMesh::makeSurfactConc()
const
478 <<
"making free-surface surfactant concentration field" <<
nl;
483 <<
"free-surface surfactant concentration field already exists"
506void Foam::interfaceTrackingFvMesh::makeBulkSurfactConc()
const
509 <<
"making volume surfactant concentration field" <<
nl;
511 if (bulkSurfactConcPtr_)
514 <<
"volume surfactant concentration field already exists"
539 bulkSurfactConc = surfactant().bulkConc();
545void Foam::interfaceTrackingFvMesh::makeSurfaceTension()
const
548 <<
"making surface tension field" <<
nl;
550 if (surfaceTensionPtr_)
553 <<
"surface tension field already exists"
567 sigma() + surfactant().dSigma(surfactantConcentration())/rho_
572void Foam::interfaceTrackingFvMesh::makeSurfactant()
const
575 <<
"making surfactant properties" <<
nl;
580 <<
"surfactant properties already exists"
585 motion().
subDict(
"surfactantProperties");
591void Foam::interfaceTrackingFvMesh::makeContactAngle()
594 <<
"making contact angle field" <<
nl;
596 if (contactAnglePtr_)
599 <<
"contact angle already exists"
614 Info<<
"Reading contact angle field" <<
endl;
633void Foam::interfaceTrackingFvMesh::updateDisplacementDirections()
635 if (normalMotionDir())
643 if (contactAnglePtr_)
645 label ngbPolyPatchID =
648 if (ngbPolyPatchID != -1)
662 .ngbPolyPatchPointNormals()
665 forAll(patchPoints, pointI)
667 pointsDisplacementDir()[patchPoints[pointI]] -=
671 & pointsDisplacementDir()[patchPoints[pointI]]
674 pointsDisplacementDir()[patchPoints[pointI]] /=
677 pointsDisplacementDir()
689 facesDisplacementDir() =
696 facesDisplacementDir()
697 *(facesDisplacementDir()&(controlPoints() - Cf))
703void Foam::interfaceTrackingFvMesh::initializeControlPointsPosition()
711 correctPointDisplacement(sweptVolCorr, displacement);
719 sweptVol[faceI] = -faces[faceI].sweptVol(
points, newPoints);
726 faceArea[faceI] = faces[faceI].unitNormal(newPoints);
733 deltaH[faceI] = sweptVol[faceI]/
734 ((faceArea[faceI] & facesDisplacementDir()[faceI]) + SMALL);
736 if (
mag(faceArea[faceI] & facesDisplacementDir()[faceI]) < SMALL)
743 <<
"Something is wrong with specified motion direction"
748 for (
const word& patchName : fixedFreeSurfacePatches_)
752 if (fixedPatchID == -1)
755 <<
"Wrong faPatch name in the fixedFreeSurfacePatches list"
756 <<
" defined in the freeSurfaceProperties dictionary"
765 deltaH[eFaces[edgeI]] *= 2.0;
769 controlPoints() += facesDisplacementDir()*deltaH;
774void Foam::interfaceTrackingFvMesh::smoothFreeSurfaceMesh()
776 Info<<
"Smoothing free surface mesh" <<
endl;
790 sweptVol[faceI] = -faces[faceI].sweptVol(
points, newPoints);
796 faceArea[faceI] = faces[faceI].unitNormal(newPoints);
801 sweptVol/(faceArea & facesDisplacementDir())
804 for (
const word& patchName : fixedFreeSurfacePatches_)
808 if (fixedPatchID == -1)
811 <<
"Wrong faPatch name fixedFreeSurfacePatches list"
812 <<
" defined in the dynamicMeshDict dictionary"
821 deltaHf[eFaces[edgeI]] *= 2.0;
825 controlPoints() += facesDisplacementDir()*deltaHf;
827 displacement = pointDisplacement();
830 refCast<velocityMotionSolver>
838 fixedValuePointPatchVectorField& fsPatchPointMeshU =
839 refCast<fixedValuePointPatchVectorField>
854void Foam::interfaceTrackingFvMesh::updateSurfaceFlux()
860void Foam::interfaceTrackingFvMesh::updateSurfactantConcentration()
862 if (!pureFreeSurface())
872 +
fam::div(Phis(), surfactantConcentration())
875 surfactant().diffusion(),
876 surfactantConcentration()
880 if (surfactant().soluble())
896 zeroGradientFaPatchScalarField::typeName
900 bulkSurfactantConcentration().boundaryField()[fsPatchIndex()];
901 Cb.correctBoundaryConditions();
906 surfactant().adsorptionCoeff()*
Cb
907 + surfactant().adsorptionCoeff()
908 *surfactant().desorptionCoeff(),
909 surfactantConcentration()
911 - surfactant().adsorptionCoeff()
912 *
Cb*surfactant().saturatedConc();
920 sigma() + surfactant().dSigma(surfactantConcentration())/rho_;
922 if (
neg(
min(surfaceTension().internalField().
field())))
925 <<
"Surface tension is negative"
932Foam::vector Foam::interfaceTrackingFvMesh::totalPressureForce()
const
942 return gSum(pressureForces);
946Foam::vector Foam::interfaceTrackingFvMesh::totalViscousForce()
const
966 U().boundaryField()[fsPatchIndex()].
snGrad()
979 return gSum(viscousForces);
983Foam::vector Foam::interfaceTrackingFvMesh::totalSurfaceTensionForce()
const
993 if (pureFreeSurface())
1000 )().internalField();
1004 surfTensionForces = surfaceTension().internalField().field()*
K*S*
n;
1007 return gSum(surfTensionForces);
1011Foam::scalar Foam::interfaceTrackingFvMesh::maxCourantNumber()
1015 if (pureFreeSurface())
1021 mesh().time().deltaT().value()/
1040 mesh().time().deltaT().value()/
1052void Foam::interfaceTrackingFvMesh::updateProperties()
1057 "transportProperties"
1066void Foam::interfaceTrackingFvMesh::correctPointDisplacement
1081 for (
const word& patchName : fixedFreeSurfacePatches_)
1095 label curFace = eFaces[edgeI];
1097 const labelList& curPoints = faces[curFace];
1099 forAll(curPoints, pointI)
1101 label curPoint = curPoints[pointI];
1102 label index = pLabels.
find(curPoint);
1115 forAll(corrPoints, pointI)
1117 label curPoint = corrPoints[pointI];
1123 label curFace =
pFaces[curPoint][faceI];
1125 label index = eFaces.
find(curFace);
1136 forAll(corrPoints, pointI)
1138 label curPoint = corrPoints[pointI];
1142 const labelList& curPointFaces = corrPointFaces[pointI];
1144 forAll(curPointFaces, faceI)
1146 const face& curFace = faces[curPointFaces[faceI]];
1148 label ptInFace = curFace.
which(curPoint);
1149 label next = curFace.
nextLabel(ptInFace);
1150 label prev = curFace.
prevLabel(ptInFace);
1154 const vector&
c = pointsDisplacementDir()[curPoint];
1156 curDisp += 2*sweptVolCorr[curPointFaces[faceI]]/((a^
b)&c);
1159 curDisp /= curPointFaces.
size();
1161 displacement[curPoint] =
1162 curDisp*pointsDisplacementDir()[curPoint];
1167 for (
const word& patchName : nonReflectingFreeSurfacePatches_)
1169 label nonReflectingPatchID =
1182 forAll(corrPoints, pointI)
1184 label curPoint = corrPoints[pointI];
1190 label curFace =
pFaces[curPoint][faceI];
1192 label index = eFaces.
find(curFace);
1204 forAll(corrPoints, pointI)
1206 label curPoint = corrPoints[pointI];
1210 const labelList& curPointFaces = corrPointFaces[pointI];
1212 forAll(curPointFaces, faceI)
1214 const face& curFace = faces[curPointFaces[faceI]];
1216 label ptInFace = curFace.
which(curPoint);
1217 label next = curFace.
nextLabel(ptInFace);
1218 label prev = curFace.
prevLabel(ptInFace);
1224 if (corrPoints.
find(next) == -1)
1239 vector c0 = displacement[p1];
1241 scalar V0 =
mag((a0^b0)&c0)/2;
1243 scalar DV = sweptVolCorr[curPointFaces[faceI]] - V0;
1245 if (corrPoints.
find(prev) != -1)
1249 (
points[next] + displacement[next])
1251 const vector&
c = pointsDisplacementDir()[curPoint];
1253 curDisp += 2*DV/((a^
b)&c);
1258 - (
points[prev] + displacement[prev]);
1260 const vector&
c = pointsDisplacementDir()[curPoint];
1262 curDisp += 2*DV/((a^
b)&c);
1266 curDisp /= curPointFaces.size();
1268 displacement[curPoint] =
1269 curDisp*pointsDisplacementDir()[curPoint];
1275void Foam::interfaceTrackingFvMesh::correctContactLinePointNormals()
1286 if (contactAnglePtr_ && correctContactLineNormals())
1288 Info<<
"Correcting contact line normals" <<
endl;
1331 refCast<const processorFaPatch>
1339 forAll(patchPointLabels, pointI)
1341 label curPoint = patchPointLabels[pointI];
1348 const labelList& curPointEdges = pointEdges[curPoint];
1350 forAll(curPointEdges, edgeI)
1352 label curEdge = curPointEdges[edgeI];
1354 if (edgeFaces[curEdge].size() == 1)
1361 label index = curEdges.
find(curEdge);
1385 vector t = edges[curEdge].vec(oldPoints);
1386 t /=
mag(t) + SMALL;
1389 pointFaces[curPoint];
1391 forAll(curPointFaces, fI)
1393 tangent.ref().field()[curPointFaces[fI]] = t;
1400 tangent.correctBoundaryConditions();
1405 label ngbPolyPatchID =
1408 if (ngbPolyPatchID != -1)
1421 - contactAnglePtr_->boundaryField()[patchI]
1436 rotationAxis /=
mag(rotationAxis) + SMALL;
1446 forAll(pointEdges, pointI)
1450 forAll(pointEdges[pointI], eI)
1454 + pointEdges[pointI][eI];
1456 vector e = edges[curEdge].vec(oldPoints);
1458 e *= (
e&rotationAxis[pointI])
1459 /
mag(
e&rotationAxis[pointI]);
1461 e /=
mag(
e) + SMALL;
1466 if (pointEdges[pointI].size() == 1)
1468 label curPoint = patchPoints[pointI];
1475 label procPatchID = -1;
1481 forAll(curPointEdges, edgeI)
1483 label curEdge = curPointEdges[edgeI];
1485 if (edgeFaces[curEdge].size() == 1)
1493 curEdges.
find(curEdge);
1512 if (procPatchID != -1)
1515 tangent.boundaryField()[procPatchID]
1516 .patchNeighbourField()()[
edgeID];
1518 t *= (t&rotationAxis[pointI])
1519 /(
mag(t&rotationAxis[pointI]) + SMALL);
1521 t /=
mag(t) + SMALL;
1527 rotationAxis[pointI] = rotAx/(
mag(rotAx) + SMALL);
1531 ngbN = ngbN*
cos(rotAngle)
1532 + rotationAxis*(rotationAxis & ngbN)*(1 -
cos(rotAngle))
1533 + (rotationAxis^ngbN)*
sin(rotAngle);
1536 forAll(patchPoints, pointI)
1538 N[patchPoints[pointI]] -=
1539 ngbN[pointI]*(ngbN[pointI]&
N[patchPoints[pointI]]);
1541 N[patchPoints[pointI]] /=
1542 mag(
N[patchPoints[pointI]]) + SMALL;
1564 fixedFreeSurfacePatches_(),
1565 nonReflectingFreeSurfacePatches_(),
1566 pointNormalsCorrectionPatches_(),
1567 normalMotionDir_(false),
1570 pureFreeSurface_(true),
1571 rigidFreeSurface_(false),
1572 correctContactLineNormals_(false),
1577 controlPointsPtr_(nullptr),
1578 motionPointsMaskPtr_(nullptr),
1579 pointsDisplacementDirPtr_(nullptr),
1580 facesDisplacementDirPtr_(nullptr),
1581 fsNetPhiPtr_(nullptr),
1583 surfactConcPtr_(nullptr),
1584 bulkSurfactConcPtr_(nullptr),
1585 surfaceTensionPtr_(nullptr),
1586 surfactantPtr_(nullptr),
1587 contactAnglePtr_(nullptr)
1670 aMeshPtr_.reset(
new faMesh(*
this));
1673 fixedFreeSurfacePatches_ =
1674 motion().get<
wordList>(
"fixedFreeSurfacePatches");
1676 pointNormalsCorrectionPatches_ =
1677 motion().
get<
wordList>(
"pointNormalsCorrectionPatches");
1679 normalMotionDir_ = motion().
get<
bool>(
"normalMotionDir");
1680 smoothing_ = motion().getOrDefault(
"smoothing",
false);
1681 pureFreeSurface_ = motion().getOrDefault(
"pureFreeSurface",
true);
1718 return *fsNetPhiPtr_;
1729 return *fsNetPhiPtr_;
1735 forAll(
Us().boundaryField(), patchI)
1739 Us().boundaryField()[patchI].
type()
1740 == calculatedFaPatchVectorField::typeName
1745 pUs =
Us().boundaryField()[patchI].patchInternalField();
1747 label ngbPolyPatchID =
1750 if (ngbPolyPatchID != -1)
1755 U().boundaryField()[ngbPolyPatchID].
type()
1756 == slipFvPatchVectorField::typeName
1760 U().boundaryField()[ngbPolyPatchID].
type()
1761 == symmetryFvPatchVectorField::typeName
1776 Us().correctBoundaryConditions();
1784 Us().ref().field() =
U().boundaryField()[fsPatchIndex()];
1792 correctUsBoundaryConditions();
1798 return *getObjectPtr<const volVectorField>(
"U");
1804 return *getObjectPtr<const volScalarField>(
"p");
1810 return *getObjectPtr<const surfaceScalarField>(
"phi");
1818 auto& SnGradU = tSnGradU.ref();
1825 -
aMesh().faceCurvatures()*(
aMesh().faceAreaNormals()&
Us())
1832 gradUs -=
n*(
n & gradUs);
1842 if (!pureFreeSurface() &&
max(
nu) > SMALL)
1844 tangentialSurfaceTensionForce =
1845 surfaceTensionGrad()().internalField();
1849 tangentialSurfaceTensionForce/(
nu + SMALL)
1861 auto& SnGradUn = tSnGradUn.ref();
1866 -
aMesh().faceCurvatures()*(
aMesh().faceAreaNormals()&
Us())
1879 auto& pressureJump = tPressureJump.ref();
1893 + 2.0*
nu*freeSurfaceSnGradUn();
1895 if (pureFreeSurface())
1897 pressureJump -= sigma().
value()*
K;
1901 pressureJump -= surfaceTension().internalField()*
K;
1904 return tPressureJump;
1910 if (!controlPointsPtr_)
1912 makeControlPoints();
1915 return *controlPointsPtr_;
1921 if (!motionPointsMaskPtr_)
1923 makeMotionPointsMask();
1926 return *motionPointsMaskPtr_;
1932 if (!pointsDisplacementDirPtr_)
1937 return *pointsDisplacementDirPtr_;
1943 if (!facesDisplacementDirPtr_)
1948 return *facesDisplacementDirPtr_;
1965 if (!surfactConcPtr_)
1970 return *surfactConcPtr_;
1977 if (!surfactConcPtr_)
1982 return *surfactConcPtr_;
1989 if (!bulkSurfactConcPtr_)
1991 makeBulkSurfactConc();
1994 return *bulkSurfactConcPtr_;
2001 if (!bulkSurfactConcPtr_)
2003 makeBulkSurfactConc();
2006 return *bulkSurfactConcPtr_;
2013 if (!surfaceTensionPtr_)
2015 makeSurfaceTension();
2018 return *surfaceTensionPtr_;
2025 if (!surfaceTensionPtr_)
2027 makeSurfaceTension();
2030 return *surfaceTensionPtr_;
2037 if (!surfactantPtr_)
2042 return *surfactantPtr_;
2053 "surfaceTensionGrad",
2065 auto& grad = tgrad.ref();
2069 grad -=
n*(
n & grad);
2070 grad.correctBoundaryConditions();
2080 if (smoothing_ && !rigidFreeSurface_)
2082 smoothFreeSurfaceMesh();
2083 clearControlPoints();
2086 updateDisplacementDirections();
2090 Info<<
"Maximal capillary Courant number: "
2091 << maxCourantNumber() <<
endl;
2095 Info<<
"Free surface curvature: min = " <<
gMin(
K)
2101 if (!rigidFreeSurface_)
2117 Info<<
"Free surface continuity error : sum local = "
2118 <<
gSum(
mag(sweptVolCorr)) <<
", global = " <<
gSum(sweptVolCorr)
2122 fsNetPhi().ref().field() = sweptVolCorr;
2128 "ddt(" +
U().
name() +
')'
2166 <<
"Unsupported temporal differencing scheme : "
2176 sweptVolCorr/(Sf*(Nf & facesDisplacementDir()))
2179 for (
const word& patchName : fixedFreeSurfacePatches_)
2181 label fixedPatchID =
2184 if (fixedPatchID == -1)
2187 <<
"Wrong faPatch name '" << patchName
2188 <<
"'in the fixedFreeSurfacePatches list"
2189 <<
" defined in dynamicMeshDict dictionary"
2198 deltaHf[eFaces[edgeI]] *= 2.0;
2202 controlPoints() += facesDisplacementDir()*deltaHf;
2204 pointField displacement(pointDisplacement());
2205 correctPointDisplacement(sweptVolCorr, displacement);
2208 refCast<velocityMotionSolver>
2216 fixedValuePointPatchVectorField& fsPatchPointMeshU =
2217 refCast<fixedValuePointPatchVectorField>
2225 fsPatchPointMeshU ==
2230 correctContactLinePointNormals();
2241 refCast<velocityMotionSolver>
2249 fixedValuePointPatchVectorField& fsPatchPointMeshU =
2250 refCast<fixedValuePointPatchVectorField>
2258 fsPatchPointMeshU ==
2266 updateSurfactantConcentration();
2278 mesh().time().timePath()/
"freeSurface",
2290 mesh().time().timePath()/
"freeSurfaceControlPoints.vtk"
2293 Info<<
"Writing free surface control points to " <<
os.
name() <<
nl;
2295 os <<
"# vtk DataFile Version 2.0" <<
nl
2296 <<
"freeSurfaceControlPoints" <<
nl
2298 <<
"DATASET POLYDATA" <<
nl;
2300 const label
nPoints = controlPoints().size();
2303 for (
const point&
p : controlPoints())
2305 os << float(
p.x()) <<
' '
2306 << float(
p.y()) <<
' '
2307 << float(
p.z()) <<
nl;
2311 for (label
id = 0;
id <
nPoints; ++id)
2313 os << 1 <<
' ' <<
id <<
nl;
CGAL::Exact_predicates_exact_constructions_kernel K
CGAL::indexedCell< K > Cb
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
vtk::internalMeshWriter writer(topoMesh, topoCells, vtk::formatType::INLINE_ASCII, runTime.path()/"blockTopology")
const uniformDimensionedVectorField & g
An indexed form of CGAL::Triangulation_cell_base_3<K> used to keep track of the Delaunay cells (tets)...
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
void correctBoundaryConditions()
Correct boundary field.
const Internal & internalField() const
Return a const-reference to the dimensioned internal field.
bool insert(const Key &key)
Insert a new entry, not overwriting existing entries.
List< Key > toc() const
The table of contents (the keys) in unsorted order.
A primitive field of type <T> with automated input and output.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Output to file stream, using an OSstream.
virtual const fileName & name() const
Read/write access to the name of the stream.
const labelListList & pointEdges() const
Return point-edge addressing.
const labelList & meshPoints() const
Return labelList of mesh points in patch.
const Field< point_type > & localPoints() const
Return pointField of points in patch.
const labelListList & pointFaces() const
Return point-face addressing.
const labelListList & edgeFaces() const
Return edge-face addressing.
const List< face_type > & localFaces() const
Return patch faces addressing into local point list.
label timeIndex() const noexcept
Return current time index.
dimensionedScalar deltaT() const
Return time step.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
std::enable_if< std::is_same< bool, TypeT >::value, bool >::type get(const label i) const
label find(const T &val, label pos=0) const
Find index of the first occurrence of the value.
void size(const label n)
Older name for setAddressableSize.
static bool & parRun() noexcept
Test if this a parallel run.
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface....
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
Generic dimensioned Type class.
const Type & value() const
Return const reference to value.
Abstract base class for geometry and/or topology changing fvMesh.
The dynamicMotionSolverFvMesh.
virtual bool update()
Update the mesh for both mesh motion and topology change.
const motionSolver & motion() const
Return the motionSolver.
const edgeScalarField & lPN() const
Return reference to PN geodesic distance.
UPtrList< const labelUList > edgeFaces() const
Return a list of edgeFaces for each patch.
label findPatchID(const word &patchName, const bool allowNotFound=true) const
Find patch index given a name, return -1 if not found.
label start() const
The start label of the edges in the faMesh edges list.
Finite area mesh (used for 2-D non-Euclidian finite area method) defined using a patch of faces on a ...
const faBoundaryMesh & boundary() const noexcept
Return constant reference to boundary mesh.
const faceList & faces() const
Return local faces.
const DimensionedField< scalar, areaMesh > & S() const
Return face areas.
bool correctPatchPointNormals(const label patchID) const
Whether point normals should be corrected for a patch.
const edgeList & edges() const noexcept
Return local edges with reordered boundary.
const uindirectPrimitivePatch & patch() const
Return constant reference to primitive patch.
const edgeVectorField & Le() const
Return edge length vectors.
const areaScalarField & faceCurvatures() const
Return face curvatures.
const pointField & points() const
Return local points.
const vectorField & pointAreaNormals() const
Return point area normals.
const areaVectorField & faceAreaNormals() const
Return face area normals.
const areaVectorField & areaCentres() const
Return face centres as areaVectorField.
const labelList & pointLabels() const
Return patch point labels.
Template specialisation for scalar faMatrix.
A face is a list of labels corresponding to mesh vertices.
label which(const label pointLabel) const
Find local index on face for the point label, same as find()
label nextLabel(const label i) const
Next vertex on face.
label prevLabel(const label i) const
Previous vertex on face.
const Time & time() const
Return the top-level database.
const surfaceVectorField & Cf() const
Return face centres as surfaceVectorField.
Second-oder Crank-Nicolson implicit ddt using the current and previous time-step fields as well as th...
Basic first-order Euler implicit/explicit ddt using only the current and previous time-step values.
Second-order backward-differencing ddt using the current and two previous time-step values.
The interfaceTrackingFvMesh.
vectorField & controlPoints()
Return control points.
void correctUsBoundaryConditions()
Correct surface velocity boundary conditions.
areaVectorField & Us()
Return free-surface velocity field.
const areaScalarField & fsNetPhi() const
Return free-surface net flux.
const surfaceScalarField & phi() const
Return constant reference to velocity field.
tmp< scalarField > freeSurfaceSnGradUn()
Return free surface normal derivative of normal velocity comp.
vectorField & pointsDisplacementDir()
Return reference to point displacement direction field.
virtual bool init(const bool doInit)
Initialise all non-demand-driven data.
const volVectorField & U() const
Return constant reference to velocity field.
void writeVTK() const
Write VTK freeSurface mesh.
void updateUs()
Update free-surface velocity field.
areaScalarField & surfaceTension()
Return surface tension field.
vectorField & facesDisplacementDir()
Return reference to control points displacement direction field.
void writeVTKControlPoints()
Write VTK freeSurface control points.
faMesh & aMesh()
Return reference to finite area mesh.
virtual bool update()
Update the mesh for both mesh motion and topology change.
* ~interfaceTrackingFvMesh()
Destructor.
tmp< areaVectorField > surfaceTensionGrad()
Return surface tension gradient.
labelList & motionPointsMask()
Return reference to motion points mask field.
volScalarField & bulkSurfactantConcentration()
Return volume surfactant concentration field.
const volScalarField & p() const
Return constant reference to pressure field.
const surfactantProperties & surfactant() const
Return surfactant properties.
areaScalarField & surfactantConcentration()
Return free-surface surfactant concentration field.
edgeScalarField & Phis()
Return free-surface fluid flux field.
tmp< scalarField > freeSurfacePressureJump()
Return free surface pressure jump.
tmp< vectorField > freeSurfaceSnGradU()
Return free surface normal derivative of velocity.
Virtual base class for mesh motion solver.
const Type & lookupObject(const word &name, const bool recursive=false) const
virtual const pointField & oldPoints() const
Return old points (mesh motion)
A simple single-phase transport model based on viscosityModel.
A class for managing temporary objects.
Abstract base class for turbulence models (RAS, LES and laminar).
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
Virtual base class for velocity motion solver.
pointVectorField & pointMotionU()
Return reference to the point motion velocity field.
Write concrete PrimitivePatch faces/points (optionally with fields) as a vtp file or a legacy vtk fil...
virtual bool writeGeometry()
Write mesh topology.
A class for handling words, derived from Foam::string.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
const volScalarField & p0
Template functions to aid in the implementation of demand driven data.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
OBJstream os(runTime.globalPath()/outputName)
const labelList nFaces(UPstream::listGatherValues< label >(aMesh.nFaces()))
compressible::turbulenceModel & turbulence
Calculate the mesh motion flux and convert fluxes from absolute to relative and back.
IOobject io("surfaceFilmProperties", mesh.time().constant(), mesh, IOobject::READ_IF_PRESENT, IOobject::NO_WRITE, false)
#define DebugInFunction
Report an information message using Foam::Info.
const dimensionedScalar a0
Bohr radius: default SI units: [m].
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
const dimensionedScalar c
Speed of light in a vacuum.
static tmp< GeometricField< Type, faePatchField, edgeMesh > > interpolate(const GeometricField< Type, faPatchField, areaMesh > &tvf, const edgeScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< GeometricField< Type, faPatchField, areaMesh > > div(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
tmp< GeometricField< Type, faPatchField, areaMesh > > edgeIntegrate(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
tmp< GeometricField< typename outerProduct< vector, Type >::type, faPatchField, areaMesh > > grad(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
tmp< faMatrix< Type > > div(const edgeScalarField &flux, const GeometricField< Type, faPatchField, areaMesh > &vf, const word &name)
tmp< faMatrix< Type > > laplacian(const GeometricField< Type, faPatchField, areaMesh > &vf)
tmp< faMatrix< Type > > ddt(const GeometricField< Type, faPatchField, areaMesh > &vf)
zeroField Sp(const Foam::zero, const GeometricField< Type, faPatchField, areaMesh > &)
A no-op source.
const std::string patch
OpenFOAM patch number as a std::string.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
@ LEGACY_ASCII
Legacy ASCII, legacyAsciiFormatter.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Type gSum(const FieldField< Field, Type > &f)
const dimensionSet dimless
Dimensionless.
List< label > labelList
A List of labels.
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
const dimensionSet dimMoles(0, 0, 0, 0, 1, 0, 0)
const dimensionSet dimArea(sqr(dimLength))
dimensionedScalar sin(const dimensionedScalar &ds)
const dimensionSet dimVelocity
tmp< GeometricField< Type, faePatchField, edgeMesh > > linearEdgeInterpolate(const GeometricField< Type, faPatchField, areaMesh > &vf)
messageStream Info
Information stream (stdout output on master, null elsewhere)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
IOField< vector > vectorIOField
vectorField with IO.
Ostream & endl(Ostream &os)
Add newline and flush stream.
constexpr scalar degToRad() noexcept
Multiplication factor for degrees to radians conversion.
GeometricField< vector, faPatchField, areaMesh > areaVectorField
GeometricField< scalar, faePatchField, edgeMesh > edgeScalarField
dimensionedScalar sqrt(const dimensionedScalar &ds)
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
errorManip< error > abort(error &err)
const dimensionSet dimForce
Field< vector > vectorField
Specialisation of Field<T> for vector.
dimensionedScalar neg(const dimensionedScalar &ds)
Type gAverage(const FieldField< Field, Type > &f)
static constexpr const zero Zero
Global zero (0)
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
Type gMin(const FieldField< Field, Type > &f)
const dimensionSet dimVolume(pow3(dimLength))
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
const dimensionSet dimDensity
void deleteDemandDrivenData(DataPtr &dataPtr)
Type gMax(const FieldField< Field, Type > &f)
Ostream & flush(Ostream &os)
Flush stream.
dimensionedScalar cos(const dimensionedScalar &ds)
constexpr char nl
The newline '\n' character (0x0a)
Info<< "Finished reading KIVA file"<< endl;cellShapeList cellShapes(nPoints);labelList cellZoning(nPoints, -1);const cellModel &hex=cellModel::ref(cellModel::HEX);labelList hexLabels(8);label activeCells=0;labelList pointMap(nPoints);forAll(pointMap, i){ pointMap[i]=i;}for(label i=0;i< nPoints;i++){ if(f[i] > 0.0) { hexLabels[0]=i;hexLabels[1]=i1tab[i];hexLabels[2]=i3tab[i1tab[i]];hexLabels[3]=i3tab[i];hexLabels[4]=i8tab[i];hexLabels[5]=i1tab[i8tab[i]];hexLabels[6]=i3tab[i1tab[i8tab[i]]];hexLabels[7]=i3tab[i8tab[i]];cellShapes[activeCells].reset(hex, hexLabels);edgeList edges=cellShapes[activeCells].edges();forAll(edges, ei) { if(edges[ei].mag(points)< SMALL) { label start=pointMap[edges[ei].start()];while(start !=pointMap[start]) { start=pointMap[start];} label end=pointMap[edges[ei].end()];while(end !=pointMap[end]) { end=pointMap[end];} label minLabel=min(start, end);pointMap[start]=pointMap[end]=minLabel;} } cellZoning[activeCells]=idreg[i];activeCells++;}}cellShapes.setSize(activeCells);cellZoning.setSize(activeCells);forAll(cellShapes, celli){ cellShape &cs=cellShapes[celli];forAll(cs, i) { cs[i]=pointMap[cs[i]];} cs.collapse();}label bcIDs[11]={-1, 0, 2, 4, -1, 5, -1, 6, 7, 8, 9};const label nBCs=12;const word *kivaPatchTypes[nBCs]={ &wallPolyPatch::typeName, &wallPolyPatch::typeName, &wallPolyPatch::typeName, &wallPolyPatch::typeName, &symmetryPolyPatch::typeName, &wedgePolyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &symmetryPolyPatch::typeName, &oldCyclicPolyPatch::typeName};enum patchTypeNames{ PISTON, VALVE, LINER, CYLINDERHEAD, AXIS, WEDGE, INFLOW, OUTFLOW, PRESIN, PRESOUT, SYMMETRYPLANE, CYCLIC};const char *kivaPatchNames[nBCs]={ "piston", "valve", "liner", "cylinderHead", "axis", "wedge", "inflow", "outflow", "presin", "presout", "symmetryPlane", "cyclic"};List< SLList< face > > pFaces[nBCs]
#define forAll(list, i)
Loop across all elements in list.
static const char *const typeName
The type name used in ensight case files.
Unit conversion functions.
const Vector< label > N(dict.get< Vector< label > >("N"))