Go to the documentation of this file.
52 bool Foam::hexCellLooper::walkHex
55 const label startFacei,
56 const label startEdgeI,
62 label facei = startFacei;
64 label edgeI = startEdgeI;
72 Pout<<
" walkHex : inserting cut onto edge:" << edgeI
78 loopWeights[cutI] = 0.5;
88 if (edgeI == startEdgeI)
98 Pout<<
"hexCellLooper::walkHex" <<
"Problem : cell:" << celli
99 <<
" collected loop:";
101 Pout<<
"loopWeights:" << loopWeights <<
endl;
111 void Foam::hexCellLooper::makeFace
120 facePoints.setSize(loop.size());
121 faceVerts.setSize(loop.size());
137 loopWeights[cutI]*v1 + (1.0-loopWeights[cutI])*v0;
146 faceVerts[cutI] = cutI;
196 loopWeights.setSize(4);
198 success = walkHex(celli, face0, edgeI, loop, loopWeights);
220 <<
"could not cut cell " << celli <<
endl;
222 fileName cutsFile(
"hexCellLooper_" +
name(celli) +
".obj");
224 Pout<<
"hexCellLooper : writing cell to " << cutsFile <<
endl;
245 label elem = loop[elemI];
247 if (loopSet.found(elem))
252 loopSet.insert(elem);
256 face faceVerts(loop.size());
259 makeFace(loop, loopWeights, faceVerts, facePoints);
261 if ((faceVerts.mag(facePoints) < SMALL) || (loop.size() < 3))
264 <<
" on points:" << facePoints <<
endl
276 const plane& cutPlane,
int debug
Static debugging option.
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
List< label > labelList
A List of labels.
vectorField pointField
pointField is a vectorField.
virtual const pointField & points() const
Return raw points.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
A class for handling file names.
virtual bool cut(const vector &refDir, const label celli, const boolList &vertIsCut, const boolList &edgeIsCut, const scalarField &edgeWeight, labelList &loop, scalarField &loopWeights) const
Create cut along circumference of celli. Gets current mesh cuts.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
prefixOSstream Pout
An Ostream wrapper for parallel output to std::cout.
virtual ~hexCellLooper()
Destructor.
Mesh consisting of general polyhedral cells.
#define forAll(list, i)
Loop across all elements in list.
Geometric class that creates a 3D plane and can return the intersection point between a line and the ...
static label edgeToEVert(const primitiveMesh &mesh, const label edgeI)
Convert edgeI to eVert.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
word name(const complex &c)
Return string representation of complex.
Implementation of cellLooper. Does pure geometric cut through cell.
Patchify triangles based on orientation w.r.t other (triangulated or triangulatable) surfaces.
Macros for easy insertion into run-time selection tables.
errorManip< error > abort(error &err)
Output to file stream, using an OSstream.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const dimensionedScalar e
Elementary charge.
A List with indirect addressing.
A face is a list of labels corresponding to mesh vertices.
Various functions to operate on Lists.
Maps a geometry to a set of cell primitives.
virtual bool cut(const vector &refDir, const label celli, const boolList &vertIsCut, const boolList &edgeIsCut, const scalarField &edgeWeight, labelList &loop, scalarField &loopWeights) const
Create cut along circumference of celli. Gets current mesh cuts.
vector point
Point is a vector.
Ostream & writeCuts(Ostream &os, const labelList &, const scalarField &) const
Write cut descriptions to Ostream.
void setSize(const label newSize)
Alias for resize(const label)
defineTypeNameAndDebug(combustionModel, 0)
#define WarningInFunction
Report a warning using Foam::Warning.
const polyMesh & mesh() const