Go to the documentation of this file.
34 template<
class EdgeOrientIntersect,
class EdgeAlphaIntersect>
39 const EdgeOrientIntersect& edgeOrientIntersect,
40 const EdgeAlphaIntersect& edgeAlphaIntersect,
41 const bool triangulate,
51 const label nCellCuts =
cellCuts.count();
60 nFaceCuts = 4*nCellCuts;
93 label unwindPoint = 0;
96 const auto unwindWalk =
97 [&](
const label failedCellId = -1) ->
void
100 dynCutPoints.
resize(unwindPoint);
103 endPoints.erase(localEndPoints);
106 if (failedCellId != -1)
108 failedCells.
insert(failedCellId);
121 localFaceLoop.
clear();
122 localEndPoints.
clear();
124 unwindPoint = dynCutPoints.size();
129 unsigned pointCutType = 0u;
131 for (
const label facei : cFace)
133 const face&
f = faces[facei];
140 if (!edgeOrientIntersect(
e))
158 label cutPointId = handledEdges.
lookup(
e, -1);
163 localEdges.
insert(
e, cutPointId);
168 cutPointId = dynCutPoints.size();
174 const scalar
alpha = edgeAlphaIntersect(
e);
180 const label endp =
e.first();
182 if (endPoints.insert(endp, cutPointId))
184 localEndPoints.
insert(endp);
189 cutPointId = endPoints[endp];
192 else if (
alpha >= (1.0 - SMALL))
196 const label endp =
e.last();
198 if (endPoints.insert(endp, cutPointId))
200 localEndPoints.
insert(endp);
205 cutPointId = endPoints[endp];
216 localEdges.
insert(
e, cutPointId);
225 const label nTargetLoop = localFaces.
size();
240 if (pointCutType == 1 || pointCutType == 2)
246 <<
"skip duplicate on-place cut for cell " << celli
247 <<
" type (" << pointCutType <<
")" <<
endl;
261 const edge refEdge = localFaces.
begin().key();
262 label nextFace = localFaces.
begin().val()[0];
265 <<
"search face " << nextFace <<
" IN " << localEdges <<
endl;
270 localFaces.
size() && loopi < 2*nTargetLoop;
279 <<
"lookup " << nextFace <<
" in " << iter.val() <<
nl;
282 const label got = iter.val().which(nextFace);
289 nextFace = iter.val()[(got?0:1)];
292 localFaceLoop.
append(localEdges[iter.key()]);
295 <<
" faces " << iter.val()
296 <<
" point " << localFaceLoop.last()
297 <<
" edge=" << iter.key() <<
" nextFace=" << nextFace
301 localFaces.
erase(iter);
315 if (nTargetLoop != localFaceLoop.size())
318 <<
"Warn expected " << nTargetLoop <<
" but got "
319 << localFaceLoop.size() <<
endl;
327 handledEdges += localEdges;
329 face f(localFaceLoop);
332 if ((
f.areaNormal(dynCutPoints) & refEdge.
vec(
points)) < 0)
340 label nTri =
f.triangles(dynCutPoints, dynCutFaces);
343 dynCutCells.
append(celli);
349 dynCutCells.
append(celli);
353 if (failedCells.
size())
356 <<
"Failed cuts for " << failedCells.
size() <<
" cells:" <<
nl
363 if (dynCutCells.empty())
365 this->storedPoints().clear();
366 this->storedFaces().clear();
371 this->storedPoints().transfer(dynCutPoints);
372 this->storedFaces().transfer(dynCutFaces);
373 meshCells_.transfer(dynCutCells);
label size() const noexcept
The number of elements in table.
A bitSet stores bits (elements with only two states) in packed internal format and supports a variety...
List< T > values(const HashTable< T, Key, Hash > &tbl, const bool doSort=false)
List of values from HashTable, optionally sorted.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
bool insert(const edge &key, const T &obj)
Copy insert a new entry, not overwriting existing entries.
Ostream & endl(Ostream &os)
Add newline and flush stream.
A HashTable with keys but without contents that is similar to std::unordered_set.
#define forAll(list, i)
Loop across all elements in list.
bool erase(const iterator &iter)
Erase an entry specified by given iterator.
DynamicList< T, SizeMin > & append(const T &val)
Append an element to the end of this list.
void clear()
Clear the addressed list, i.e. set the size to zero.
#define forAllIters(container, iter)
Iterate across all elements in the container object.
const T & lookup(const Key &key, const T &deflt) const
Return hashed entry if it exists, or return the given default.
List< Key > sortedToc() const
The table of contents (the keys) in sorted order.
void walkCellCuts(const primitiveMesh &mesh, const bitSet &cellCuts, const EdgeOrientIntersect &edgeOrientIntersect, const EdgeAlphaIntersect &edgeAlphaIntersect, const bool triangulate, label nFaceCuts=0)
Walk cell cuts to create faces.
void resize(const label nElem)
Alter addressable list size.
iterator begin()
iterator set to the beginning of the HashTable
#define DebugInfo
Report an information message using Foam::Info.
FlatOutput< Container > flatOutput(const Container &obj, label len=0)
Global flatOutput function.
void clear()
Clear all entries from table.
vector vec(const UList< point > &pts) const
Return the vector (end - start)
const dimensionedScalar e
Elementary charge.
bool insert(const Key &key)
Insert a new entry, not overwriting existing entries.
A face is a list of labels corresponding to mesh vertices.
const volScalarField & p0
bool found(const Key &key) const
Return true if hashed entry is found in table.
#define WarningInFunction
Report a warning using Foam::Warning.
A cell is defined as a list of faces with extra functionality.
Description of cuts across cells.
Cell-face mesh analysis engine.