52 const face&
f = pFaces[fI];
54 label oCI = pOwner[fI];
56 const point& oCc = pC[oCI];
62 scalar thisBaseMinTetQuality = VGREAT;
64 const point& tetBasePt = pPts[f[faceBasePtI]];
66 for (
label tetPtI = 1; tetPtI < f.size() - 1; tetPtI++)
68 label facePtI = (tetPtI + faceBasePtI) % f.
size();
69 label otherFacePtI = f.fcIndex(facePtI);
73 label ptAI = f[facePtI];
74 label ptBI = f[otherFacePtI];
76 const point& pA = pPts[ptAI];
77 const point& pB = pPts[ptBI];
81 tetQualities[0] = tet.
quality();
86 label ptAI = f[otherFacePtI];
87 label ptBI = f[facePtI];
89 const point& pA = pPts[ptAI];
90 const point& pB = pPts[ptBI];
94 tetQualities[1] = tet.
quality();
97 if (
min(tetQualities) < thisBaseMinTetQuality)
99 thisBaseMinTetQuality =
min(tetQualities);
103 if (thisBaseMinTetQuality > tol)
123 return findSharedBasePoint
147 const face&
f = pFaces[fI];
149 label cI = pOwner[fI];
151 bool own = (pOwner[fI] == cI);
153 const point& cC = pC[cI];
157 scalar thisBaseMinTetQuality = VGREAT;
159 const point& tetBasePt = pPts[f[faceBasePtI]];
161 for (
label tetPtI = 1; tetPtI < f.size() - 1; tetPtI++)
163 label facePtI = (tetPtI + faceBasePtI) % f.
size();
164 label otherFacePtI = f.fcIndex(facePtI);
172 ptBI = f[otherFacePtI];
176 ptAI = f[otherFacePtI];
180 const point& pA = pPts[ptAI];
181 const point& pB = pPts[ptBI];
185 scalar tetQuality = tet.
quality();
187 if (tetQuality < thisBaseMinTetQuality)
189 thisBaseMinTetQuality = tetQuality;
193 if (thisBaseMinTetQuality > tol)
222 for (
label fI = 0; fI < nInternalFaces; fI++)
224 tetBasePtIs[fI] = findSharedBasePoint(mesh, fI, tol, report);
229 for(
label faceI = nInternalFaces; faceI < mesh.
nFaces(); faceI++)
231 neighbourCellCentres[faceI - nInternalFaces] =
235 syncTools::swapBoundaryFacePositions(mesh, neighbourCellCentres);
242 mesh.
nFaces() - nInternalFaces,
248 label fI = nInternalFaces, bFI = 0;
256 if (patches[patchI].coupled())
259 refCast<const coupledPolyPatch>(patches[patchI]);
263 boundaryFaceTetBasePtIs[bFI] = findSharedBasePoint
267 neighbourCellCentres[bFI],
276 boundaryFaceTetBasePtIs[bFI] = -2;
281 boundaryFaceTetBasePtIs[bFI] = findBasePoint
294 syncTools::syncBoundaryFaceList
297 boundaryFaceTetBasePtIs,
303 label fI = nInternalFaces, bFI = 0;
308 label& bFTetBasePtI = boundaryFaceTetBasePtIs[bFI];
310 if (bFTetBasePtI == -2)
315 "polyMeshTetDecomposition::findFaceBasePts" 317 "const polyMesh& mesh, " 321 ) <<
"Coupled face base point exchange failure for face " 326 if (bFTetBasePtI < 1)
336 if (patches[patchI].coupled())
339 refCast<const coupledPolyPatch>(patches[patchI]);
365 bFTetBasePtI = mesh.
faces()[fI].
size() - bFTetBasePtI;
397 syncTools::swapBoundaryFacePositions(mesh, neiCc);
403 label nErrorTets = 0;
407 const face&
f = fcs[faceI];
434 const face& f = fcs[faceI];
458 if (findSharedBasePoint(mesh, faceI, tol, report) == -1)
472 if (patches[patchI].coupled())
496 if (findBasePoint(mesh, faceI, tol, report) == -1)
515 Info<<
" ***Error in face tets: " 516 << nErrorTets <<
" faces with low quality or negative volume " 517 <<
"decomposition tets." <<
endl;
541 static label nWarnings = 0;
542 static const label maxWarnings = 100;
547 const face&
f = pFaces[fI];
553 bool own = (pOwner[fI] == cI);
557 if (tetBasePtI == -1)
559 if (nWarnings < maxWarnings)
564 "polyMeshTetDecomposition::faceTetIndices" 570 ) <<
"No base point for face " << fI <<
", " << f
571 <<
", produces a valid tet decomposition." 575 if (nWarnings == maxWarnings)
577 Warning<<
"Suppressing any further warnings." <<
endl;
584 for (
label tetPtI = 1; tetPtI < f.
size() - 1; tetPtI++)
588 label facePtI = (tetPtI + tetBasePtI) % f.
size();
591 faceTetIs.
cell() = cI;
593 faceTetIs.
face() = fI;
600 faceTetIs.
facePtB() = otherFacePtI;
604 faceTetIs.
facePtA() = otherFacePtI;
608 faceTetIs.
tetPt() = tetPtI;
623 static label nWarnings = 0;
624 static const label maxWarnings = 100;
629 if (tetBasePtI == -1)
631 if (nWarnings < maxWarnings)
636 "polyMeshTetDecomposition::triangleTetIndices" 643 ) <<
"No base point for face " << fI <<
", " << f
644 <<
", produces a valid tet decomposition." 648 if (nWarnings == maxWarnings)
650 Warning<<
"Suppressing any further warnings." <<
endl;
659 label facePtI = (tetPtI + tetBasePtI) % f.
size();
662 faceTetIs.
cell() = cI;
664 faceTetIs.
face() = fI;
671 faceTetIs.
facePtB() = otherFacePtI;
675 faceTetIs.
facePtA() = otherFacePtI;
679 faceTetIs.
tetPt() = tetPtI;
694 const cell& thisCell = pCells[cI];
700 nTets += pFaces[thisCell[cFI]].
size() - 2;
707 label fI = thisCell[cFI];
709 cellTets.
append(faceTetIndices(mesh, fI, cI));
726 const cell& thisCell = pCells[cI];
733 label fI = thisCell[cFI];
734 const face&
f = pFaces[fI];
736 for (
label tetPtI = 1; tetPtI < f.
size() - 1; tetPtI++)
753 tetContainingPt = faceTetIs;
758 if (tetContainingPt.
cell() != -1)
764 return tetContainingPt;
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
label cell() const
Return the cell.
static tetIndices triangleTetIndices(const polyMesh &mesh, label fI, label cI, const label tetPtI)
Return the tet decomposition of the given triangle of the given face.
static tetIndices findTet(const polyMesh &mesh, label cI, const point &pt)
Find the tet decomposition of the cell containing the given point.
static label findSharedBasePoint(const polyMesh &mesh, label fI, const point &nCc, scalar tol, bool report=false)
Find the first suitable base point to use for a minimum.
A cell is defined as a list of faces with extra functionality.
static labelList findFaceBasePts(const polyMesh &mesh, scalar tol=minTetQuality, bool report=false)
Find a suitable base point for each face for decomposition.
const labelList & tetBasePtIs() const
Return the tetBasePtIs.
label faceBasePt() const
Return the face base point.
static label findBasePoint(const polyMesh &mesh, label fI, scalar tol, bool report=false)
Find the base point to use for a minimum triangle.
static const scalar minTetQuality
Minimum tetrahedron quality.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
static bool checkFaceTets(const polyMesh &mesh, scalar tol=minTetQuality, const bool report=false, labelHashSet *setPtr=NULL)
Check face-decomposition tet volume.
void size(const label)
Override size to be inconsistent with allocated storage.
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
label fcIndex(const label i) const
Return the forward circular index, i.e. the next index.
const cellList & cells() const
bool inside(const point &pt) const
Return true if point is inside tetrahedron.
const vectorField & cellCentres() const
static List< tetIndices > faceTetIndices(const polyMesh &mesh, label fI, label cI)
Return the tet decomposition of the given face, with.
label nextLabel(const label i) const
Next vertex on face.
tetrahedron< point, const point & > tetPointRef
scalar quality() const
Return quality: Ratio of tetrahedron and circum-sphere.
tetPointRef tet(const polyMesh &mesh) const
Return the geometry corresponding to this tet from the.
A face is a list of labels corresponding to mesh vertices.
Ostream & endl(Ostream &os)
Add newline and flush stream.
#define WarningIn(functionName)
Report a warning using Foam::Warning.
virtual const pointField & points() const
Return raw points.
errorManip< error > abort(error &err)
virtual const labelList & faceOwner() const
Return face owner.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
label nInternalFaces() const
const labelList & patchID() const
Per boundary face label the patch index.
Mesh consisting of general polyhedral cells.
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
label facePtA() const
Return face point A.
virtual bool owner() const =0
Does this side own the patch ?
A List obtained as a section of another List.
virtual const faceList & faces() const
Return raw faces.
label tetPt() const
Return the characterising tetPtI.
Info<< "Finished reading KIVA file"<< endl;cellShapeList cellShapes(nPoints);labelList cellZoning(nPoints,-1);const cellModel &hex=*(cellModeller::lookup("hex"));labelList hexLabels(8);label activeCells=0;labelList pointMap(nPoints);forAll(pointMap, i){pointMap[i]=i;}for(label i=0;i< nPoints;i++){if(f[i] > 0.0){hexLabels[0]=i;hexLabels[1]=i1tab[i];hexLabels[2]=i3tab[i1tab[i]];hexLabels[3]=i3tab[i];hexLabels[4]=i8tab[i];hexLabels[5]=i1tab[i8tab[i]];hexLabels[6]=i3tab[i1tab[i8tab[i]]];hexLabels[7]=i3tab[i8tab[i]];cellShapes[activeCells]=cellShape(hex, hexLabels);edgeList edges=cellShapes[activeCells].edges();forAll(edges, ei){if(edges[ei].mag(points)< SMALL){label start=pointMap[edges[ei].start()];while(start!=pointMap[start]){start=pointMap[start];}label end=pointMap[edges[ei].end()];while(end!=pointMap[end]){end=pointMap[end];}label minLabel=min(start, end);pointMap[start]=pointMap[end]=minLabel;}}cellZoning[activeCells]=idreg[i];activeCells++;}}cellShapes.setSize(activeCells);cellZoning.setSize(activeCells);forAll(cellShapes, celli){cellShape &cs=cellShapes[celli];forAll(cs, i){cs[i]=pointMap[cs[i]];}cs.collapse();}label bcIDs[11]={-1, 0, 2, 4,-1, 5,-1, 6, 7, 8, 9};const label nBCs=12;const word *kivaPatchTypes[nBCs]={&wallPolyPatch::typeName,&wallPolyPatch::typeName,&wallPolyPatch::typeName,&wallPolyPatch::typeName,&symmetryPolyPatch::typeName,&wedgePolyPatch::typeName,&polyPatch::typeName,&polyPatch::typeName,&polyPatch::typeName,&polyPatch::typeName,&symmetryPolyPatch::typeName,&oldCyclicPolyPatch::typeName};enum patchTypeNames{PISTON, VALVE, LINER, CYLINDERHEAD, AXIS, WEDGE, INFLOW, OUTFLOW, PRESIN, PRESOUT, SYMMETRYPLANE, CYCLIC};const char *kivaPatchNames[nBCs]={"piston","valve","liner","cylinderHead","axis","wedge","inflow","outflow","presin","presout","symmetryPlane","cyclic"};List< SLList< face > > pFaces[nBCs]
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
The coupledPolyPatch is an abstract base class for patches that couple regions of the computational d...
virtual const labelList & faceNeighbour() const
Return face neighbour.
Storage and named access for the indices of a tet which is part of the decomposition of a cell...
static List< tetIndices > cellTetIndices(const polyMesh &mesh, label cI)
Return the tet decomposition of the given cell, see.
label facePtB() const
Return face point B.
bool isInternalFace(const label faceIndex) const
Return true if given face label is internal to the mesh.
const vectorField & faceCentres() const
label face() const
Return the face.
bool insert(const Key &key)
Insert a new entry.