70 static const scalar edgeTol = 1
e-3;
75 Info<<
"Writing " << msg <<
" (" <<
cells.
size() <<
") to cellSet " 90 if (
mag(normal &
vector(1, 0, 0)) > 1-edgeTol)
94 else if (
mag(normal &
vector(0, 1, 0)) > 1-edgeTol)
98 else if (
mag(normal &
vector(0, 0, 1)) > 1-edgeTol)
117 scalar maxX = -GREAT;
121 scalar maxY = -GREAT;
125 scalar maxZ = -GREAT;
128 scalar minOther = GREAT;
129 scalar maxOther = -GREAT;
135 const edge&
e = edges[edgei];
139 scalar eMag =
mag(eVec);
143 if (
mag(eVec &
x) > 1-edgeTol)
145 minX =
min(minX, eMag);
146 maxX =
max(maxX, eMag);
149 else if (
mag(eVec &
y) > 1-edgeTol)
151 minY =
min(minY, eMag);
152 maxY =
max(maxY, eMag);
155 else if (
mag(eVec & z) > 1-edgeTol)
157 minZ =
min(minZ, eMag);
158 maxZ =
max(maxZ, eMag);
163 minOther =
min(minOther, eMag);
164 maxOther =
max(maxOther, eMag);
169 <<
"Mesh edge statistics:" <<
nl 170 <<
" x aligned : number:" << nX <<
"\tminLen:" << minX
171 <<
"\tmaxLen:" << maxX <<
nl 172 <<
" y aligned : number:" << nY <<
"\tminLen:" << minY
173 <<
"\tmaxLen:" << maxY <<
nl 174 <<
" z aligned : number:" << nZ <<
"\tminLen:" << minZ
175 <<
"\tmaxLen:" << maxZ <<
nl 176 <<
" other : number:" <<
mesh.
nEdges() - nX - nY - nZ
177 <<
"\tminLen:" << minOther
178 <<
"\tmaxLen:" << maxOther <<
nl <<
endl;
180 if (excludeCmpt == 0)
182 return min(minY,
min(minZ, minOther));
184 else if (excludeCmpt == 1)
186 return min(minX,
min(minZ, minOther));
188 else if (excludeCmpt == 2)
190 return min(minX,
min(minY, minOther));
194 return min(minX,
min(minY,
min(minZ, minOther)));
225 emptyPolyPatch::typeName
250 Info<<
"Created patch oldInternalFaces at " << patchi <<
endl;
254 Info<<
"Reusing patch oldInternalFaces at " << patchi <<
endl;
263 void selectCurvatureCells
268 const scalar nearDist,
269 const scalar curvature,
301 void addCutNeighbours
304 const bool selectInside,
305 const bool selectOutside,
315 for (
const label celli : cutCells)
321 const label facei = cFaces[i];
332 if (selectInside && inside.
found(nbr))
334 addCutFaces.insert(nbr);
336 else if (selectOutside && outside.
found(nbr))
338 addCutFaces.insert(nbr);
344 Info<<
" Selected an additional " << addCutFaces.size()
345 <<
" neighbours of cutCells to refine" <<
endl;
347 for (
const label facei : addCutFaces)
349 cutCells.insert(facei);
358 bool limitRefinementLevel
361 const label limitDiff,
370 if (!excludeCells.
found(celli))
376 label nbr = cCells[i];
378 if (!excludeCells.
found(nbr))
380 if (refLevel[celli] - refLevel[nbr] >= limitDiff)
383 <<
"Level difference between neighbouring cells " 384 << celli <<
" and " << nbr
385 <<
" greater than or equal to " << limitDiff <<
endl 386 <<
"refLevels:" << refLevel[celli] <<
' ' 397 for (
const label celli : cutCells)
404 const label nbr = cCells[i];
406 if (!excludeCells.
found(nbr) && !cutCells.found(nbr))
408 if (refLevel[celli] + 1 - refLevel[nbr] >= limitDiff)
410 addCutCells.insert(nbr);
416 if (addCutCells.size())
420 Info<<
"Added an additional " << addCutCells.size() <<
" cells" 421 <<
" to satisfy 1:" << limitDiff <<
" refinement level" 424 for (
const label celli : addCutCells)
426 cutCells.insert(celli);
431 Info<<
"Added no additional cells" 432 <<
" to satisfy 1:" << limitDiff <<
" refinement level" 465 for (label celli = oldCells; celli <
mesh.
nCells(); celli++)
472 forAll(addedCells, oldCelli)
474 const labelList& added = addedCells[oldCelli];
480 label masterLevel = ++refLevel[oldCelli];
484 refLevel[added[i]] = masterLevel;
495 const label writeMesh,
507 <<
" Removing:" << cellLabels.size() <<
" cells" <<
endl;
510 labelList exposedFaces(cellRemover.getExposedFaces(cellLabels));
513 cellRemover.setRefinement
528 if (morphMap().hasMotionPoints())
534 cellRemover.updateMesh(morphMap());
537 const labelList& cellMap = morphMap().cellMap();
543 newRefLevel[i] = refLevel[cellMap[i]];
577 const bool selectCut,
578 const bool selectInside,
579 const bool selectOutside,
581 const label nCutLayers,
608 selected.
addSet(cutCells);
619 Info<<
"Determined cell status:" <<
endl 620 <<
" inside :" << inside.size() <<
endl 621 <<
" outside :" << outside.size() <<
endl 622 <<
" cutCells:" << cutCells.
size() <<
endl 623 <<
" selected:" << selected.
size() <<
endl 626 writeSet(inside,
"inside cells");
627 writeSet(outside,
"outside cells");
628 writeSet(cutCells,
"cut cells");
629 writeSet(selected,
"selected cells");
634 int main(
int argc,
char *argv[])
638 "Refine cells near to a surface" 647 label newPatchi = addPatch(
mesh,
"oldInternalFaces");
654 Info<<
"Checking for motionProperties\n" <<
endl;
675 if (motionProperties.get<
bool>(
"twoDMotion"))
686 "snappyRefineMeshDict",
696 const label nCutLayers(refineDict.
get<label>(
"nCutLayers"));
697 const label cellLimit(refineDict.
get<label>(
"cellLimit"));
698 const bool selectCut(refineDict.
get<
bool>(
"selectCut"));
699 const bool selectInside(refineDict.
get<
bool>(
"selectInside"));
700 const bool selectOutside(refineDict.
get<
bool>(
"selectOutside"));
701 const bool selectHanging(refineDict.
get<
bool>(
"selectHanging"));
703 const scalar minEdgeLen(refineDict.
get<scalar>(
"minEdgeLen"));
704 const scalar maxEdgeLen(refineDict.
get<scalar>(
"maxEdgeLen"));
705 const scalar curvature(refineDict.
get<scalar>(
"curvature"));
706 const scalar curvDist(refineDict.
get<scalar>(
"curvatureDistance"));
708 const label refinementLimit(refineDict.
get<label>(
"splitLevel"));
709 const bool writeMesh(refineDict.
get<
bool>(
"writeMesh"));
711 Info<<
"Cells to be used for meshing (0=false, 1=true):" <<
nl 712 <<
" cells cut by surface : " << selectCut <<
nl 713 <<
" cells inside of surface : " << selectInside <<
nl 714 <<
" cells outside of surface : " << selectOutside <<
nl 715 <<
" hanging cells : " << selectHanging <<
nl 719 if (nCutLayers > 0 && selectInside)
722 <<
"Illogical settings : Both nCutLayers>0 and selectInside true." 724 <<
"This would mean that inside cells get removed but should be" 725 <<
" included in final mesh" <<
endl;
738 forAll(outsidePts, outsidei)
740 const point& outsidePoint = outsidePts[outsidei];
742 if (queryMesh.findCell(outsidePoint, -1,
false) == -1)
745 <<
"outsidePoint " << outsidePoint
746 <<
" is not inside any cell" 768 label maxLevel =
max(refLevel);
772 Info<<
"Read existing refinement level from file " 774 <<
" min level : " <<
min(refLevel) <<
nl 775 <<
" max level : " << maxLevel <<
nl 780 Info<<
"Created zero refinement level in file " 789 direction normalDir(getNormalDir(correct2DPtr));
790 scalar meshMinEdgeLen = getEdgeStats(
mesh, normalDir);
792 while (meshMinEdgeLen > minEdgeLen)
819 Info<<
" Selected " << cutCells.
name() <<
" with " 820 << cutCells.
size() <<
" cells" <<
endl;
822 if ((curvDist > 0) && (meshMinEdgeLen < maxEdgeLen))
839 Info<<
" Selected " << curveCells.name() <<
" with " 840 << curveCells.size() <<
" cells" <<
endl;
859 cutCells.
subset(curveCells);
861 Info<<
" Removed cells not on surface curvature. Selected " 870 subsetMesh(
mesh, writeMesh, newPatchi, inside, cutCells, refLevel);
890 <<
" Selected for refinement :" << cutCells.
size() <<
nl 893 if (cutCells.
empty())
895 Info<<
"Stopping refining since 0 cells selected to be refined ..." 902 Info<<
"Stopping refining since cell limit reached ..." <<
nl 931 meshMinEdgeLen = getEdgeStats(
mesh, normalDir);
967 Info<<
"Detected " << hanging.
size() <<
" hanging cells" 968 <<
" (cells with all points on" 969 <<
" outside of cellSet 'selected').\nThese would probably result" 970 <<
" in flattened cells when snapping the mesh to the surface" 973 Info<<
"Refining " << hanging.
size() <<
" hanging cells" <<
nl 990 doRefinement(
mesh, refineDict, hanging, refLevel);
1001 else if (!writeMesh)
label findPatchID(const word &patchName, const bool allowNotFound=true) const
Find patch index given a name, return -1 if not found.
ITstream & lookup(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return an entry data stream. FatalIOError if not found, or not a stream. ...
Various (local, not parallel) searches on polyMesh; uses (demand driven) octree to search...
static void addNote(const string ¬e)
Add extra notes for the usage information.
void size(const label n)
Older name for setAddressableSize.
A class for handling file names.
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
void transfer(List< T > &list)
Transfer the contents of the argument List into this list and annul the argument list.
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
bool found(const Key &key) const
Same as contains()
virtual const labelList & faceNeighbour() const
Return face neighbour.
const word & name() const noexcept
Return the object name.
Cell-face mesh analysis engine.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Add elements to current set.
constexpr char nl
The newline '\n' character (0x0a)
Does multiple pass refinement to refine cells in multiple directions.
Given list of cells to remove, insert all the topology changes.
Ostream & endl(Ostream &os)
Add newline and flush stream.
virtual void updateMesh(const mapPolyMesh &morphMap)
Update any stored data for new labels.
virtual void movePoints(const pointField &)
Move points, returns volumes swept by faces in motion.
static unsigned int defaultPrecision() noexcept
Return the default precision.
static void noParallel()
Remove the parallel options.
A bounding box defined in terms of min/max extrema points.
const cellList & cells() const
Ignore writing from objectRegistry::writeObject()
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a T. FatalIOError if not found, or if the number of tokens is incorrect.
const Time & time() const
Return the top-level database.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
fileName objectPath() const
The complete path + object name.
bool isInternalFace(const label faceIndex) const noexcept
Return true if given face label is internal to the mesh.
static labelHashSet getHangingCells(const primitiveMesh &mesh, const labelHashSet &internalCells)
Get cells using points on 'outside' only.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
virtual const pointField & points() const
Return raw points.
#define forAll(list, i)
Loop across all elements in list.
HashSet< label, Hash< label > > labelHashSet
A HashSet of labels, uses label hasher.
label size() const noexcept
The number of elements in table.
Helper class to search on triSurface.
vectorField pointField
pointField is a vectorField.
const dimensionedScalar e
Elementary charge.
void setSize(const label n)
Alias for resize()
virtual void addSet(const topoSet &set)
Add elements present in set.
An edge is a list of two vertex labels. This can correspond to a directed graph edge or an edge on a ...
Class applies a two-dimensional correction to mesh motion point field.
const polyBoundaryMesh & boundaryMesh() const noexcept
Return boundary mesh.
A class for handling words, derived from Foam::string.
static word defaultRegion
Return the default region name.
label size() const noexcept
The number of entries in the list.
const word & system() const noexcept
Return system name.
virtual const labelList & faceOwner() const
Return face owner.
Reading is optional [identical to LAZY_READ].
virtual bool write(const bool writeOnProc=true) const
Write using setting from DB.
virtual void subset(const topoSet &set)
Subset contents. Only elements present in both sets remain.
label nInternalFaces() const noexcept
Number of internal faces.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
const vectorField & cellCentres() const
errorManip< error > abort(error &err)
static void getSurfaceSets(const polyMesh &mesh, const fileName &surfName, const triSurface &surf, const triSurfaceSearch &querySurf, const pointField &outsidePts, const label nCutLayers, labelHashSet &inside, labelHashSet &outside, labelHashSet &cut)
Divide cells into cut,inside and outside.
label nEdges() const
Number of mesh edges.
virtual bool write(const bool writeOnProc=true) const
Write mesh using IO settings from time.
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO...
void removeBoundary()
Remove boundary patches.
bool empty() const noexcept
True if the hash table is empty.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
static word timeName(const scalar t, const int precision=precision_)
Return a time name for the given scalar time value formatted with the given precision.
const word & constant() const noexcept
Return constant name.
const vector & planeNormal() const
Return plane normal.
Empty front and back plane patch. Used for 2-D geometries.
const triSurface & surface() const
Return reference to the surface.
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
A topoSetCellSource to select cells based on relation to a surface given by an external file...
#define WarningInFunction
Report a warning using Foam::Warning.
string & expand(const bool allowEmpty=false)
Inplace expand initial tags, tildes, and all occurrences of environment variables as per stringOps::e...
label nCells() const noexcept
Number of mesh cells.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers...
A collection of cell labels.
Direct mesh changes based on v1.3 polyTopoChange syntax.
const polyBoundaryMesh & patches
Automatically write from objectRegistry::writeObject()
messageStream Info
Information stream (stdout output on master, null elsewhere)
Pointer management similar to std::unique_ptr, with some additional methods and type checking...
Mesh consisting of general polyhedral cells.
List< Key > toc() const
The table of contents (the keys) in unsorted order.
List< label > labelList
A List of labels.
A patch is a list of labels that address the faces in the global face list.
Triangulated surface description with patch information.
void addPatches(polyPatchList &plist, const bool validBoundary=true)
Add boundary patches.
Defines the attributes of an object for which implicit objectRegistry management is supported...
const labelListList & cellCells() const
uindirectPrimitivePatch pp(UIndirectList< face >(mesh.faces(), faceLabels), mesh.points())
static constexpr const zero Zero
Global zero (0)