51bool Foam::hexCellLooper::walkHex
54 const label startFacei,
55 const label startEdgeI,
58 scalarField& loopWeights
61 label facei = startFacei;
63 label edgeI = startEdgeI;
71 Pout<<
" walkHex : inserting cut onto edge:" << edgeI
77 loopWeights[cutI] = 0.5;
87 if (edgeI == startEdgeI)
97 Pout<<
"hexCellLooper::walkHex" <<
"Problem : cell:" << celli
98 <<
" collected loop:";
100 Pout<<
"loopWeights:" << loopWeights <<
endl;
110void Foam::hexCellLooper::makeFace
119 facePoints.setSize(loop.size());
120 faceVerts.setSize(loop.size());
124 label
cut = loop[cutI];
128 label edgeI = getEdge(
cut);
136 loopWeights[cutI]*v1 + (1.0-loopWeights[cutI])*v0;
140 label vertI = getVertex(
cut);
145 faceVerts[cutI] = cutI;
190 success = walkHex(celli, face0, edgeI, loop, loopWeights);
212 <<
"could not cut cell " << celli <<
endl;
214 fileName cutsFile(
"hexCellLooper_" +
name(celli) +
".obj");
216 Pout<<
"hexCellLooper : writing cell to " << cutsFile <<
endl;
237 label elem = loop[elemI];
239 if (loopSet.
found(elem))
251 makeFace(loop, loopWeights, faceVerts, facePoints);
253 if ((faceVerts.
mag(facePoints) < SMALL) || (loop.
size() < 3))
256 <<
" on points:" << facePoints <<
endl
268 const plane& cutPlane,
Various functions to operate on Lists.
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
bool insert(const Key &key)
Insert a new entry, not overwriting existing entries.
bool found(const Key &key) const
Return true if hashed entry is found in table.
void setSize(const label n)
Alias for resize()
Output to file stream, using an OSstream.
A List with indirect addressing. Like IndirectList but does not store addressing.
bool empty() const noexcept
True if the UList is empty (ie, size() is zero)
iterator end() noexcept
Return an iterator to end traversing the UList.
void size(const label n)
Older name for setAddressableSize.
Abstract base class. Concrete implementations know how to cut a cell (i.e. determine a loop around th...
Maps a geometry to a set of cell primitives.
Ostream & writeCuts(Ostream &os, const labelList &, const scalarField &) const
Write cut descriptions to Ostream.
static label edgeToEVert(const primitiveMesh &mesh, const label edgeI)
Convert edgeI to eVert.
const polyMesh & mesh() const
A face is a list of labels corresponding to mesh vertices.
scalar mag(const UList< point > &p) const
Magnitude of face area.
A class for handling file names.
Implementation of cellLooper. Does pure geometric cut through cell.
Implementation of cellLooper.
Geometric class that creates a 3D plane and can return the intersection point between a line and the ...
Mesh consisting of general polyhedral cells.
virtual const pointField & points() const
Return raw points.
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
Patchify triangles based on orientation w.r.t other (triangulated or triangulatable) surfaces.
A class for handling words, derived from Foam::string.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
#define WarningInFunction
Report a warning using Foam::Warning.
List< label > labelList
A List of labels.
vectorField pointField
pointField is a vectorField.
vector point
Point is a vector.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
errorManip< error > abort(error &err)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
#define forAll(list, i)
Loop across all elements in list.