Commit 35fecf7e authored by Franjo's avatar Franjo
Browse files

Face-edges calculation tolerant of topological flaws

parent 815ff5f4
......@@ -29,6 +29,8 @@ Description
#include "VRWGraphSMPModifier.H"
#include "demandDrivenData.H"
#include <set>
# ifdef USE_OMP
#include <omp.h>
# endif
......@@ -70,49 +72,42 @@ void triSurfAddressing::calculateEdges() const
# ifdef USE_OMP
# pragma omp for schedule(static)
# endif
forAll(facets_, fI)
forAll(pFacets, pI)
{
const labelledTri& tri = facets_[fI];
std::set<std::pair<label, label> > edgesAtPoint;
forAll(tri, pI)
forAllRow(pFacets, pI, pfI)
{
const edge fe(tri[pI], tri[(pI+1)%3]);
const label s = fe.start();
const label e = fe.end();
DynList<label> edgeFaces;
bool store(true);
const label triI = pFacets(pI, pfI);
const labelledTri& tri = facets_[triI];
//- find all faces attached to this edge
//- store the edge in case the face faceI is the face
//- with the smallest label
forAllRow(pFacets, s, pfI)
label pos(-1);
forAll(tri, i)
{
const label ofI = pFacets(s, pfI);
const labelledTri& of = facets_[ofI];
label pos(-1);
for(label i=0;i<3;++i)
if( of[i] == e )
{
pos = i;
break;
}
if( pos < 0 )
continue;
if( ofI < fI )
if( tri[i] == pI )
{
store = false;
break;
pos = i;
if( tri[(pos+1)%3] >= pI )
edgesAtPoint.insert
(
std::make_pair(pI, tri[(pos+1)%3])
);
if( tri[(pos+2)%3] >= pI )
edgesAtPoint.insert
(
std::make_pair(pI, tri[(pos+2)%3])
);
}
edgeFaces.append(ofI);
}
if( store )
edgesHelper.append(fe);
}
std::set<std::pair<label, label> >::const_iterator it;
for(it=edgesAtPoint.begin();it!=edgesAtPoint.end();++it)
edgesHelper.append(edge(it->first, it->second));
}
//- this enables other threads to see the number of edges
......@@ -158,56 +153,58 @@ void triSurfAddressing::calculateFacetEdges() const
const edgeLongList& edges = this->edges();
const VRWGraph& pointFaces = this->pointFacets();
facetEdgesPtr_ = new VRWGraph(facets_.size());
facetEdgesPtr_ = new VRWGraph(facets_.size(), 3, -1);
VRWGraph& faceEdges = *facetEdgesPtr_;
labelList nfe(facets_.size());
# ifdef USE_OMP
const label nThreads = 3 * omp_get_num_procs();
# pragma omp parallel num_threads(nThreads)
# pragma omp parallel for schedule(dynamic, 100)
# endif
forAll(edges, edgeI)
{
# ifdef USE_OMP
# pragma omp for schedule(static)
# endif
forAll(facets_, fI)
nfe[fI] = 3;
# ifdef USE_OMP
# pragma omp barrier
# pragma omp master
# endif
VRWGraphSMPModifier(faceEdges).setSizeAndRowSize(nfe);
# ifdef USE_OMP
# pragma omp barrier
const edge ee = edges[edgeI];
const label pI = ee.start();
# pragma omp for schedule(static)
# endif
forAll(edges, edgeI)
forAllRow(pointFaces, pI, pfI)
{
const edge ee = edges[edgeI];
const label pI = ee.start();
const label fI = pointFaces(pI, pfI);
forAllRow(pointFaces, pI, pfI)
const labelledTri& tri = facets_[fI];
forAll(tri, eI)
{
const label fI = pointFaces(pI, pfI);
const edge e(tri[eI], tri[(eI+1)%3]);
const labelledTri& tri = facets_[fI];
forAll(tri, eI)
if( e == ee )
{
const edge e(tri[eI], tri[(eI+1)%3]);
if( e == ee )
{
faceEdges(fI, eI) = edgeI;
break;
}
faceEdges(fI, eI) = edgeI;
}
}
}
}
# ifdef DEBUGTriSurfAddressing
forAll(faceEdges, triI)
{
forAllRow(faceEdges, triI, feI)
{
if( facets_[triI][feI] < 0 || facets_[triI][feI] >= points_.size() )
FatalErrorIn
(
"void triSurfAddressing::calculateFacetEdges() const"
) << "Invalid entry in triangle " << triI
<< " " << facets_[triI] << abort(FatalError);
const label edgeI = faceEdges(triI, feI);
if( edgeI < 0 || edgeI >= edges.size() )
FatalErrorIn
(
"void triSurfAddressing::calculateFacetEdges() const"
) << "Invalid entry in face " << triI << " "
<< facets_[triI] << " edges "
<< faceEdges[triI] << abort(FatalError);
}
}
# endif
}
void triSurfAddressing::calculateEdgeFacets() const
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment