zhyang-dev/openfoam
1
0
Fork 0
mirror of https://develop.openfoam.com/Development/openfoam.git synced 2025-11-28 03:28:01 +00:00
Code Issues Actions Packages Projects Releases Wiki Activity

Adding filter limit scalars to the Delaunay cells. After initial face

Browse Source 
check, identifying faces that have bad quality, identifying the points
that support the faces, and 'limiting' them.

Wrong faces are also those that are attached to a cell with 1-3 faces.

At each iteration of filtering, the Delaunay cells are all reindexed.
Now a separate function.

Currently 'limiting' means: do not collapse any face that has any
vertex that has been limited at all and do not move a vertex to the
boundary during the surface smooth if it has been limited.  This is
the most agressive form of limiting and leaves faces in the mesh that
might be safely removed with an less aggressive filter.

Limiting of all points connected to a cell which has a limited
point could be used (comment out code).

Removed some redundant code from indexedCell.
This commit is contained in:
graham
2010-01-13 16:12:05 +00:00
parent 6d4cb53c56
commit cf0c71789a
6 changed files with 321 additions and 110 deletions
Download Patch File Download Diff File Expand all files Collapse all files

40
src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
View File

@ -1030,6 +1030,46 @@ Foam::face Foam::conformalVoronoiMesh::buildDualFace
}
Foam::scalar Foam::conformalVoronoiMesh::minFilterLimit
(
const Triangulation::Finite_edges_iterator& eit
) const
{
Cell_circulator ccStart = incident_cells(*eit);
Cell_circulator cc = ccStart;
scalar minFilterLimit = GREAT;
do
{
if (cc->cellIndex() < 0)
{
Cell_handle c = eit->first;
Vertex_handle vA = c->vertex(eit->second);
Vertex_handle vB = c->vertex(eit->third);
FatalErrorIn("Foam::conformalVoronoiMesh::buildDualFace")
<< "Dual face uses circumcenter defined by a "
<< "Delaunay tetrahedron with no internal "
<< "or boundary points. Defining Delaunay edge ends: "
<< topoint(vA->point()) << " "
<< topoint(vB->point()) << nl
<< exit(FatalError);
}
if (cc->filterLimit() < minFilterLimit)
{
minFilterLimit = cc->filterLimit();
}
cc++;
} while (cc != ccStart);
return minFilterLimit;
}
bool Foam::conformalVoronoiMesh::ownerAndNeighbour
(
Vertex_handle vA,

12
src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
View File

@ -357,6 +357,14 @@ private:
const Triangulation::Finite_edges_iterator& eit
) const;
//- Finds the minimum filterLimit of the dual vertices
// (Delaunay cells) that form the dual face produced by the
// supplied edge
scalar minFilterLimit
(
const Triangulation::Finite_edges_iterator& eit
) const;
//- Determines the owner and neighbour labels for dual cells
// corresponding to the dual face formed by the supplied
// Delaunay vertices. If the dual face is a boundary face
@ -552,6 +560,10 @@ private:
// elements damage the mesh quality, allowing backtracking.
label checkPolyMeshQuality(const pointField& pts) const;
//- Index all of the the Delaunay cells and calculate their
//- dual points
void indexDualVertices(pointField& pts);
//- Re-index all of the the Delaunay cells
void reindexDualVertices(const Map<label>& dualPtIndexMap);

271
src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C
View File

