64 static const scalar edgeTol = 1
e-3;
69 Info<<
"Writing " << msg <<
" (" <<
cells.
size() <<
") to cellSet "
84 if (
mag(normal &
vector(1, 0, 0)) > 1-edgeTol)
88 else if (
mag(normal &
vector(0, 1, 0)) > 1-edgeTol)
92 else if (
mag(normal &
vector(0, 0, 1)) > 1-edgeTol)
111 scalar maxX = -great;
115 scalar maxY = -great;
119 scalar maxZ = -great;
122 scalar minOther = great;
123 scalar maxOther = -great;
129 const edge&
e = edges[edgeI];
133 scalar eMag =
mag(eVec);
137 if (
mag(eVec &
x) > 1-edgeTol)
139 minX =
min(minX, eMag);
140 maxX =
max(maxX, eMag);
143 else if (
mag(eVec &
y) > 1-edgeTol)
145 minY =
min(minY, eMag);
146 maxY =
max(maxY, eMag);
149 else if (
mag(eVec & z) > 1-edgeTol)
151 minZ =
min(minZ, eMag);
152 maxZ =
max(maxZ, eMag);
157 minOther =
min(minOther, eMag);
158 maxOther =
max(maxOther, eMag);
163 <<
"Mesh edge statistics:" <<
nl
164 <<
" x aligned : number:" << nX <<
"\tminLen:" << minX
165 <<
"\tmaxLen:" << maxX <<
nl
166 <<
" y aligned : number:" << nY <<
"\tminLen:" << minY
167 <<
"\tmaxLen:" << maxY <<
nl
168 <<
" z aligned : number:" << nZ <<
"\tminLen:" << minZ
169 <<
"\tmaxLen:" << maxZ <<
nl
170 <<
" other : number:" << mesh.
nEdges() - nX - nY - nZ
171 <<
"\tminLen:" << minOther
172 <<
"\tmaxLen:" << maxOther <<
nl <<
endl;
174 if (excludeCmpt == 0)
176 return min(minY,
min(minZ, minOther));
178 else if (excludeCmpt == 1)
180 return min(minX,
min(minZ, minOther));
182 else if (excludeCmpt == 2)
184 return min(minX,
min(minY, minOther));
188 return min(minX,
min(minY,
min(minZ, minOther)));
215 emptyPolyPatch::typeName
248 void selectCurvatureCells
253 const scalar nearDist,
254 const scalar curvature,
286 void addCutNeighbours
289 const bool selectInside,
290 const bool selectOutside,
302 const label celli = iter.key();
307 const label facei = cFaces[i];
318 if (selectInside && inside.
found(nbr))
320 addCutFaces.insert(nbr);
322 else if (selectOutside && outside.
found(nbr))
324 addCutFaces.insert(nbr);
330 Info<<
" Selected an additional " << addCutFaces.size()
331 <<
" neighbours of cutCells to refine" <<
endl;
335 cutCells.
insert(iter.key());
344 bool limitRefinementLevel
347 const label limitDiff,
356 if (!excludeCells.
found(celli))
362 label nbr = cCells[i];
364 if (!excludeCells.
found(nbr))
366 if (refLevel[celli] - refLevel[nbr] >= limitDiff)
369 <<
"Level difference between neighbouring cells "
370 << celli <<
" and " << nbr
371 <<
" greater than or equal to " << limitDiff <<
endl
372 <<
"refLevels:" << refLevel[celli] <<
' '
386 const label celli = iter.key();
391 const label nbr = cCells[i];
393 if (!excludeCells.
found(nbr) && !cutCells.
found(nbr))
395 if (refLevel[celli] + 1 - refLevel[nbr] >= limitDiff)
397 addCutCells.insert(nbr);
403 if (addCutCells.size())
407 Info<<
"Added an additional " << addCutCells.size() <<
" cells"
408 <<
" to satisfy 1:" << limitDiff <<
" refinement level"
413 cutCells.
insert(iter.key());
419 Info<<
"Added no additional cells"
420 <<
" to satisfy 1:" << limitDiff <<
" refinement level"
454 for (
label celli = oldCells; celli < mesh.
nCells(); celli++)
461 forAll(addedCells, oldCelli)
463 const labelList& added = addedCells[oldCelli];
469 label masterLevel = ++refLevel[oldCelli];
473 refLevel[added[i]] = masterLevel;
484 const label writeMesh,
495 Info<<
"Mesh has:" << mesh.
nCells() <<
" cells."
496 <<
" Removing:" << cellLabels.size() <<
" cells" <<
endl;
499 labelList exposedFaces(cellRemover.getExposedFaces(cellLabels));
502 cellRemover.setRefinement
518 cellRemover.topoChange(map());
521 const labelList& cellMap = map().cellMap();
527 newRefLevel[i] = refLevel[cellMap[i]];
531 refLevel.transfer(newRefLevel);
535 Info<<
"Writing refined mesh to time " << runTime.
name() <<
nl
561 const bool selectCut,
562 const bool selectInside,
563 const bool selectOutside,
565 const label nCutLayers,
592 selected.
addSet(cutCells);
603 Info<<
"Determined cell status:" <<
endl
604 <<
" inside :" << inside.
size() <<
endl
605 <<
" outside :" << outside.
size() <<
endl
606 <<
" cutCells:" << cutCells.
size() <<
endl
607 <<
" selected:" << selected.
size() <<
endl
610 writeSet(inside,
"inside cells");
611 writeSet(outside,
"outside cells");
612 writeSet(cutCells,
"cut cells");
613 writeSet(selected,
"selected cells");
618 int main(
int argc,
char *argv[])
627 label newPatchi = addPatch(mesh,
"oldInternalFaces");
634 Info<<
"Checking for motionProperties\n" <<
endl;
648 if (motionObj.headerOk())
650 Info<<
"Reading " << runTime.
constant() /
"motionProperties"
655 Switch twoDMotion(motionProperties.lookup(
"twoDMotion"));
668 "autoRefineMeshDict",
685 scalar minEdgeLen(refineDict.
lookup<scalar>(
"minEdgeLen"));
686 scalar maxEdgeLen(refineDict.
lookup<scalar>(
"maxEdgeLen"));
687 scalar curvature(refineDict.
lookup<scalar>(
"curvature"));
688 scalar curvDist(refineDict.
lookup<scalar>(
"curvatureDistance"));
693 Info<<
"Cells to be used for meshing (0=false, 1=true):" <<
nl
694 <<
" cells cut by surface : " << selectCut <<
nl
695 <<
" cells inside of surface : " << selectInside <<
nl
696 <<
" cells outside of surface : " << selectOutside <<
nl
697 <<
" hanging cells : " << selectHanging <<
nl
701 if (nCutLayers > 0 && selectInside)
704 <<
"Illogical settings : Both nCutLayers>0 and selectInside true."
706 <<
"This would mean that inside cells get removed but should be"
707 <<
" included in final mesh" <<
endl;
720 forAll(outsidePts, outsideI)
722 const point& outsidePoint = outsidePts[outsideI];
724 if (queryMesh.findCell(outsidePoint, -1,
false) == -1)
727 <<
"outsidePoint " << outsidePoint
728 <<
" is not inside any cell"
754 Info<<
"Read existing refinement level from file "
756 <<
" min level : " <<
min(refLevel) <<
nl
757 <<
" max level : " << maxLevel <<
nl
762 Info<<
"Created zero refinement level in file "
771 direction normalDir(getNormalDir(correct2DPtr));
772 scalar meshMinEdgeLen = getEdgeStats(mesh, normalDir);
774 while (meshMinEdgeLen > minEdgeLen)
801 Info<<
" Selected " << cutCells.
name() <<
" with "
802 << cutCells.
size() <<
" cells" <<
endl;
804 if ((curvDist > 0) && (meshMinEdgeLen < maxEdgeLen))
821 Info<<
" Selected " << curveCells.name() <<
" with "
822 << curveCells.size() <<
" cells" <<
endl;
841 cutCells.
subset(curveCells);
843 Info<<
" Removed cells not on surface curvature. Selected "
852 subsetMesh(mesh, writeMesh, newPatchi, inside, cutCells, refLevel);
872 <<
" Selected for refinement :" << cutCells.
size() <<
nl
875 if (cutCells.
empty())
877 Info<<
"Stopping refining since 0 cells selected to be refined ..."
882 if ((mesh.
nCells() + 8*cutCells.
size()) > cellLimit)
884 Info<<
"Stopping refining since cell limit reached ..." <<
nl
885 <<
"Would refine from " << mesh.
nCells()
886 <<
" to " << mesh.
nCells() + 8*cutCells.
size() <<
" cells."
904 Info<<
" Writing refined mesh to time " << runTime.
name()
913 meshMinEdgeLen = getEdgeStats(mesh, normalDir);
949 Info<<
"Detected " << hanging.
size() <<
" hanging cells"
950 <<
" (cells with all points on"
951 <<
" outside of cellSet 'selected').\nThese would probably result"
952 <<
" in flattened cells when snapping the mesh to the surface"
955 Info<<
"Refining " << hanging.
size() <<
" hanging cells" <<
nl
972 doRefinement(mesh, refineDict, hanging, refLevel);
974 Info<<
"Writing refined mesh to time " << runTime.
name() <<
nl
985 Info<<
"Writing refined mesh to time " << runTime.
name() <<
nl
#define forAll(list, i)
Loop across all elements in list.
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
static twoDPointCorrector & New(const word &name, const polyMesh &mesh)
Construct and return the named DemandDrivenMeshObject.
bool insert(const Key &key)
Insert a new entry.
List< Key > toc() const
Return the table of contents.
label size() const
Return number of elements in table.
bool empty() const
Return true if the hash table is empty.
bool found(const Key &) const
Return true if hashedEntry is found in table.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
fileName relativeObjectPath() const
Return complete relativePath + object name.
const word & name() const
Return name.
static unsigned int defaultPrecision()
Return the default precision.
void size(const label)
Override size to be inconsistent with allocated storage.
void setSize(const label)
Reset size of List.
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
static const word & constant()
Return constant name.
static const word & system()
Return system name.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
label size() const
Return the number of elements in the UList.
label size() const
Return the number of elements in the UPtrList.
static void noParallel()
Remove the parallel options.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
A bounding box defined in terms of the points at its extremities.
A collection of cell labels.
virtual void topoChange(const polyTopoChangeMap &map)
Update any stored data for new labels.
A list of keyword definitions, which are a keyword followed by any number of values (e....
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
const word & name() const
Return const reference to name.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Empty front and back plane patch. Used for 2-D geometries.
A class for handling file names.
Various (local, not parallel) searches on polyMesh; uses (demand driven) octree to search.
Does multiple pass refinement to refine cells in multiple directions.
const Time & time() const
Return time.
label findIndex(const word &patchName) const
Find patch index given a name.
Mesh consisting of general polyhedral cells.
static word defaultRegion
Return the default region name.
virtual const labelList & faceOwner() const
Return face owner.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
virtual const labelList & faceNeighbour() const
Return face neighbour.
void addPatches(const List< polyPatch * > &, const bool validBoundary=true)
Add boundary patches.
void removeBoundary()
Remove boundary patches.
A patch is a list of labels that address the faces in the global face list.
label start() const
Return start label of this patch in the polyMesh face list.
virtual autoPtr< polyPatch > clone(const polyBoundaryMesh &bm) const
Construct and return a clone, resetting the boundary mesh.
Direct mesh changes based on v1.3 polyTopoChange syntax.
Cell-face mesh analysis engine.
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
const vectorField & cellCentres() const
const labelListList & cellCells() const
label nInternalFaces() const
virtual const pointField & points() const =0
Return mesh points.
bool isInternalFace(const label faceIndex) const
Return true if given face label is internal to the mesh.
const cellList & cells() const
virtual bool write(const bool write=true) const
Write using setting from DB.
Given list of cells to remove insert all the topology changes.
string & expand(const bool allowEmpty=false)
Expand initial tildes and all occurrences of environment variables.
static labelHashSet getHangingCells(const primitiveMesh &mesh, const labelHashSet &internalCells)
Get cells using points on 'outside' only.
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.
A topoSetSource to select cells based on relation to surface.
virtual void addSet(const topoSet &set)
Add elements present in set.
virtual void subset(const topoSet &set)
Subset contents. Only elements present in both sets remain.
Helper class to search on triSurface.
const triSurface & surface() const
Return reference to the surface.
Triangulated surface description with patch information.
Class applies a two-dimensional correction to mesh motion point field.
const vector & planeNormal() const
Return plane normal.
Templated form of IOobject providing type information for file reading and header type checking.
A class for handling words, derived from string.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
int main(int argc, char *argv[])
const fvPatchList & patches
#define WarningInFunction
Report a warning using Foam::Warning.
errorManipArg< error, int > exit(error &err, const int errNo=1)
List< label > labelList
A List of labels.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Ostream & endl(Ostream &os)
Add newline and flush stream.
errorManip< error > abort(error &err)
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
Vector< scalar > vector
A scalar version of the templated Vector.
dimensioned< scalar > mag(const dimensioned< Type > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys.