diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
index cde15d06621414b1d4eb3c1604e9ea198d15a42b..f9648dc1e7749c6cd78f444c5a872d8aefbc6fe2 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
@@ -536,6 +536,11 @@ private:
const Delaunay::Finite_edges_iterator& eit
) const;
+ boolList dualFaceBoundaryPoints
+ (
+ const Delaunay::Finite_edges_iterator& eit
+ ) const;
+
//- Finds the maximum filterCount of the dual vertices
// (Delaunay cells) that form the dual face produced by the
// supplied edge
@@ -842,6 +847,13 @@ private:
const Delaunay::Finite_facets_iterator& fit
) const;
+ //- Merge adjacent edges that are not attached to other faces
+ label mergeNearlyParallelEdges
+ (
+ const pointField& pts,
+ const scalar maxCosAngle
+ );
+
//- Merge vertices that are very close together
void mergeCloseDualVertices
(
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C
index 7cc8957de09059bfcccf443f1ca48f886cdba376..85b4ba7b27af55cdb62988913931f71464cbc992 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C
@@ -212,6 +212,20 @@ void Foam::conformalVoronoiMesh::calcDualMesh
Info<< nl << "Collapsing unnecessary faces" << endl;
collapseFaces(points, boundaryPts, deferredCollapseFaces);
+
+ const scalar maxCosAngle
+ = cos(degToRad(cvMeshControls().edgeMergeAngle()));
+
+ Info<< nl << "Merging adjacent edges which have an angle "
+ << "of greater than "
+ << cvMeshControls().edgeMergeAngle() << ": " << endl;
+
+ label nRemovedEdges =
+ mergeNearlyParallelEdges(points, maxCosAngle);
+
+ reduce(nRemovedEdges, sumOp<label>());
+
+ Info<< " Merged " << nRemovedEdges << " edges" << endl;
}
labelHashSet wrongFaces = checkPolyMeshQuality(points);
@@ -652,6 +666,113 @@ Foam::label Foam::conformalVoronoiMesh::mergeCloseDualVertices
}
+Foam::label Foam::conformalVoronoiMesh::mergeNearlyParallelEdges
+(
+ const pointField& pts,
+ const scalar maxCosAngle
+)
+{
+ List<HashSet<label> > pointFaceCount(number_of_cells());
+ labelList pointNeighbour(number_of_cells(), -1);
+
+ Map<label> dualPtIndexMap;
+
+ for
+ (
+ Delaunay::Finite_edges_iterator eit = finite_edges_begin();
+ eit != finite_edges_end();
+ ++eit
+ )
+ {
+ Cell_handle c = eit->first;
+ Vertex_handle vA = c->vertex(eit->second);
+ Vertex_handle vB = c->vertex(eit->third);
+
+ if (isBoundaryDualFace(eit))
+ {
+ const face f = buildDualFace(eit);
+
+ forAll(f, pI)
+ {
+ const label pIndex = f[pI];
+
+ const label prevPointI = f.prevLabel(pI);
+ const label nextPointI = f.nextLabel(pI);
+
+ pointFaceCount[pIndex].insert(prevPointI);
+ pointFaceCount[pIndex].insert(nextPointI);
+ pointNeighbour[pIndex] = nextPointI;
+ }
+ }
+ else if
+ (
+ vA->internalOrBoundaryPoint()
+ || vB->internalOrBoundaryPoint()
+ )
+ {
+ const face f = buildDualFace(eit);
+ const boolList faceBoundaryPoints = dualFaceBoundaryPoints(eit);
+
+ forAll(f, pI)
+ {
+ const label pIndex = f[pI];
+
+ const label prevPointI = f.prevLabel(pI);
+ const label nextPointI = f.nextLabel(pI);
+
+ pointFaceCount[pIndex].insert(prevPointI);
+ pointFaceCount[pIndex].insert(nextPointI);
+
+ if (faceBoundaryPoints[pI] == false)
+ {
+ pointNeighbour[pIndex] = nextPointI;
+ }
+ else
+ {
+ if (faceBoundaryPoints[prevPointI] == true)
+ {
+ pointNeighbour[pIndex] = prevPointI;
+ }
+ else if (faceBoundaryPoints[nextPointI] == true)
+ {
+ pointNeighbour[pIndex] = nextPointI;
+ }
+ else
+ {
+ pointNeighbour[pIndex] = pIndex;
+ }
+ }
+ }
+ }
+ }
+
+ forAll(pointFaceCount, pI)
+ {
+ if (pointFaceCount[pI].size() == 2)
+ {
+ List<vector> edges(2, vector(0, 0, 0));
+
+ label count = 0;
+ forAllConstIter(HashSet<label>, pointFaceCount[pI], iter)
+ {
+ edges[count] = pts[pI] - pts[iter.key()];
+ edges[count] /= mag(edges[count]) + VSMALL;
+ count++;
+ }
+
+ if (mag(edges[0] & edges[1]) > maxCosAngle)
+ {
+ dualPtIndexMap.insert(pI, pointNeighbour[pI]);
+ }
+ }
+ }
+
+ reindexDualVertices(dualPtIndexMap);
+
+ return dualPtIndexMap.size();
+}
+
+
void Foam::conformalVoronoiMesh::smoothSurface
(
pointField& pts,
@@ -696,7 +817,6 @@ void Foam::conformalVoronoiMesh::smoothSurface
} while (nCollapsedFaces > 0);
// Force all points of boundary faces to be on the surface
-
for
(
Delaunay::Finite_cells_iterator cit = finite_cells_begin();
@@ -917,6 +1037,7 @@ void Foam::conformalVoronoiMesh::collapseFaces
mergeCloseDualVertices(pts, boundaryPts);
+
if (nCollapsedFaces > 0)
{
Info<< " Collapsed " << nCollapsedFaces << " faces" << endl;
@@ -932,6 +1053,7 @@ void Foam::conformalVoronoiMesh::collapseFaces
}
} while (nCollapsedFaces > 0);
+
}
@@ -1064,6 +1186,25 @@ Foam::conformalVoronoiMesh::collapseFace
bool allowEarlyCollapseToPoint = true;
+
+// // Quick exit
+// label smallEdges = 0;
+// const edgeList& fEdges = f.edges();
+// forAll(fEdges, eI)
+// {
+// const edge& e = fEdges[eI];
+//
+// if (e.mag(pts) < 0.2*targetFaceSize)
+// {
+// smallEdges++;
+// }
+// }
+// if (smallEdges == 0)
+// {
+// return fcmNone;
+// }
+
+
// if (maxFC > cvMeshControls().filterCountSkipThreshold() - 3)
// {
// limitToQuadsOrTris = true;
@@ -1071,6 +1212,7 @@ Foam::conformalVoronoiMesh::collapseFace
// allowEarlyCollapseToPoint = false;
// }
+
collapseSizeLimitCoeff *= pow
(
cvMeshControls().filterErrorReductionCoeff(),
@@ -1158,6 +1300,91 @@ Foam::conformalVoronoiMesh::collapseFace
}
}
+
+ scalar maxDist = 0;
+ scalar minDist = GREAT;
+//
+// if (f.size() <= 3)
+// {
+// const edgeList& fEdges = f.edges();
+//
+// forAll(fEdges, eI)
+// {
+// const edge& e = fEdges[eI];
+// const scalar d = e.mag(pts);
+//
+// if (d > maxDist)
+// {
+// maxDist = d;
+// collapseAxis = e.vec(pts);
+// }
+// else if (d < minDist && d != 0)
+// {
+// minDist = d;
+// }
+// }
+// }
+// else
+// {
+// forAll(f, pI)
+// {
+// for (label i = pI + 1; i < f.size(); ++i)
+// {
+// if
+// (
+// f[i] != f.nextLabel(pI)
+// && f[i] != f.prevLabel(pI)
+// )
+// {
+// scalar d = mag(pts[f[pI]] - pts[f[i]]);
+//
+// if (d > maxDist)
+// {
+// maxDist = d;
+// collapseAxis = pts[f[pI]] - pts[f[i]];
+// }
+// else if (d < minDist && d != 0)
+// {
+// minDist = d;
+// }
+// }
+// }
+// }
+// }
+//
+// const edgeList& fEdges = f.edges();
+//
+// scalar perimeter = 0;
+//
+// forAll(fEdges, eI)
+// {
+// const edge& e = fEdges[eI];
+// const scalar d = e.mag(pts);
+//
+// perimeter += d;
+//
+//// collapseAxis += e.vec(pts);
+//
+// if (d > maxDist)
+// {
+// collapseAxis = e.vec(pts);
+// maxDist = d;
+// }
+// else if (d < minDist && d != 0)
+// {
+// minDist = d;
+// }
+// }
+//
+// collapseAxis /= mag(collapseAxis);
+//
+//// Info<< f.size() << " " << minDist << " " << maxDist << " | "
+//// << collapseAxis << endl;
+//
+// aspectRatio = maxDist/minDist;
+//
+// aspectRatio = min(aspectRatio, sqr(perimeter)/(16.0*fA));
+
if (magSqr(collapseAxis) < VSMALL)
{
WarningIn
@@ -1170,9 +1397,9 @@ Foam::conformalVoronoiMesh::collapseFace
// Output face and collapse axis for visualisation
Pout<< "# Aspect ratio = " << aspectRatio << nl
- << "# inertia = " << J << nl
- << "# determinant = " << detJ << nl
- << "# eigenvalues = " << eigenValues(J) << nl
+// << "# inertia = " << J << nl
+// << "# determinant = " << detJ << nl
+// << "# eigenvalues = " << eigenValues(J) << nl
<< "# collapseAxis = " << collapseAxis << nl
<< "# facePts = " << facePts << nl
<< endl;
@@ -1280,8 +1507,6 @@ Foam::conformalVoronoiMesh::collapseFace
{
scalar guardFraction = cvMeshControls().edgeCollapseGuardFraction();
- cvMeshControls().maxCollapseFaceToPointSideLengthCoeff();
-
if
(
allowEarlyCollapseToPoint
@@ -1337,6 +1562,8 @@ Foam::conformalVoronoiMesh::collapseFace
Foam::point collapseToPt =
collapseAxis*(sum(dNeg)/dNeg.size() - dShift) + fC;
+// DynamicList<label> faceBoundaryPts(f.size());
+
forAll(facePtsNeg, fPtI)
{
if (boundaryPts[facePtsNeg[fPtI]] == true)
@@ -1351,9 +1578,32 @@ Foam::conformalVoronoiMesh::collapseFace
collapseToPt = pts[collapseToPtI];
break;
+
+// faceBoundaryPts.append(facePtsNeg[fPtI]);
}
}
+// if (!faceBoundaryPts.empty())
+// {
+// if (faceBoundaryPts.size() == 2)
+// {
+// collapseToPtI = faceBoundaryPts[0];
+//
+// collapseToPt =
+// 0.5*(pts[faceBoundaryPts[0]] + pts[faceBoundaryPts[1]]);
+// }
+// else if (faceBoundaryPts.size() < f.size())
+// {
+// face bFace(faceBoundaryPts);
+//
+// collapseToPtI = faceBoundaryPts.first();
+//
+// collapseToPt = bFace.centre(pts);
+// }
+// }
+//
+// faceBoundaryPts.clear();
+
// ...otherwise arbitrarily choosing the most distant
// point as the index to collapse to.
@@ -1384,9 +1634,30 @@ Foam::conformalVoronoiMesh::collapseFace
collapseToPt = pts[collapseToPtI];
break;
+
+// faceBoundaryPts.append(facePtsNeg[fPtI]);
}
}
+// if (!faceBoundaryPts.empty())
+// {
+// if (faceBoundaryPts.size() == 2)
+// {
+// collapseToPtI = faceBoundaryPts[0];
+//
+// collapseToPt =
+// 0.5*(pts[faceBoundaryPts[0]] + pts[faceBoundaryPts[1]]);
+// }
+// else if (faceBoundaryPts.size() < f.size())
+// {
+// face bFace(faceBoundaryPts);
+//
+// collapseToPtI = faceBoundaryPts.first();
+//
+// collapseToPt = bFace.centre(pts);
+// }
+// }
+
// ...otherwise arbitrarily choosing the most distant
// point as the index to collapse to.
@@ -1406,6 +1677,8 @@ Foam::conformalVoronoiMesh::collapseFace
Foam::point collapseToPt = fC;
+ DynamicList<label> faceBoundaryPts(f.size());
+
forAll(facePts, fPtI)
{
if (boundaryPts[facePts[fPtI]] == true)
@@ -1415,11 +1688,32 @@ Foam::conformalVoronoiMesh::collapseFace
// use the first boundary point encountered if
// there are multiple boundary points.
- collapseToPtI = facePts[fPtI];
+// collapseToPtI = facePts[fPtI];
+//
+// collapseToPt = pts[collapseToPtI];
+//
+// break;
- collapseToPt = pts[collapseToPtI];
+ faceBoundaryPts.append(facePts[fPtI]);
+ }
+ }
- break;
+ if (!faceBoundaryPts.empty())
+ {
+ if (faceBoundaryPts.size() == 2)
+ {
+ collapseToPtI = faceBoundaryPts[0];
+
+ collapseToPt =
+ 0.5*(pts[faceBoundaryPts[0]] + pts[faceBoundaryPts[1]]);
+ }
+ else if (faceBoundaryPts.size() < f.size())
+ {
+ face bFace(faceBoundaryPts);
+
+ collapseToPtI = faceBoundaryPts.first();
+
+ collapseToPt = bFace.centre(pts);
}
}
@@ -1858,7 +2152,10 @@ Foam::labelHashSet Foam::conformalVoronoiMesh::findOffsetPatchFaces
}
}
- offsetBoundaryCells.write();
+ if (cvMeshControls().objOutput())
+ {
+ offsetBoundaryCells.write();
+ }
return offsetBoundaryCells;
}
@@ -2092,13 +2389,7 @@ void Foam::conformalVoronoiMesh::indexDualVertices
pts[dualVertI] = cit->dual();
- if
- (
- !cit->vertex(0)->internalOrBoundaryPoint()
- || !cit->vertex(1)->internalOrBoundaryPoint()
- || !cit->vertex(2)->internalOrBoundaryPoint()
- || !cit->vertex(3)->internalOrBoundaryPoint()
- )
+ if (cit->boundaryDualVertex())
{
// This is a boundary dual vertex
boundaryPts[dualVertI] = true;
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H
index 1e24b7d8ef29cf030425647cfc70fc03f50bf26a..43d8c077fc785724ad83b37ac8e48a7a76f77ba9 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H
@@ -366,6 +366,49 @@ inline bool Foam::conformalVoronoiMesh::isBoundaryDualFace
}
+inline Foam::List<bool> Foam::conformalVoronoiMesh::dualFaceBoundaryPoints
+(
+ const Delaunay::Finite_edges_iterator& eit
+) const
+{
+ Cell_circulator ccStart = incident_cells(*eit);
+ Cell_circulator cc1 = ccStart;
+ Cell_circulator cc2 = cc1;
+
+ // Advance the second circulator so that it always stays on the next
+ // cell around the edge;
+ cc2++;
+
+ DynamicList<bool> tmpFaceBoundaryPoints;
+
+ do
+ {
+ label cc1I = cc1->cellIndex();
+
+ label cc2I = cc2->cellIndex();
+
+ if (cc1I != cc2I)
+ {
+ if (cc1->boundaryDualVertex())
+ {
+ tmpFaceBoundaryPoints.append(true);
+ }
+ else
+ {
+ tmpFaceBoundaryPoints.append(false);
+ }
+ }
+
+ cc1++;
+
+ cc2++;
+
+ } while (cc1 != ccStart);
+
+ return tmpFaceBoundaryPoints;
+}
+
+
inline Foam::List<Foam::label> Foam::conformalVoronoiMesh::processorsAttached
(
const Delaunay::Finite_facets_iterator& fit