@ -76,53 +76,28 @@ void Foam::conformalVoronoiMesh::calcDualMesh
}
}
// Indexing Delaunay cells, which are the dual vertices
indexDualVertices(points);
label dualVertI = 0;
points.setSize(number_of_cells());
for
(
Triangulation::Finite_cells_iterator cit = finite_cells_begin();
cit != finite_cells_end();
++cit
)
{
if
(
cit->vertex(0)->internalOrBoundaryPoint()
|| cit->vertex(1)->internalOrBoundaryPoint()
|| cit->vertex(2)->internalOrBoundaryPoint()
|| cit->vertex(3)->internalOrBoundaryPoint()
)
{
cit->cellIndex() = dualVertI;
points[dualVertI] = topoint(dual(cit));
dualVertI++;
}
else
{
cit->cellIndex() = -1;
}
// No-risk face filtering to get rid of zero area faces and
// establish if the mesh can be produced at all to the
// specified criteria
Info<< nl << " Merging close points" << endl;
// There is no guarantee that a merge of close points is
// no-risk, but it seems to work using 1e-4 as the mergeClosenessCoeff
mergeCloseDualVertices(points);
}
points.setSize(dualVertI);
// No-risk face filtering to get rid of zero area faces and
// establish if the mesh can be produced at all to the
// specified criteria
Info<< nl << " Merging close points" << endl;
mergeCloseDualVertices(points);
label nInitialBadQualityFaces = checkPolyMeshQuality(points);
Info<< "Initial check before face collapse, found "
<< nInitialBadQualityFaces << " bad quality faces"
<< endl;
HashSet<labelPair, labelPair::Hash<> > deferredCollapseFaces;
if (nInitialBadQualityFaces > 0)
{
Info<< "A mesh could not be produced to satisfy the specified quality "
@ -132,24 +107,41 @@ void Foam::conformalVoronoiMesh::calcDualMesh
<< "be applied in areas where problems are occurring."
<< endl;
}
HashSet<labelPair, labelPair::Hash<> > deferredCollapseFaces;
else
{
// Risky and undo-able face filtering to reduce the face count
// as much as possible staying within the specified criteria
label nBadQualityFaces = 0;
Info<< nl << " Smoothing surface" << endl;
do
{
// Reindexing the Delaunay cells and regenerating the
// points resets the mesh to the starting condition.
smoothSurface(points);
indexDualVertices(points);
Info<< nl << " Collapsing unnecessary faces" << endl;
{
Info<< nl << " Merging close points" << endl;
collapseFaces(points, deferredCollapseFaces);
mergeCloseDualVertices(points);
}
label nBadQualityFaces = checkPolyMeshQuality(points);
{
// Risky and undo-able face filtering to reduce the face count
// as much as possible, staying within the specified criteria
Info<< "Found " << nBadQualityFaces << " bad quality faces" << endl;
Info<< nl << " Smoothing surface" << endl;
smoothSurface(points);
Info<< nl << " Collapsing unnecessary faces" << endl;
collapseFaces(points, deferredCollapseFaces);
}
nBadQualityFaces = checkPolyMeshQuality(points);
Info<< "Found " << nBadQualityFaces << " bad quality faces" << endl;
} while (nBadQualityFaces > 0);
}
// Final dual face and owner neighbour construction
@ -485,6 +477,12 @@ void Foam::conformalVoronoiMesh::smoothSurface(pointField& pts)
{
label ptI = cit->cellIndex();
if (cit->filterLimit() < 1.0)
{
// This vertex has been limited, skip
continue;
}
// Only cells with indices > -1 are valid
if (ptI > -1)
{
@ -580,6 +578,14 @@ Foam::label Foam::conformalVoronoiMesh::smoothSurfaceDualFaces
continue;
}
scalar minFL = minFilterLimit(eit);
if (minFL < 1.0)
{
// A vertex on this face has been limited, skip
continue;
}
pointIndexHit surfHit;
label hitSurface;
@ -749,6 +755,14 @@ Foam::label Foam::conformalVoronoiMesh::collapseFaces
continue;
}
scalar minFL = minFilterLimit(eit);
if (minFL < 1.0)
{
// A vertex on this face has been limited, skip
continue;
}
scalar targetFaceSize = averageAnyCellSize(vA, vB);
faceCollapseMode mode = collapseFace
@ -1256,31 +1270,138 @@ Foam::label Foam::conformalVoronoiMesh::checkPolyMeshQuality
motionSmoother::checkMesh
(
false,
false, // report
pMesh,
cvMeshControls().cvMeshDict().subDict("meshQualityControls"),
wrongFaces
);
label nInvalidPolyhedra = 0;
const cellList& cells = pMesh.cells();
forAll(cells, cI)
{
if (cells[cI].size() < 4)
if (cells[cI].size() < 4 && cells[cI].size() > 0)
{
Info<< "cell " << cI
<< " has " << cells[cI].size() << " faces."
<< endl;
// Info<< "cell " << cI << " " << cells[cI]
// << " has " << cells[cI].size() << " faces."
// << endl;
nInvalidPolyhedra++;
forAll(cells[cI], cFI)
{
wrongFaces.insert(cells[cI][cFI]);
}
}
}
Info<< "Cells with fewer than 4 faces : "
<< nInvalidPolyhedra << endl;
PackedBoolList ptToBeLimited(pts.size(), false);
forAllConstIter(labelHashSet, wrongFaces, iter)
{
const face f = pMesh.faces()[iter.key()];
forAll(f, fPtI)
{
ptToBeLimited[f[fPtI]] = true;
}
}
// // Limit connected cells
// labelHashSet limitCells(pMesh.nCells()/100);
// const labelListList& ptCells = pMesh.pointCells();
// forAllConstIter(labelHashSet, wrongFaces, iter)
// {
// label faceI = iter.key();
// const face f = pMesh.faces()[iter.key()];
// Info<< faceI << " " << pMesh.faces()[faceI] << endl;
// forAll(f, fPtI)
// {
// label ptI = f[fPtI];
// const labelList& pC = ptCells[ptI];
// forAll(pC, pCI)
// {
// limitCells.insert(pC[pCI]);
// }
// }
// }
// const labelListList& cellPts = pMesh.cellPoints();
// forAllConstIter(labelHashSet, limitCells, iter)
// {
// label cellI = iter.key();
// const labelList& cP = cellPts[cellI];
// forAll(cP, cPI)
// {
// ptToBeLimited[cP[cPI]] = true;
// }
// }
// Limit Delaunay cell filter values
for
(
Triangulation::Finite_cells_iterator cit = finite_cells_begin();
cit != finite_cells_end();
++cit
)
{
label cI = cit->cellIndex();
if (cI >= 0)
{
if (ptToBeLimited[cI] == true)
{
cit->filterLimit() *= 0.75;
}
}
}
// Determine the minimum minFilterLimit
scalar minMinFilterLimit = GREAT;
for
(
Triangulation::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
(
vA->internalOrBoundaryPoint()
|| vB->internalOrBoundaryPoint()
)
{
scalar minFL = minFilterLimit(eit);
if (minFL < minMinFilterLimit)
{
minMinFilterLimit = minFL;
}
}
}
Info<< "minMinFilterLimit " << minMinFilterLimit << endl;
return wrongFaces.size();
// For parallel running:
@ -1292,6 +1413,46 @@ Foam::label Foam::conformalVoronoiMesh::checkPolyMeshQuality
}
void Foam::conformalVoronoiMesh::indexDualVertices
(
pointField& pts
)
{
// Indexing Delaunay cells, which are the dual vertices
label dualVertI = 0;
pts.setSize(number_of_cells());
for
(
Triangulation::Finite_cells_iterator cit = finite_cells_begin();
cit != finite_cells_end();
++cit
)
{
if
(
cit->vertex(0)->internalOrBoundaryPoint()
|| cit->vertex(1)->internalOrBoundaryPoint()
|| cit->vertex(2)->internalOrBoundaryPoint()
|| cit->vertex(3)->internalOrBoundaryPoint()
)
{
cit->cellIndex() = dualVertI;
pts[dualVertI] = topoint(dual(cit));
dualVertI++;
}
else
{
cit->cellIndex() = -1;
}
}
pts.setSize(dualVertI);
}
void Foam::conformalVoronoiMesh::reindexDualVertices
(
const Map<label>& dualPtIndexMap

60
src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell.H
View File

@ -36,6 +36,7 @@ Description
#include <CGAL/Triangulation_3.h>
#include "indexedVertex.H"
#include "scalar.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -43,7 +44,7 @@ namespace CGAL
{
/*---------------------------------------------------------------------------*\
Class indexedCell Declaration
Class indexedCell Declaration
\*---------------------------------------------------------------------------*/
template<class Gt, class Cb=CGAL::Triangulation_cell_base_3<Gt> >
@ -54,20 +55,15 @@ class indexedCell
// Private data
//- The index for this tetrahedral cell
// -1: tetrahedral and changed and associated data needs updating
// >=0: index of cells, set by external update algorithm
int index_;
//- The applicable fraction of the face filtering mechanism
// controls
Foam::scalar filterLimit_;
public:
enum cellTypes
{
ctUnchanged = 0,
ctChanged = -1,
ctSaveChanged = -2
};
typedef typename Cb::Vertex_handle Vertex_handle;
typedef typename Cb::Cell_handle Cell_handle;
@ -82,18 +78,22 @@ public:
indexedCell()
:
Cb(),
index_(ctChanged)
index_(-1),
filterLimit_(1.0)
{}
indexedCell
(
Vertex_handle v0, Vertex_handle v1, Vertex_handle v2, Vertex_handle v3
)
:
Cb(v0, v1, v2, v3),
index_(ctChanged)
index_(-1),
filterLimit_(1.0)
{}
indexedCell
(
Vertex_handle v0,
@ -107,39 +107,33 @@ public:
)
:
Cb(v0, v1, v2, v3, n0, n1, n2, n3),
index_(ctChanged)
index_(-1),
filterLimit_(1.0)
{}
void set_vertex(int i, Vertex_handle v)
{
index_ = ctChanged;
Cb::set_vertex(i, v);
}
void set_vertices()
{
index_ = ctChanged;
Cb::set_vertices();
}
void set_vertices
(
Vertex_handle v0, Vertex_handle v1, Vertex_handle v2, Vertex_handle v3
)
{
index_ = ctChanged;
Cb::set_vertices(v0, v1, v2, v3);
}
int& cellIndex()
{
return index_;
}
int cellIndex() const
{
return index_;
}
inline Foam::scalar& filterLimit()
{
return filterLimit_;
}
inline Foam::scalar filterLimit() const
{
return filterLimit_;
}
};

28
src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedVertex.H
View File

@ -35,7 +35,6 @@ Description
#define indexedVertex_H
#include <CGAL/Triangulation_3.h>
#include "List.H"
#include "tensor.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -44,7 +43,7 @@ namespace CGAL
{
/*---------------------------------------------------------------------------*\
Class indexedVertex Declaration
Class indexedVertex Declaration
\*---------------------------------------------------------------------------*/
template<class Gt, class Vb=CGAL::Triangulation_vertex_base_3<Gt> >
@ -91,6 +90,7 @@ public:
typedef indexedVertex<Gt,Vb2> Other;
};
indexedVertex()
:
Vb(),
@ -100,6 +100,7 @@ public:
targetCellSize_(0.0)
{}
indexedVertex(const Point& p)
:
Vb(p),
@ -109,6 +110,7 @@ public:
targetCellSize_(0.0)
{}
indexedVertex(const Point& p, Cell_handle f)
:
Vb(f, p),
@ -118,6 +120,7 @@ public:
targetCellSize_(0.0)
{}
indexedVertex(Cell_handle f)
:
Vb(f),
@ -133,70 +136,83 @@ public:
return index_;
}
int index() const
{
return index_;
}
int& type()
{
return type_;
}
int type() const
{
return type_;
}
inline Foam::tensor& alignment()
{
return alignment_;
}
inline const Foam::tensor& alignment() const
{
return alignment_;
}
inline Foam::scalar& targetCellSize()
{
return targetCellSize_;
}
inline Foam::scalar targetCellSize() const
{
return targetCellSize_;
}
inline bool uninitialised() const
{
return type_ == ptInternalPoint && index_ == ptInternalPoint;
}
//- Is point a far-point
inline bool farPoint() const
{
return type_ == ptFarPoint;
}
//- Is point internal, i.e. not on boundary
inline bool internalPoint() const
{
return type_ <= ptInternalPoint;
}
//- Set the point to be internal
inline void setInternal()
{
type_ = ptInternalPoint;
}
//- Is point internal and near the boundary
inline bool nearBoundary() const
{
return type_ == ptNearBoundaryPoint;
}
//- Set the point to be near the boundary
inline void setNearBoundary()
{
@ -209,12 +225,14 @@ public:
return type_ == ptMirrorPoint;
}
//- Either master or slave of pointPair.
inline bool pairPoint() const
{
return type_ >= 0;
}
//- Master of a pointPair is the lowest numbered one.
inline bool ppMaster() const
{
@ -228,6 +246,7 @@ public:
}
}
//- Slave of a pointPair is the highest numbered one.
inline bool ppSlave() const
{
@ -241,18 +260,21 @@ public:
}
}
//- Either original internal point or master of pointPair.
inline bool internalOrBoundaryPoint() const
{
return internalPoint() || ppMaster();
}
//- Is point near the boundary or part of the boundary definition
inline bool nearOrOnBoundary() const
{
return pairPoint() || mirrorPoint() || nearBoundary();
}
//- Do the two given vertices consitute a boundary point-pair
inline friend bool pointPair
(
@ -263,6 +285,7 @@ public:
return v0.index_ == v1.type_ || v1.index_ == v0.type_;
}
//- Do the three given vertices consitute a boundary triangle
inline friend bool boundaryTriangle
(
@ -276,6 +299,7 @@ public:
|| (v2.pairPoint() && pointPair(v0, v1));
}
//- Do the three given vertices consitute an outside triangle
inline friend bool outsideTriangle
(

20
src/mesh/conformalVoronoiMesh/log
View File

@ -1,20 +0,0 @@
/*---------------------------------------------------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: dev |
| \\ / A nd | Web: www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
Build : dev-79ee164b7a49
Exec : cvMesh
Date : Jan 12 2010
Time : 19:49:15
Host : rockall
PID : 14798
Case : /home/graham/OpenFOAM/OpenFOAM-dev/src/mesh/conformalVoronoiMesh
nProcs : 1
SigFpe : Enabling floating point exception trapping (FOAM_SIGFPE).
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time