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.first()].push_uniq(regi);
50 boundaryPointRegions[
e.second()].push_uniq(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();
74 if (surfMesh.isInternalEdge(edgeI))
79 const edge&
e = edges[edgeI];
83 const DynamicList<label>& startPtRegions =
84 boundaryPtRegions[surfPtToBoundaryPt[meshPoints[
e.first()]]];
86 const DynamicList<label>& endPtRegions =
87 boundaryPtRegions[surfPtToBoundaryPt[meshPoints[
e.second()]]];
89 if (startPtRegions.size() > 1 && endPtRegions.size() > 1)
91 region = startPtRegions[0];
94 <<
"Both points in edge are in different regions." 95 <<
" Assigning edge to region " << region
98 else if (startPtRegions.size() > 1 || endPtRegions.size() > 1)
102 startPtRegions.size() > 1
109 startPtRegions[0] == endPtRegions[0]
110 && startPtRegions[0] != -1
113 region = startPtRegions[0];
118 newMapEdgesRegion.insert(
e, region);
119 patchSizes[region]++;
123 mapEdgesRegion.transfer(newMapEdgesRegion);
129 Foam::shortEdgeFilter2D::shortEdgeFilter2D
132 const dictionary&
dict 135 shortEdgeFilterFactor_(
dict.
get<scalar>(
"shortEdgeFilterFactor")),
136 edgeAttachedToBoundaryFactor_
138 dict.getOrDefault<scalar>(
"edgeAttachedToBoundaryFactor", 2.0)
166 OFstream str(
"indirectPatchEdges.obj");
169 Info<<
"Writing indirectPatchEdges to " << str.name() <<
endl;
173 const edge&
e = iter.key();
187 ms_ = MeshedSurface<face>(std::move(
points), std::move(faces));
189 Info<<
"Meshed surface stats before edge filtering :" <<
endl;
190 ms_.writeStats(
Info);
196 ms_.write(
"MeshedSurface_preFilter.obj");
207 const edgeList& edges = ms_.edges();
208 const faceList& faces = ms_.surfFaces();
209 const labelList& meshPoints = ms_.meshPoints();
210 const labelList& boundaryPoints = ms_.boundaryPoints();
213 label nPointsToRemove = 0;
215 labelList pointsToRemove(ms_.points().size(), -1);
218 labelList newFaceVertexCount(faces.size(), -1);
221 newFaceVertexCount[facei] = faces[facei].size();
226 List<DynamicList<label>> boundaryPointRegions
231 assignBoundaryPointRegions(boundaryPointRegions);
235 Info<<
" Marking edges attached to boundaries." <<
endl;
236 boolList edgeAttachedToBoundary(edges.size(),
false);
239 const edge&
e = edges[edgeI];
240 const label startVertex =
e.start();
241 const label endVertex =
e.end();
243 forAll(boundaryPoints, bPoint)
247 boundaryPoints[bPoint] == startVertex
248 || boundaryPoints[bPoint] == endVertex
251 edgeAttachedToBoundary[edgeI] =
true;
258 const edge&
e = edges[edgeI];
261 const label startVertex =
e.start();
262 const label endVertex =
e.end();
267 points[meshPoints[startVertex]]
268 -
points[meshPoints[endVertex]]
271 if (edgeAttachedToBoundary[edgeI])
273 edgeLength *= edgeAttachedToBoundaryFactor_;
276 scalar shortEdgeFilterValue = 0.0;
278 const labelList& psEdges = ms_.pointEdges()[startVertex];
279 const labelList& peEdges = ms_.pointEdges()[endVertex];
283 const edge& psE = edges[psEdges[psEdgeI]];
284 if (edgeI != psEdges[psEdgeI])
286 shortEdgeFilterValue +=
289 points[meshPoints[psE.start()]]
290 -
points[meshPoints[psE.end()]]
297 const edge& peE = edges[peEdges[peEdgeI]];
298 if (edgeI != peEdges[peEdgeI])
300 shortEdgeFilterValue +=
303 points[meshPoints[peE.start()]]
304 -
points[meshPoints[peE.end()]]
309 shortEdgeFilterValue *=
310 shortEdgeFilterFactor_
311 /(psEdges.size() + peEdges.size() - 2);
313 edge lookupInPatchEdge
315 meshPoints[startVertex],
316 meshPoints[endVertex]
321 edgeLength < shortEdgeFilterValue
322 || indirectPatchEdge_.found(lookupInPatchEdge)
325 bool flagDegenerateFace =
false;
330 const face&
f = ms_.localFaces()[
pFaces[pFacei]];
334 if (
f[fp] == endVertex)
337 if (newFaceVertexCount[
pFaces[pFacei]] < 4)
339 flagDegenerateFace =
true;
343 newFaceVertexCount[
pFaces[pFacei]]--;
350 if (newFaceVertexCount[
pFaces[pFacei]] < 3)
352 flagDegenerateFace =
true;
361 pointsToRemove[meshPoints[startVertex]] == -1
362 && pointsToRemove[meshPoints[endVertex]] == -1
363 && !flagDegenerateFace
366 const DynamicList<label>& startVertexRegions =
367 boundaryPointRegions[meshPoints[startVertex]];
368 const DynamicList<label>& endVertexRegions =
369 boundaryPointRegions[meshPoints[endVertex]];
371 if (startVertexRegions.size() && endVertexRegions.size())
373 if (startVertexRegions.size() > 1)
375 pointsToRemove[meshPoints[endVertex]] =
376 meshPoints[startVertex];
380 pointsToRemove[meshPoints[startVertex]] =
381 meshPoints[endVertex];
384 else if (startVertexRegions.size())
386 pointsToRemove[meshPoints[endVertex]] =
387 meshPoints[startVertex];
391 pointsToRemove[meshPoints[startVertex]] =
392 meshPoints[endVertex];
400 label totalNewPoints =
points.size() - nPointsToRemove;
404 label numberRemoved = 0;
407 labelList newPtToOldPt(totalNewPoints, -1);
412 if (pointsToRemove[pointi] == -1)
414 newPoints[pointi - numberRemoved] =
points[pointi];
415 newPointNumbers[pointi] = pointi - numberRemoved;
416 newPtToOldPt[pointi - numberRemoved] = pointi;
429 label newFaceSize = 0;
434 const face&
f = faces[facei];
437 newFace.setSize(
f.size());
442 label pointi =
f[fp];
444 if (pointsToRemove[pointi] == -1)
446 newFace[newFaceSize++] = newPointNumbers[pointi];
450 label
newPointi = pointsToRemove[pointi];
460 label totalChain = 0;
461 for (label nChain = 0; nChain <= totalChain; ++nChain)
463 if (newPointNumbers[pChain] != -1)
465 newFace[newFaceSize++] = newPointNumbers[pChain];
466 newPointNumbers[pointi] = newPointNumbers[pChain];
467 maxChain =
max(totalChain, maxChain);
472 <<
"Point " << pChain
473 <<
" marked for deletion as well as point " 475 <<
" Incrementing maxChain by 1 from " 476 << totalChain <<
" to " << totalChain + 1
480 pChain = pointsToRemove[pChain];
487 newPointNumbers[pointi] = newPointNumbers[
newPointi];
493 newFace.setSize(newFaceSize);
495 if (newFace.size() > 2)
497 newFaces[newFacei++] = face(newFace);
502 <<
"Only " << newFace.size() <<
" in face " << facei
507 newFaces.setSize(newFacei);
509 MeshedSurface<face> fMesh(std::move(newPoints), std::move(newFaces));
514 boundaryPointRegions,
520 forAll(newPointNumbers, pointi)
522 if (newPointNumbers[pointi] == -1)
525 << pointi <<
" will be deleted and " << newPointNumbers[pointi]
526 <<
", so it will not be replaced. " 527 <<
"This will cause edges to be deleted." <<
endl;
535 Info<<
" Maximum number of chained collapses = " << maxChain <<
endl;
544 os <<
"Short Edge Filtering Information:" <<
nl 545 <<
" shortEdgeFilterFactor: " << shortEdgeFilterFactor_ <<
nl 546 <<
" edgeAttachedToBoundaryFactor: " << edgeAttachedToBoundaryFactor_
549 forAll(patchNames_, patchi)
551 os <<
" Patch " << patchNames_[patchi]
552 <<
", size " << patchSizes_[patchi] <<
endl;
555 os <<
" There are " << mapEdgesRegion_.size()
556 <<
" boundary edges." <<
endl;
558 os <<
" Mesh Info:" <<
nl 559 <<
" Points: " << ms_.nPoints() <<
nl 560 <<
" Faces: " << ms_.size() <<
nl 561 <<
" Edges: " << ms_.nEdges() <<
nl 562 <<
" Internal: " << ms_.nInternalEdges() <<
nl 563 <<
" 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
void resize_nocopy(const label len)
Adjust allocated size of list without necessarily.
#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 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)