ForkMaster/netgenplugin
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

0021140: EDF 1759 SMESH: Netgen1D2D fails on subshape

Browse Source 
fix regression with seam edges
This commit is contained in:
eap 2011-05-31 14:53:31 +00:00
parent e32634a346
commit f12c83a835
1 changed files with 102 additions and 47 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

149
src/NETGENPlugin/NETGENPlugin_Mesher.cxx
View File

@ -368,6 +368,42 @@ void NETGENPlugin_Mesher::PrepareOCCgeometry(netgen::OCCGeometry& occgeo,
occgeo.face_maxh_modified = 0; occgeo.face_maxh_modified = 0;
#endif #endif
// { // set netgen::mparam.minh
// TopLoc_Location loc;
// int i1, i2, i3;
// const int* pi[4] = { &i1, &i2, &i3, &i1 };
// double maxh = 1e100;
// for ( int i = 0; i < occgeo.fmap.Extent(); ++i )
// {
// Handle(Poly_Triangulation) triangulation =
// BRep_Tool::Triangulation ( TopoDS::Face( occgeo.fmap(i+1) ), loc);
// if ( triangulation.IsNull() ) continue;
// const TColgp_Array1OfPnt& points = triangulation->Nodes();
// const Poly_Array1OfTriangle& trias = triangulation->Triangles();
// for ( int iT = trias.Lower(); iT <= trias.Upper(); ++iT )
// {
// trias(iT).Get( i1, i2, i3 );
// for ( int j = 0; j < 3; ++j )
// {
// double dist2 = points(*pi[j]).SquareDistance( points( *pi[j+1] ));
// if ( dist2 < maxh )
// maxh = dist2;
// }
// }
// }
// maxh = sqrt( maxh );
// if ( maxh > 0.5 * occgeo.boundingbox.Diam() ) // no or too rough triangulation
// {
// netgen::mparam.minh = occgeo.boundingbox.Diam()*1e-24;
// cout << "DEFAULT mparams.minh = " <<netgen::mparam.minh << endl;
// }
// else
// {
// netgen::mparam.minh = maxh * 2;
// cout << "TRIANGULATION mparams.minh = " <<netgen::mparam.minh << endl;
// }
// }
} }
namespace namespace
@ -404,10 +440,11 @@ namespace
*/ */
//================================================================================ //================================================================================
list< TopoDS_Edge > getConnectedEdges( const TopoDS_Edge& edge, list< TopoDS_Edge > getConnectedEdges( const TopoDS_Edge& edge,
const TopoDS_Face& face, const TopoDS_Face& face,
const set< SMESH_subMesh* > & computedSM, const set< SMESH_subMesh* > & computedSM,
const SMESH_MesherHelper& helper ) const SMESH_MesherHelper& helper,
map< SMESH_subMesh*, set< int > >& addedEdgeSM2Faces)
{ {
// get ordered EDGEs // get ordered EDGEs
TopoDS_Vertex v1; TopoDS_Vertex v1;
@ -431,48 +468,59 @@ namespace
edges.push_back( e ); edges.push_back( e );
return edges; return edges;
} }
// find not computed or not connected EDGEs around <edge>
// get all computed EDGEs connected to <edge>
list< TopoDS_Edge >::iterator eItBack = eItFwd, ePrev; list< TopoDS_Edge >::iterator eItBack = eItFwd, ePrev;
TopoDS_Vertex vCommon; TopoDS_Vertex vCommon;
TopTools_MapOfShape eAdded; // map used not to add a seam edge twice to <edges>
eAdded.Add( edge );
// put edges after <edge> at <edges> head // put edges before <edge> to <edges> back
while ( edges.back() != *eItFwd ) while ( edges.begin() != eItFwd )
edges.splice( edges.begin(), edges, --edges.end() ); edges.splice( edges.end(), edges, edges.begin() );
// search forward
ePrev = eItFwd;
while ( ++eItFwd != edges.end() )
{
SMESH_subMesh* sm = helper.GetMesh()->GetSubMesh( *eItFwd );
bool connected = TopExp::CommonVertex( *ePrev, *eItFwd, vCommon );
bool computed = sm->IsMeshComputed();
bool added = addedEdgeSM2Faces[ sm ].count( helper.GetSubShapeID() );
bool doubled = !eAdded.Add( *eItFwd );
bool orientOK = (( ePrev ->Orientation() < TopAbs_INTERNAL ) ==
( eItFwd->Orientation() < TopAbs_INTERNAL ) );
if ( !connected || !computed || !orientOK || added || doubled )
{
// stop advancement; move edges from tail to head
while ( edges.back() != *ePrev )
edges.splice( edges.begin(), edges, --edges.end() );
break;
}
ePrev = eItFwd;
}
// search backward // search backward
while ( eItBack != edges.begin() ) while ( eItBack != edges.begin() )
{ {
ePrev = eItBack; ePrev = eItBack;
--eItBack; --eItBack;
SMESH_subMesh* sm = helper.GetMesh()->GetSubMesh( *eItBack );
bool connected = TopExp::CommonVertex( *ePrev, *eItBack, vCommon ); bool connected = TopExp::CommonVertex( *ePrev, *eItBack, vCommon );
bool computed = helper.GetMesh()->GetSubMesh( *eItBack )->IsMeshComputed(); bool computed = sm->IsMeshComputed();
bool added = addedEdgeSM2Faces[ sm ].count( helper.GetSubShapeID() );
bool doubled = !eAdded.Add( *eItBack );
bool orientOK = (( ePrev ->Orientation() < TopAbs_INTERNAL ) == bool orientOK = (( ePrev ->Orientation() < TopAbs_INTERNAL ) ==
( eItBack->Orientation() < TopAbs_INTERNAL ) ); ( eItBack->Orientation() < TopAbs_INTERNAL ) );
if ( !connected || !computed || !orientOK) if ( !connected || !computed || !orientOK || added || doubled)
{ {
// move edges from head to tail // stop advancement
while ( edges.begin() != eItBack ) edges.erase( edges.begin(), ePrev );
edges.splice( edges.end(), edges, edges.begin() );
edges.erase( eItBack );
break; break;
} }
} }
// search forward
ePrev = eItFwd;
while ( ++eItFwd != edges.end() )
{
bool connected = TopExp::CommonVertex( *ePrev, *eItFwd, vCommon );
bool computed = helper.GetMesh()->GetSubMesh( *eItFwd )->IsMeshComputed();
bool orientOK = (( ePrev ->Orientation() < TopAbs_INTERNAL ) ==
( eItFwd->Orientation() < TopAbs_INTERNAL ) );
if ( !connected || !computed || !orientOK )
{
edges.erase( eItFwd, edges.end() );
break;
}
ePrev = eItFwd;
}
if ( edges.front() != edges.back() ) if ( edges.front() != edges.back() )
{ {
// assure that the 1st vertex is meshed // assure that the 1st vertex is meshed
@ -508,7 +556,7 @@ bool NETGENPlugin_Mesher::fillNgMesh(const netgen::OCCGeometry& occgeom,
SMESH_MesherHelper helper (*_mesh); SMESH_MesherHelper helper (*_mesh);
int faceID = occgeom.fmap.Extent(); int faceNgID = occgeom.fmap.Extent();
list< SMESH_subMesh* >::const_iterator smIt, smEnd = meshedSM.end(); list< SMESH_subMesh* >::const_iterator smIt, smEnd = meshedSM.end();
for ( smIt = meshedSM.begin(); smIt != smEnd; ++smIt ) for ( smIt = meshedSM.begin(); smIt != smEnd; ++smIt )
@ -532,19 +580,24 @@ bool NETGENPlugin_Mesher::fillNgMesh(const netgen::OCCGeometry& occgeom,
PShapeIteratorPtr fIt = helper.GetAncestors( geomEdge, *sm->GetFather(), TopAbs_FACE ); PShapeIteratorPtr fIt = helper.GetAncestors( geomEdge, *sm->GetFather(), TopAbs_FACE );
while ( const TopoDS_Shape * anc = fIt->next() ) while ( const TopoDS_Shape * anc = fIt->next() )
{ {
int faceID = occgeom.fmap.FindIndex( *anc ); faceNgID = occgeom.fmap.FindIndex( *anc );
if ( faceID < 1 ) if ( faceNgID < 1 )
continue; // meshed face continue; // meshed face
if ( visitedEdgeSM2Faces[ sm ].count( faceID ))
int faceSMDSId = helper.GetMeshDS()->ShapeToIndex( *anc );
if ( visitedEdgeSM2Faces[ sm ].count( faceSMDSId ))
continue; // already treated EDGE continue; // already treated EDGE
TopoDS_Face face = TopoDS::Face( occgeom.fmap( faceID )); TopoDS_Face face = TopoDS::Face( occgeom.fmap( faceNgID ));
if ( face.Orientation() >= TopAbs_INTERNAL ) if ( face.Orientation() >= TopAbs_INTERNAL )
face.Orientation( TopAbs_FORWARD ); // issue 0020676 face.Orientation( TopAbs_FORWARD ); // issue 0020676
// get all meshed EDGEs of the FACE connected to geomEdge (issue 0021140) // get all meshed EDGEs of the FACE connected to geomEdge (issue 0021140)
helper.SetSubShape( face ); helper.SetSubShape( face );
list< TopoDS_Edge > edges = getConnectedEdges( geomEdge, face, computedSM, helper ); list< TopoDS_Edge > edges = getConnectedEdges( geomEdge, face, computedSM, helper,
visitedEdgeSM2Faces );
if ( edges.empty() )
continue; // wrong ancestor?
// find out orientation of <edges> within <face> // find out orientation of <edges> within <face>
TopoDS_Edge eNotSeam = edges.front(); TopoDS_Edge eNotSeam = edges.front();
@ -564,9 +617,9 @@ bool NETGENPlugin_Mesher::fillNgMesh(const netgen::OCCGeometry& occgeom,
const vector<UVPtStruct>& points = fSide.GetUVPtStruct(); const vector<UVPtStruct>& points = fSide.GetUVPtStruct();
int i, nbSeg = fSide.NbSegments(); int i, nbSeg = fSide.NbSegments();
// remembre EDGEs of fSide to treat only once // remember EDGEs of fSide to treat only once
for ( int iE = 0; iE < fSide.NbEdges(); ++iE ) for ( int iE = 0; iE < fSide.NbEdges(); ++iE )
visitedEdgeSM2Faces[ helper.GetMesh()->GetSubMesh( fSide.Edge(iE )) ].insert( faceID ); visitedEdgeSM2Faces[ helper.GetMesh()->GetSubMesh( fSide.Edge(iE )) ].insert(faceSMDSId);
double otherSeamParam = 0; double otherSeamParam = 0;
bool isSeam = false; bool isSeam = false;
@ -585,10 +638,12 @@ bool NETGENPlugin_Mesher::fillNgMesh(const netgen::OCCGeometry& occgeom,
isSeam = false; isSeam = false;
if ( helper.IsRealSeam( p1.node->getshapeId() )) if ( helper.IsRealSeam( p1.node->getshapeId() ))
{ {
geomEdge = fSide.Edge( fSide.EdgeIndex( 0.5 * ( p1.normParam + p2.normParam ))); TopoDS_Edge e = fSide.Edge( fSide.EdgeIndex( 0.5 * ( p1.normParam + p2.normParam )));
isSeam = helper.IsRealSeam( geomEdge ); isSeam = helper.IsRealSeam( e );
if ( isSeam ) if ( isSeam )
{
otherSeamParam = helper.GetOtherParam( helper.GetPeriodicIndex() & 1 ? p2.u : p2.v ); otherSeamParam = helper.GetOtherParam( helper.GetPeriodicIndex() & 1 ? p2.u : p2.v );
}
} }
} }
netgen::Segment seg; netgen::Segment seg;
@ -608,7 +663,7 @@ bool NETGENPlugin_Mesher::fillNgMesh(const netgen::OCCGeometry& occgeom,
//seg.epgeominfo[ 0 ].edgenr = seg.epgeominfo[ 1 ].edgenr = occgeom.emap.FindIndex( geomEdge ); //seg.epgeominfo[ 0 ].edgenr = seg.epgeominfo[ 1 ].edgenr = occgeom.emap.FindIndex( geomEdge );
//seg.epgeominfo[ iEnd ].edgenr = edgeID; // = geom.emap.FindIndex(edge); //seg.epgeominfo[ iEnd ].edgenr = edgeID; // = geom.emap.FindIndex(edge);
seg.si = faceID; // = geom.fmap.FindIndex (face); seg.si = faceNgID; // = geom.fmap.FindIndex (face);
seg.edgenr = ngMesh.GetNSeg() + 1; // segment id seg.edgenr = ngMesh.GetNSeg() + 1; // segment id
ngMesh.AddSegment (seg); ngMesh.AddSegment (seg);
@ -681,11 +736,11 @@ bool NETGENPlugin_Mesher::fillNgMesh(const netgen::OCCGeometry& occgeom,
if ( solidID1 && id != solidID1 ) solidID2 = id; if ( solidID1 && id != solidID1 ) solidID2 = id;
else solidID1 = id; else solidID1 = id;
} }
faceID++; faceNgID++;
//_faceDescriptors[ faceID ].first = solidID1; //_faceDescriptors[ faceNgID ].first = solidID1;
//_faceDescriptors[ faceID ].second = solidID2; //_faceDescriptors[ faceNgID ].second = solidID2;
// Add ng face descriptors of meshed faces // Add ng face descriptors of meshed faces
ngMesh.AddFaceDescriptor (netgen::FaceDescriptor(faceID, solidID1, solidID2, 0)); ngMesh.AddFaceDescriptor (netgen::FaceDescriptor(faceNgID, solidID1, solidID2, 0));
// Orient the face correctly in solidID1 (issue 0020206) // Orient the face correctly in solidID1 (issue 0020206)
bool reverse = false; bool reverse = false;
@ -700,7 +755,7 @@ bool NETGENPlugin_Mesher::fillNgMesh(const netgen::OCCGeometry& occgeom,
// Add surface elements // Add surface elements
netgen::Element2d tri(3); netgen::Element2d tri(3);
tri.SetIndex ( faceID ); tri.SetIndex ( faceNgID );
#ifdef DUMP_TRIANGLES #ifdef DUMP_TRIANGLES
@ -1657,8 +1712,8 @@ bool NETGENPlugin_Mesher::Compute()
int key = (*it).first; int key = (*it).first;
double val = (*it).second; double val = (*it).second;
const TopoDS_Shape& shape = ShapesWithLocalSize.FindKey(key); const TopoDS_Shape& shape = ShapesWithLocalSize.FindKey(key);
int faceID = occgeo.fmap.FindIndex(shape); int faceNgID = occgeo.fmap.FindIndex(shape);
occgeo.SetFaceMaxH(faceID, val); occgeo.SetFaceMaxH(faceNgID, val);
} }
} }
#endif #endif