39 void Foam::shortEdgeFilter2D::assignBoundaryPointRegions
41 List<DynamicList<label>>& boundaryPointRegions
46 const edge&
e = iter.key();
47 const label regi = iter.val();
49 boundaryPointRegions[
e.start()].appendUniq(regi);
50 boundaryPointRegions[
e.end()].appendUniq(regi);
55 void Foam::shortEdgeFilter2D::updateEdgeRegionMap
57 const MeshedSurface<face>& surfMesh,
58 const List<DynamicList<label>>& boundaryPtRegions,
60 EdgeMap<label>& mapEdgesRegion,
64 EdgeMap<label> newMapEdgesRegion(mapEdgesRegion.size());
66 const edgeList& edges = surfMesh.edges();
67 const labelList& meshPoints = surfMesh.meshPoints();
69 patchSizes.
setSize(patchNames_.size(), 0);
74 if (surfMesh.isInternalEdge(edgeI))
79 const edge&
e = edges[edgeI];
81 const label startI = meshPoints[
e[0]];
82 const label endI = meshPoints[
e[1]];
86 const DynamicList<label> startPtRegions =
87 boundaryPtRegions[surfPtToBoundaryPt[startI]];
88 const DynamicList<label> endPtRegions =
89 boundaryPtRegions[surfPtToBoundaryPt[endI]];
91 if (startPtRegions.size() > 1 && endPtRegions.size() > 1)
93 region = startPtRegions[0];
96 <<
"Both points in edge are in different regions." 97 <<
" Assigning edge to region " << region
100 else if (startPtRegions.size() > 1 || endPtRegions.size() > 1)
104 startPtRegions.size() > 1
111 startPtRegions[0] == endPtRegions[0]
112 && startPtRegions[0] != -1
115 region = startPtRegions[0];
120 newMapEdgesRegion.insert(
e, region);
121 patchSizes[region]++;
125 mapEdgesRegion.transfer(newMapEdgesRegion);
131 Foam::shortEdgeFilter2D::shortEdgeFilter2D
134 const dictionary&
dict 137 shortEdgeFilterFactor_(
dict.
get<scalar>(
"shortEdgeFilterFactor")),
138 edgeAttachedToBoundaryFactor_
140 dict.getOrDefault<scalar>(
"edgeAttachedToBoundaryFactor", 2.0)
168 OFstream str(
"indirectPatchEdges.obj");
171 Info<<
"Writing indirectPatchEdges to " << str.name() <<
endl;
175 const edge&
e = iter.key();
189 ms_ = MeshedSurface<face>(std::move(
points), std::move(faces));
191 Info<<
"Meshed surface stats before edge filtering :" <<
endl;
192 ms_.writeStats(
Info);
198 ms_.write(
"MeshedSurface_preFilter.obj");
209 const edgeList& edges = ms_.edges();
210 const faceList& faces = ms_.surfFaces();
211 const labelList& meshPoints = ms_.meshPoints();
212 const labelList& boundaryPoints = ms_.boundaryPoints();
215 label nPointsToRemove = 0;
217 labelList pointsToRemove(ms_.points().size(), -1);
220 labelList newFaceVertexCount(faces.size(), -1);
223 newFaceVertexCount[facei] = faces[facei].size();
228 List<DynamicList<label>> boundaryPointRegions
233 assignBoundaryPointRegions(boundaryPointRegions);
237 Info<<
" Marking edges attached to boundaries." <<
endl;
238 boolList edgeAttachedToBoundary(edges.size(),
false);
241 const edge&
e = edges[edgeI];
242 const label startVertex =
e.start();
243 const label endVertex =
e.end();
245 forAll(boundaryPoints, bPoint)
249 boundaryPoints[bPoint] == startVertex
250 || boundaryPoints[bPoint] == endVertex
253 edgeAttachedToBoundary[edgeI] =
true;
260 const edge&
e = edges[edgeI];
263 const label startVertex =
e.start();
264 const label endVertex =
e.end();
269 points[meshPoints[startVertex]]
270 -
points[meshPoints[endVertex]]
273 if (edgeAttachedToBoundary[edgeI])
275 edgeLength *= edgeAttachedToBoundaryFactor_;
278 scalar shortEdgeFilterValue = 0.0;
280 const labelList& psEdges = ms_.pointEdges()[startVertex];
281 const labelList& peEdges = ms_.pointEdges()[endVertex];
285 const edge& psE = edges[psEdges[psEdgeI]];
286 if (edgeI != psEdges[psEdgeI])
288 shortEdgeFilterValue +=
291 points[meshPoints[psE.start()]]
292 -
points[meshPoints[psE.end()]]
299 const edge& peE = edges[peEdges[peEdgeI]];
300 if (edgeI != peEdges[peEdgeI])
302 shortEdgeFilterValue +=
305 points[meshPoints[peE.start()]]
306 -
points[meshPoints[peE.end()]]
311 shortEdgeFilterValue *=
312 shortEdgeFilterFactor_
313 /(psEdges.size() + peEdges.size() - 2);
315 edge lookupInPatchEdge
317 meshPoints[startVertex],
318 meshPoints[endVertex]
323 edgeLength < shortEdgeFilterValue
324 || indirectPatchEdge_.found(lookupInPatchEdge)
327 bool flagDegenerateFace =
false;
332 const face&
f = ms_.localFaces()[
pFaces[pFacei]];
336 if (
f[fp] == endVertex)
339 if (newFaceVertexCount[
pFaces[pFacei]] < 4)
341 flagDegenerateFace =
true;
345 newFaceVertexCount[
pFaces[pFacei]]--;
352 if (newFaceVertexCount[
pFaces[pFacei]] < 3)
354 flagDegenerateFace =
true;
363 pointsToRemove[meshPoints[startVertex]] == -1
364 && pointsToRemove[meshPoints[endVertex]] == -1
365 && !flagDegenerateFace
368 const DynamicList<label>& startVertexRegions =
369 boundaryPointRegions[meshPoints[startVertex]];
370 const DynamicList<label>& endVertexRegions =
371 boundaryPointRegions[meshPoints[endVertex]];
373 if (startVertexRegions.size() && endVertexRegions.size())
375 if (startVertexRegions.size() > 1)
377 pointsToRemove[meshPoints[endVertex]] =
378 meshPoints[startVertex];
382 pointsToRemove[meshPoints[startVertex]] =
383 meshPoints[endVertex];
386 else if (startVertexRegions.size())
388 pointsToRemove[meshPoints[endVertex]] =
389 meshPoints[startVertex];
393 pointsToRemove[meshPoints[startVertex]] =
394 meshPoints[endVertex];
402 label totalNewPoints =
points.size() - nPointsToRemove;
406 label numberRemoved = 0;
409 labelList newPtToOldPt(totalNewPoints, -1);
414 if (pointsToRemove[pointi] == -1)
416 newPoints[pointi - numberRemoved] =
points[pointi];
417 newPointNumbers[pointi] = pointi - numberRemoved;
418 newPtToOldPt[pointi - numberRemoved] = pointi;
431 label newFaceSize = 0;
436 const face&
f = faces[facei];
439 newFace.setSize(
f.size());
444 label pointi =
f[fp];
446 if (pointsToRemove[pointi] == -1)
448 newFace[newFaceSize++] = newPointNumbers[pointi];
452 label
newPointi = pointsToRemove[pointi];
462 label totalChain = 0;
463 for (label nChain = 0; nChain <= totalChain; ++nChain)
465 if (newPointNumbers[pChain] != -1)
467 newFace[newFaceSize++] = newPointNumbers[pChain];
468 newPointNumbers[pointi] = newPointNumbers[pChain];
469 maxChain =
max(totalChain, maxChain);
474 <<
"Point " << pChain
475 <<
" marked for deletion as well as point " 477 <<
" Incrementing maxChain by 1 from " 478 << totalChain <<
" to " << totalChain + 1
482 pChain = pointsToRemove[pChain];
489 newPointNumbers[pointi] = newPointNumbers[
newPointi];
495 newFace.setSize(newFaceSize);
497 if (newFace.size() > 2)
499 newFaces[newFacei++] = face(newFace);
504 <<
"Only " << newFace.size() <<
" in face " << facei
509 newFaces.setSize(newFacei);
511 MeshedSurface<face> fMesh(std::move(newPoints), std::move(newFaces));
516 boundaryPointRegions,
522 forAll(newPointNumbers, pointi)
524 if (newPointNumbers[pointi] == -1)
527 << pointi <<
" will be deleted and " << newPointNumbers[pointi]
528 <<
", so it will not be replaced. " 529 <<
"This will cause edges to be deleted." <<
endl;
537 Info<<
" Maximum number of chained collapses = " << maxChain <<
endl;
546 os <<
"Short Edge Filtering Information:" <<
nl 547 <<
" shortEdgeFilterFactor: " << shortEdgeFilterFactor_ <<
nl 548 <<
" edgeAttachedToBoundaryFactor: " << edgeAttachedToBoundaryFactor_
551 forAll(patchNames_, patchi)
553 os <<
" Patch " << patchNames_[patchi]
554 <<
", size " << patchSizes_[patchi] <<
endl;
557 os <<
" There are " << mapEdgesRegion_.size()
558 <<
" boundary edges." <<
endl;
560 os <<
" Mesh Info:" <<
nl 561 <<
" Points: " << ms_.nPoints() <<
nl 562 <<
" Faces: " << ms_.size() <<
nl 563 <<
" Edges: " << ms_.nEdges() <<
nl 564 <<
" Internal: " << ms_.nInternalEdges() <<
nl 565 <<
" External: " << ms_.nEdges() - ms_.nInternalEdges()
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
List of values generated by applying the access operation to each list item.
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
List< edge > edgeList
List of edge.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
constexpr char nl
The newline '\n' character (0x0a)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Foam::point toPoint3D(const point2D &) const
#define forAll(list, i)
Loop across all elements in list.
List< face > faceList
List of faces.
unsigned int count(const UList< bool > &bools, const bool val=true)
Count number of 'true' entries.
vectorField pointField
pointField is a vectorField.
const dimensionedScalar e
Elementary charge.
void setSize(const label n)
Alias for resize()
void calcDual(point2DField &dualPoints, faceList &dualFaces, wordList &patchNames, labelList &patchSizes, EdgeMap< label > &mapEdgesRegion, EdgeMap< label > &indirectPatchEdge) const
Calculates dual points (circumcentres of tets) and faces.
int debug
Static debugging option.
OBJstream os(runTime.globalPath()/outputName)
defineTypeNameAndDebug(combustionModel, 0)
Info<< "Finished reading KIVA file"<< endl;cellShapeList cellShapes(nPoints);labelList cellZoning(nPoints, -1);const cellModel &hex=cellModel::ref(cellModel::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].reset(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]
vector2DField point2DField
point2DField is a vector2DField.
#define WarningInFunction
Report a warning using Foam::Warning.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Conformal-Voronoi 2D automatic mesher with grid or read initial points and point position relaxation ...
List< label > labelList
A List of labels.
void writeInfo(Ostream &os)
List< bool > boolList
A List of bools.
forAllConstIters(mixture.phases(), phase)
static constexpr const zero Zero
Global zero (0)