-// Copyright (C) 2007-2008 CEA/DEN, EDF R&D, OPEN CASCADE
+// Copyright (C) 2007-2010 CEA/DEN, EDF R&D, OPEN CASCADE
//
// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
//
// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
//
+
//=============================================================================
// File : NETGENPlugin_NETGEN_3D.cxx
// Moved here from SMESH_NETGEN_3D.cxx
Netgen include files
*/
+#define OCCGEOMETRY
+#include <occgeom.hpp>
namespace nglib {
#include <nglib.h>
}
*/
//=============================================================================
-bool NETGENPlugin_NETGEN_3D::CheckHypothesis
- (SMESH_Mesh& aMesh,
- const TopoDS_Shape& aShape,
- SMESH_Hypothesis::Hypothesis_Status& aStatus)
+bool NETGENPlugin_NETGEN_3D::CheckHypothesis (SMESH_Mesh& aMesh,
+ const TopoDS_Shape& aShape,
+ Hypothesis_Status& aStatus)
{
MESSAGE("NETGENPlugin_NETGEN_3D::CheckHypothesis");
return isOk;
}
-namespace
-{
- //================================================================================
- /*!
- * \brief It correctly initializes netgen library at constructor and
- * correctly finishes using netgen library at destructor
- */
- //================================================================================
-
- struct TNetgenLibWrapper
- {
- Ng_Mesh* _ngMesh;
- TNetgenLibWrapper()
- {
- Ng_Init();
- _ngMesh = Ng_NewMesh();
- }
- ~TNetgenLibWrapper()
- {
- Ng_DeleteMesh( _ngMesh );
- Ng_Exit();
- NETGENPlugin_Mesher::RemoveTmpFiles();
- }
- };
-
- //================================================================================
- /*!
- * \brief Find mesh faces on non-internal geom faces sharing internal edge
- * some nodes of which are to be doubled to make the second border of the "crack"
- */
- //================================================================================
-
- void findBorders( const set<int>& internalShapeIds,
- SMESH_MesherHelper& helper,
- TIDSortedElemSet & borderElems,
- set<int> & borderFaceIds )
- {
- SMESH_Mesh* mesh = helper.GetMesh();
- SMESHDS_Mesh* meshDS = helper.GetMeshDS();
-
- // loop on internal geom edges
- set<int>::const_iterator intShapeId = internalShapeIds.begin();
- for ( ; intShapeId != internalShapeIds.end(); ++intShapeId )
- {
- const TopoDS_Shape& s = meshDS->IndexToShape( *intShapeId );
- if ( s.ShapeType() != TopAbs_EDGE ) continue;
-
- // get internal and non-internal geom faces sharing the internal edge <s>
- int intFace = 0;
- set<int>::iterator bordFace = borderFaceIds.end();
- TopTools_ListIteratorOfListOfShape ancIt( mesh->GetAncestors( s ));
- for ( ; ancIt.More(); ancIt.Next() )
- {
- if ( ancIt.Value().ShapeType() != TopAbs_FACE ) continue;
- int faceID = meshDS->ShapeToIndex( ancIt.Value() );
- if ( internalShapeIds.count( faceID ))
- intFace = faceID;
- else
- bordFace = borderFaceIds.insert( faceID ).first;
- }
- if ( bordFace == borderFaceIds.end() || !intFace ) continue;
-
- // get all links of mesh faces on internal geom face sharing nodes on edge <s>
- set< SMESH_OrientedLink > links; //!< links of faces on internal geom face
- TIDSortedElemSet suspectFaces; //!< mesh faces on border geom faces
- SMESHDS_SubMesh* intFaceSM = meshDS->MeshElements( intFace );
- if ( !intFaceSM || intFaceSM->NbElements() == 0 ) continue;
- SMESH_subMeshIteratorPtr smIt = mesh->GetSubMesh( s )->getDependsOnIterator(true,true);
- while ( smIt->more() )
- {
- SMESHDS_SubMesh* sm = smIt->next()->GetSubMeshDS();
- if ( !sm ) continue;
- SMDS_NodeIteratorPtr nIt = sm->GetNodes();
- while ( nIt->more() )
- {
- const SMDS_MeshNode* nOnEdge = nIt->next();
- SMDS_ElemIteratorPtr fIt = nOnEdge->GetInverseElementIterator(SMDSAbs_Face);
- while ( fIt->more() )
- {
- const SMDS_MeshElement* f = fIt->next();
- if ( intFaceSM->Contains( f ))
- {
- int nbNodes = f->NbNodes() / ( f->IsQuadratic() ? 2 : 1 );
- for ( int i = 0; i < nbNodes; ++i )
- links.insert( SMESH_OrientedLink( f->GetNode(i), f->GetNode((i+1)%nbNodes)));
- }
- else
- {
- suspectFaces.insert( f );
- }
- }
- }
- }
- // <suspectFaces> having link with same orientation as mesh faces on
- // the internal geom face are <borderElems>.
- // Some of them have only one node on edge s, we collect them to decide
- // later by links of found <borderElems>
- TIDSortedElemSet posponedFaces;
- set< SMESH_OrientedLink > borderLinks;
- TIDSortedElemSet::iterator fIt = suspectFaces.begin();
- for ( ; fIt != suspectFaces.end(); ++fIt )
- {
- const SMDS_MeshElement* f = *fIt;
- bool linkFound = false, isBorder = false;
- list< SMESH_OrientedLink > faceLinks;
- int nbNodes = f->NbNodes() / ( f->IsQuadratic() ? 2 : 1 );
- for ( int i = 0; i < nbNodes; ++i )
- {
- SMESH_OrientedLink link( f->GetNode(i), f->GetNode((i+1)%nbNodes));
- faceLinks.push_back( link );
- if ( !linkFound )
- {
- set< SMESH_OrientedLink >::iterator foundLink = links.find( link );
- if ( foundLink != links.end() )
- {
- linkFound= true;
- isBorder = ( foundLink->_reversed == link._reversed );
- if ( !isBorder ) break;
- }
- }
- }
- if ( isBorder )
- {
- borderElems.insert( f );
- borderLinks.insert( faceLinks.begin(), faceLinks.end() );
- }
- else if ( !linkFound )
- {
- posponedFaces.insert( f );
- }
- }
- // decide on posponedFaces
- for ( fIt = posponedFaces.begin(); fIt != posponedFaces.end(); ++fIt )
- {
- const SMDS_MeshElement* f = *fIt;
- int nbNodes = f->NbNodes() / ( f->IsQuadratic() ? 2 : 1 );
- for ( int i = 0; i < nbNodes; ++i )
- {
- SMESH_OrientedLink link( f->GetNode(i), f->GetNode((i+1)%nbNodes));
- set< SMESH_OrientedLink >::iterator foundLink = borderLinks.find( link );
- if ( foundLink != borderLinks.end() )
- {
- if ( foundLink->_reversed != link._reversed )
- borderElems.insert( f );
- break;
- }
- }
- }
- }
- }
-}
-
//=============================================================================
/*!
*Here we are going to use the NETGEN mesher
SMESH_MesherHelper helper(aMesh);
bool _quadraticMesh = helper.IsQuadraticSubMesh(aShape);
+ helper.SetElementsOnShape( true );
int Netgen_NbOfNodes = 0;
int Netgen_param2ndOrder = 0;
int Netgen_triangle[3];
int Netgen_tetrahedron[4];
- TNetgenLibWrapper ngLib;
+ NETGENPlugin_NetgenLibWrapper ngLib;
Ng_Mesh * Netgen_mesh = ngLib._ngMesh;
- // maps of 1) ordinary nodes and 2) doubled nodes on internal shapes
- typedef map< const SMDS_MeshNode*, int, TIDCompare > TNodeToIDMap;
- typedef TNodeToIDMap::value_type TN2ID;
- TNodeToIDMap nodeToNetgenID[2];
-
+ // vector of nodes in which node index == netgen ID
+ vector< const SMDS_MeshNode* > nodeVec;
{
const int invalid_ID = -1;
SMESH::Controls::Area areaControl;
SMESH::Controls::TSequenceOfXYZ nodesCoords;
- // --------------------------------------------------------------------
- // Issue 0020676 (StudyFiss_bugNetgen3D.hdf). Pb with internal face.
- // Find internal geom faces, edges and vertices.
- // Nodes and faces built on the found internal shapes
- // will be doubled in Netgen input to make two borders of the "crack".
- // --------------------------------------------------------------------
-
- set<int> internalShapeIds;
- set<int> borderFaceIds; //!< non-internal geom faces sharing internal edge
-
- // mesh faces on <borderFaceIds>, having nodes on internal edge that are
- // to be replaced by doubled nodes
- TIDSortedElemSet borderElems;
+ // maps nodes to ng ID
+ typedef map< const SMDS_MeshNode*, int, TIDCompare > TNodeToIDMap;
+ typedef TNodeToIDMap::value_type TN2ID;
+ TNodeToIDMap nodeToNetgenID;
- // find "internal" faces and edges
- TopExp_Explorer exFa(aShape,TopAbs_FACE), exEd, exVe;
- for ( ; exFa.More(); exFa.Next())
- {
- if ( exFa.Current().Orientation() == TopAbs_INTERNAL )
- {
- internalShapeIds.insert( meshDS->ShapeToIndex( exFa.Current() ));
- for ( exEd.Init( exFa.Current(), TopAbs_EDGE ); exEd.More(); exEd.Next())
- if ( helper.NbAncestors( exEd.Current(), aMesh, TopAbs_FACE ) > 1 )
- internalShapeIds.insert( meshDS->ShapeToIndex( exEd.Current() ));
- }
- }
- if ( !internalShapeIds.empty() )
- {
- // find internal vertices,
- // we consider vertex internal if it is shared by more than one internal edge
- TopTools_ListIteratorOfListOfShape ancIt;
- set<int>::iterator intShapeId = internalShapeIds.begin();
- for ( ; intShapeId != internalShapeIds.end(); ++intShapeId )
- {
- const TopoDS_Shape& s = meshDS->IndexToShape( *intShapeId );
- if ( s.ShapeType() != TopAbs_EDGE ) continue;
- for ( exVe.Init( s, TopAbs_VERTEX ); exVe.More(); exVe.Next())
- {
- set<int> internalEdges;
- for ( ancIt.Initialize( aMesh.GetAncestors( exVe.Current() ));
- ancIt.More(); ancIt.Next() )
- {
- if ( ancIt.Value().ShapeType() != TopAbs_EDGE ) continue;
- int edgeID = meshDS->ShapeToIndex( ancIt.Value() );
- if ( internalShapeIds.count( edgeID ))
- internalEdges.insert( edgeID );
- }
- if ( internalEdges.size() > 1 )
- internalShapeIds.insert( meshDS->ShapeToIndex( exVe.Current() ));
- }
- }
- // find border shapes and mesh elements
- findBorders( internalShapeIds, helper, borderElems, borderFaceIds );
- }
+ // find internal shapes
+ NETGENPlugin_Internals internals( aMesh, aShape, /*is3D=*/true );
// ---------------------------------
// Feed the Netgen with surface mesh
if ( aMesh.NbQuadrangles() > 0 )
Adaptor.Compute(aMesh,aShape);
- for ( exFa.ReInit(); exFa.More(); exFa.Next())
+ for ( TopExp_Explorer exFa( aShape, TopAbs_FACE ); exFa.More(); exFa.Next())
{
const TopoDS_Shape& aShapeFace = exFa.Current();
int faceID = meshDS->ShapeToIndex( aShapeFace );
- bool isInternalFace = internalShapeIds.count( faceID );
- bool isBorderFace = borderFaceIds.count( faceID );
+ bool isInternalFace = internals.isInternalShape( faceID );
bool isRev = false;
if ( checkReverse && !isInternalFace &&
helper.NbAncestors(aShapeFace, aMesh, aShape.ShapeType()) > 1 )
}
// Add nodes of triangles and triangles them-selves to netgen mesh
- // a triangle on internal face is added twice,
- // on border face, once but with doubled nodes
- bool isBorder = ( isBorderFace && borderElems.count( elem ));
- int nbDblLoops = ( isInternalFace && isTraingle || isBorder ) ? 2 : 1;
-
for ( int i = 0; i < trias.size(); ++i )
{
- bool reverse = isRev;
- for ( int isDblF = isBorder; isDblF < nbDblLoops; ++isDblF, reverse = !reverse )
+ // add three nodes of triangle
+ bool hasDegen = false;
+ for ( int iN = 0; iN < 3; ++iN )
{
- // add three nodes of triangle
- bool hasDegen = false;
- for ( int iN = 0; iN < 3; ++iN )
+ const SMDS_MeshNode* node = trias[i]->GetNode( iN );
+ int shapeID = node->GetPosition()->GetShapeId();
+ if ( node->GetPosition()->GetTypeOfPosition() == SMDS_TOP_EDGE &&
+ helper.IsDegenShape( shapeID ))
+ {
+ // ignore all nodes on degeneraged edge and use node on its vertex instead
+ TopoDS_Shape vertex = TopoDS_Iterator( meshDS->IndexToShape( shapeID )).Value();
+ node = SMESH_Algo::VertexNode( TopoDS::Vertex( vertex ), meshDS );
+ hasDegen = true;
+ }
+ int& ngID = nodeToNetgenID.insert(TN2ID( node, invalid_ID )).first->second;
+ if ( ngID == invalid_ID )
{
- const SMDS_MeshNode* node = trias[i]->GetNode( iN );
- int shapeID = node->GetPosition()->GetShapeId();
- if ( node->GetPosition()->GetTypeOfPosition() == SMDS_TOP_EDGE &&
- helper.IsDegenShape( shapeID ))
- {
- // ignore all nodes on degeneraged edge and use node on its vertex instead
- TopoDS_Shape vertex = TopoDS_Iterator( meshDS->IndexToShape( shapeID )).Value();
- node = SMESH_Algo::VertexNode( TopoDS::Vertex( vertex ), meshDS );
- hasDegen = true;
- }
- bool isDblN = isDblF && internalShapeIds.count( shapeID );
- int& ngID = nodeToNetgenID[isDblN].insert(TN2ID( node, invalid_ID )).first->second;
- if ( ngID == invalid_ID )
- {
- ngID = ++Netgen_NbOfNodes;
- Netgen_point [ 0 ] = node->X();
- Netgen_point [ 1 ] = node->Y();
- Netgen_point [ 2 ] = node->Z();
- Ng_AddPoint(Netgen_mesh, Netgen_point);
- }
- Netgen_triangle[ iN ] = ngID;
+ ngID = ++Netgen_NbOfNodes;
+ Netgen_point [ 0 ] = node->X();
+ Netgen_point [ 1 ] = node->Y();
+ Netgen_point [ 2 ] = node->Z();
+ Ng_AddPoint(Netgen_mesh, Netgen_point);
}
- // add triangle
- if ( hasDegen && (Netgen_triangle[0] == Netgen_triangle[1] ||
- Netgen_triangle[0] == Netgen_triangle[2] ||
- Netgen_triangle[2] == Netgen_triangle[1] ))
- break;
- if ( reverse )
- swap( Netgen_triangle[1], Netgen_triangle[2] );
+ Netgen_triangle[ isRev ? 2-iN : iN ] = ngID;
+ }
+ // add triangle
+ if ( hasDegen && (Netgen_triangle[0] == Netgen_triangle[1] ||
+ Netgen_triangle[0] == Netgen_triangle[2] ||
+ Netgen_triangle[2] == Netgen_triangle[1] ))
+ continue;
+
+ Ng_AddSurfaceElement(Netgen_mesh, NG_TRIG, Netgen_triangle);
+ if ( isInternalFace && isTraingle )
+ {
+ swap( Netgen_triangle[1], Netgen_triangle[2] );
Ng_AddSurfaceElement(Netgen_mesh, NG_TRIG, Netgen_triangle);
}
}
} // loop on elements on a face
} // loop on faces of a SOLID or SHELL
+ // insert old nodes into nodeVec
+ nodeVec.resize( nodeToNetgenID.size() + 1, 0 );
+ TNodeToIDMap::iterator n_id = nodeToNetgenID.begin();
+ for ( ; n_id != nodeToNetgenID.end(); ++n_id )
+ nodeVec[ n_id->second ] = n_id->first;
+ nodeToNetgenID.clear();
+
+ if ( internals.hasInternalVertexInSolid() )
+ {
+ netgen::OCCGeometry occgeo;
+ NETGENPlugin_Mesher::addIntVerticesInSolids( occgeo,
+ (netgen::Mesh&) *Netgen_mesh,
+ nodeVec,
+ internals);
+ Netgen_NbOfNodes = Ng_GetNP(Netgen_mesh);
+ }
}
// -------------------------
}
int Netgen_NbOfNodesNew = Ng_GetNP(Netgen_mesh);
-
- int Netgen_NbOfTetra = Ng_GetNE(Netgen_mesh);
+ int Netgen_NbOfTetra = Ng_GetNE(Netgen_mesh);
MESSAGE("End of Volume Mesh Generation. status=" << status <<
", nb new nodes: " << Netgen_NbOfNodesNew - Netgen_NbOfNodes <<
// Feed back the SMESHDS with the generated Nodes and Volume Elements
// -------------------------------------------------------------------
- // vector of nodes in which node index == netgen ID
- vector< const SMDS_MeshNode* > nodeVec ( Netgen_NbOfNodesNew + 1 );
- // insert old nodes into nodeVec
- for ( int isDbl = 0; isDbl < 2; ++isDbl )
- {
- TNodeToIDMap::iterator n_id = nodeToNetgenID[isDbl].begin();
- for ( ; n_id != nodeToNetgenID[isDbl].end(); ++n_id )
- nodeVec[ n_id->second ] = n_id->first;
- nodeToNetgenID[isDbl].clear();
- }
if ( status == NG_VOLUME_FAILURE )
{
SMESH_ComputeErrorPtr err = NETGENPlugin_Mesher::readErrors(nodeVec);
if ( isOK )
{
// create and insert new nodes into nodeVec
+ nodeVec.resize( Netgen_NbOfNodesNew + 1, 0 );
int nodeIndex = Netgen_NbOfNodes + 1;
- int shapeID = meshDS->ShapeToIndex( aShape );
for ( ; nodeIndex <= Netgen_NbOfNodesNew; ++nodeIndex )
{
Ng_GetPoint( Netgen_mesh, nodeIndex, Netgen_point );
- SMDS_MeshNode * node = meshDS->AddNode(Netgen_point[0],
- Netgen_point[1],
- Netgen_point[2]);
- meshDS->SetNodeInVolume(node, shapeID);
- nodeVec.at(nodeIndex) = node;
+ nodeVec.at(nodeIndex) = helper.AddNode(Netgen_point[0], Netgen_point[1], Netgen_point[2]);
}
// create tetrahedrons
for ( int elemIndex = 1; elemIndex <= Netgen_NbOfTetra; ++elemIndex )
{
Ng_GetVolumeElement(Netgen_mesh, elemIndex, Netgen_tetrahedron);
- SMDS_MeshVolume * elt = helper.AddVolume (nodeVec.at( Netgen_tetrahedron[0] ),
- nodeVec.at( Netgen_tetrahedron[1] ),
- nodeVec.at( Netgen_tetrahedron[2] ),
- nodeVec.at( Netgen_tetrahedron[3] ));
- meshDS->SetMeshElementOnShape(elt, shapeID );
+ helper.AddVolume (nodeVec.at( Netgen_tetrahedron[0] ),
+ nodeVec.at( Netgen_tetrahedron[1] ),
+ nodeVec.at( Netgen_tetrahedron[2] ),
+ nodeVec.at( Netgen_tetrahedron[3] ));
}
}
// using adaptor
const list<const SMDS_FaceOfNodes*>* faces = Adaptor.GetTriangles(elem);
if(faces==0)
- return error( COMPERR_BAD_INPUT_MESH,
- SMESH_Comment("No triangles in adaptor for element ")<<elem->GetID());
+ continue; // Issue 0020682. There already can be 3d mesh
trias.assign( faces->begin(), faces->end() );
}
else {
int Netgen_triangle[3];
int Netgen_tetrahedron[4];
- TNetgenLibWrapper ngLib;
+ NETGENPlugin_NetgenLibWrapper ngLib;
Ng_Mesh * Netgen_mesh = ngLib._ngMesh;
// set nodes and remember thier netgen IDs
