From b5527bc7f9c2c15c426415a7558d2a608a477f23 Mon Sep 17 00:00:00 2001 From: eap Date: Thu, 30 Dec 2010 09:32:03 +0000 Subject: [PATCH] 0020832: EDF 1359 SMESH : Automatic meshing of boundary layers for NETGEN_3D to work after StdMeshers_ViscousLayers + opt-hypos="ViscousLayers" --- resources/NETGENPlugin.xml | 2 + src/NETGENPlugin/NETGENPlugin_NETGEN_3D.cxx | 262 ++++++++------------ src/NETGENPlugin/NETGENPlugin_NETGEN_3D.hxx | 7 +- 3 files changed, 114 insertions(+), 157 deletions(-) diff --git a/resources/NETGENPlugin.xml b/resources/NETGENPlugin.xml index 0d3ec83..92657a4 100644 --- a/resources/NETGENPlugin.xml +++ b/resources/NETGENPlugin.xml @@ -55,6 +55,7 @@ label-id="Tetrahedron (Netgen)" icon-id="mesh_algo_tetra.png" hypos="MaxElementVolume" + opt-hypos="ViscousLayers" need-geom="false" input="TRIA,QUAD" dim="3"/> @@ -70,6 +71,7 @@ label-id="Netgen 1D-2D-3D" icon-id="mesh_algo_netgen_2d3d.png" hypos="NETGEN_Parameters, NETGEN_SimpleParameters_3D" + opt-hypos="ViscousLayers" dim="3" support-submeshes="true" /> diff --git a/src/NETGENPlugin/NETGENPlugin_NETGEN_3D.cxx b/src/NETGENPlugin/NETGENPlugin_NETGEN_3D.cxx index fb0a1bd..d6faefb 100644 --- a/src/NETGENPlugin/NETGENPlugin_NETGEN_3D.cxx +++ b/src/NETGENPlugin/NETGENPlugin_NETGEN_3D.cxx @@ -39,9 +39,11 @@ #include "SMESH_ControlsDef.hxx" #include "SMESH_Gen.hxx" #include "SMESH_Mesh.hxx" -#include "SMESH_MesherHelper.hxx" #include "SMESH_MeshEditor.hxx" +#include "SMESH_MesherHelper.hxx" +#include "StdMeshers_MaxElementVolume.hxx" #include "StdMeshers_QuadToTriaAdaptor.hxx" +#include "StdMeshers_ViscousLayers.hxx" #include #include @@ -78,14 +80,14 @@ using namespace std; */ //============================================================================= -NETGENPlugin_NETGEN_3D::NETGENPlugin_NETGEN_3D(int hypId, int studyId, - SMESH_Gen* gen) +NETGENPlugin_NETGEN_3D::NETGENPlugin_NETGEN_3D(int hypId, int studyId, SMESH_Gen* gen) : SMESH_3D_Algo(hypId, studyId, gen) { MESSAGE("NETGENPlugin_NETGEN_3D::NETGENPlugin_NETGEN_3D"); _name = "NETGEN_3D"; _shapeType = (1 << TopAbs_SHELL) | (1 << TopAbs_SOLID);// 1 bit /shape type _compatibleHypothesis.push_back("MaxElementVolume"); + _compatibleHypothesis.push_back("ViscousLayers"); _maxElementVolume = 0.; @@ -118,39 +120,35 @@ bool NETGENPlugin_NETGEN_3D::CheckHypothesis (SMESH_Mesh& aMesh, MESSAGE("NETGENPlugin_NETGEN_3D::CheckHypothesis"); _hypMaxElementVolume = NULL; + _viscousLayersHyp = NULL; _maxElementVolume = DBL_MAX; - list::const_iterator itl; - const SMESHDS_Hypothesis* theHyp; - - const list& hyps = GetUsedHypothesis(aMesh, aShape); - int nbHyp = hyps.size(); - if (!nbHyp) + const list& hyps = + GetUsedHypothesis(aMesh, aShape, /*ignoreAuxiliary=*/false); + list ::const_iterator h = hyps.begin(); + if ( h == hyps.end()) { aStatus = SMESH_Hypothesis::HYP_OK; - //aStatus = SMESH_Hypothesis::HYP_MISSING; return true; // can work with no hypothesis } - itl = hyps.begin(); - theHyp = (*itl); // use only the first hypothesis - - string hypName = theHyp->GetName(); - - bool isOk = false; - - if (hypName == "MaxElementVolume") + aStatus = HYP_OK; + for ( ; h != hyps.end(); ++h ) { - _hypMaxElementVolume = static_cast (theHyp); - ASSERT(_hypMaxElementVolume); - _maxElementVolume = _hypMaxElementVolume->GetMaxVolume(); - isOk =true; - aStatus = SMESH_Hypothesis::HYP_OK; + if ( !_hypMaxElementVolume ) + _hypMaxElementVolume = dynamic_cast< const StdMeshers_MaxElementVolume*> ( *h ); + if ( !_viscousLayersHyp ) + _viscousLayersHyp = dynamic_cast< const StdMeshers_ViscousLayers*> ( *h ); + + if ( *h != _hypMaxElementVolume && + *h != _viscousLayersHyp ) + aStatus = HYP_INCOMPATIBLE; } - else - aStatus = SMESH_Hypothesis::HYP_INCOMPATIBLE; - return isOk; + if ( _hypMaxElementVolume ) + _maxElementVolume = _hypMaxElementVolume->GetMaxVolume(); + + return aStatus == HYP_OK; } //============================================================================= @@ -187,9 +185,6 @@ bool NETGENPlugin_NETGEN_3D::Compute(SMESH_Mesh& aMesh, { const int invalid_ID = -1; - SMESH::Controls::Area areaControl; - SMESH::Controls::TSequenceOfXYZ nodesCoords; - // maps nodes to ng ID typedef map< const SMDS_MeshNode*, int, TIDCompare > TNodeToIDMap; typedef TNodeToIDMap::value_type TN2ID; @@ -205,9 +200,19 @@ bool NETGENPlugin_NETGEN_3D::Compute(SMESH_Mesh& aMesh, TopAbs_ShapeEnum mainType = aMesh.GetShapeToMesh().ShapeType(); bool checkReverse = ( mainType == TopAbs_COMPOUND || mainType == TopAbs_COMPSOLID ); - StdMeshers_QuadToTriaAdaptor Adaptor; + StdMeshers_ProxyMesh::Ptr proxyMesh( new StdMeshers_ProxyMesh( aMesh )); + if ( _viscousLayersHyp ) + { + proxyMesh = _viscousLayersHyp->Compute( aMesh, aShape ); + if ( !proxyMesh ) + return false; + } if ( aMesh.NbQuadrangles() > 0 ) - Adaptor.Compute(aMesh,aShape); + { + StdMeshers_QuadToTriaAdaptor* Adaptor = new StdMeshers_QuadToTriaAdaptor; + Adaptor->Compute(aMesh,aShape,proxyMesh.get()); + proxyMesh.reset( Adaptor ); + } for ( TopExp_Explorer exFa( aShape, TopAbs_FACE ); exFa.More(); exFa.Next()) { @@ -221,7 +226,7 @@ bool NETGENPlugin_NETGEN_3D::Compute(SMESH_Mesh& aMesh, // so we use it as less as possible isRev = SMESH_Algo::IsReversedSubMesh( TopoDS::Face(aShapeFace), meshDS ); - const SMESHDS_SubMesh * aSubMeshDSFace = meshDS->MeshElements( aShapeFace ); + const SMESHDS_SubMesh * aSubMeshDSFace = proxyMesh->GetSubMesh( aShapeFace ); if ( !aSubMeshDSFace ) continue; SMDS_ElemIteratorPtr iteratorElem = aSubMeshDSFace->GetElements(); while ( iteratorElem->more() ) // loop on elements on a geom face @@ -230,72 +235,49 @@ bool NETGENPlugin_NETGEN_3D::Compute(SMESH_Mesh& aMesh, const SMDS_MeshElement* elem = iteratorElem->next(); if ( !elem ) return error( COMPERR_BAD_INPUT_MESH, "Null element encounters"); - vector< const SMDS_MeshElement* > trias; - bool isTraingle = ( elem->NbNodes() == ( elem->IsQuadratic() ? 6 : 3 )); - if ( !isTraingle ) - { - // use adaptor to convert quadrangle face into triangles - const list* faces = Adaptor.GetTriangles(elem); - if(faces==0) - return error( COMPERR_BAD_INPUT_MESH, - SMESH_Comment("No triangles in adaptor for element ")<GetID()); - trias.assign( faces->begin(), faces->end() ); - } - else - { - trias.push_back( elem ); - } + if ( elem->NbCornerNodes() != 3 ) + return error( COMPERR_BAD_INPUT_MESH, "Not triangle element encounters"); + // Add nodes of triangles and triangles them-selves to netgen mesh - for ( int i = 0; i < trias.size(); ++i ) + // 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 = elem->GetNode( iN ); + const int shapeID = node->GetPosition()->GetShapeId(); + if ( node->GetPosition()->GetTypeOfPosition() == SMDS_TOP_EDGE && + helper.IsDegenShape( shapeID )) { - 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 ) - { - 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[ isRev ? 2-iN : iN ] = ngID; + // 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; } - // 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 ) + int& ngID = nodeToNetgenID.insert(TN2ID( node, invalid_ID )).first->second; + if ( ngID == invalid_ID ) { - swap( Netgen_triangle[1], Netgen_triangle[2] ); - Ng_AddSurfaceElement(Netgen_mesh, NG_TRIG, Netgen_triangle); + 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[ isRev ? 2-iN : iN ] = ngID; } -#ifdef _DEBUG_ - // check if a trainge is degenerated - areaControl.GetPoints( elem, nodesCoords ); - double area = areaControl.GetValue( nodesCoords ); - if ( area <= DBL_MIN ) { - MESSAGE( "Warning: Degenerated " << elem ); + // 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 && !proxyMesh->IsTemporary( elem )) + { + swap( Netgen_triangle[1], Netgen_triangle[2] ); + Ng_AddSurfaceElement(Netgen_mesh, NG_TRIG, Netgen_triangle); } -#endif } // loop on elements on a face } // loop on faces of a SOLID or SHELL @@ -409,54 +391,14 @@ bool NETGENPlugin_NETGEN_3D::Compute(SMESH_Mesh& aMesh, MESSAGE("NETGENPlugin_NETGEN_3D::Compute with maxElmentsize = " << _maxElementVolume); const int invalid_ID = -1; bool _quadraticMesh = false; - typedef map< const SMDS_MeshNode*, int, TIDCompare > TNodeToIDMap; - TNodeToIDMap nodeToNetgenID; - list< const SMDS_MeshElement* > triangles; - SMESHDS_Mesh* MeshDS = aHelper->GetMeshDS(); SMESH_MesherHelper::MType MeshType = aHelper->IsQuadraticMesh(); - + if(MeshType == SMESH_MesherHelper::COMP) return error( COMPERR_BAD_INPUT_MESH, SMESH_Comment("Mesh with linear and quadratic elements given.")); else if (MeshType == SMESH_MesherHelper::QUADRATIC) _quadraticMesh = true; - - StdMeshers_QuadToTriaAdaptor Adaptor; - if ( aMesh.NbQuadrangles() > 0 ) - Adaptor.Compute(aMesh); - - SMDS_FaceIteratorPtr fIt = MeshDS->facesIterator(/*idInceasingOrder=*/true); - while( fIt->more()) - { - // check element - const SMDS_MeshElement* elem = fIt->next(); - if ( !elem ) - return error( COMPERR_BAD_INPUT_MESH, "Null element encounters"); - - vector< const SMDS_MeshElement* > trias; - bool isTraingle = ( elem->NbCornerNodes() == 3 ); - if ( !isTraingle ) { - // using adaptor - const list* faces = Adaptor.GetTriangles(elem); - if(faces==0) - continue; // Issue 0020682. There already can be 3d mesh - trias.assign( faces->begin(), faces->end() ); - } - else { - trias.push_back( elem ); - } - for ( int i = 0; i < trias.size(); ++i ) - { - triangles.push_back( trias[i] ); - for ( int iN = 0; iN < 3; ++iN ) - { - const SMDS_MeshNode* node = trias[i]->GetNode( iN ); - // put elem nodes to nodeToNetgenID map - nodeToNetgenID.insert( make_pair( node, invalid_ID )); - } - } - } // --------------------------------- // Feed the Netgen with surface mesh @@ -474,41 +416,53 @@ bool NETGENPlugin_NETGEN_3D::Compute(SMESH_Mesh& aMesh, NETGENPlugin_NetgenLibWrapper ngLib; Ng_Mesh * Netgen_mesh = ngLib._ngMesh; - // set nodes and remember thier netgen IDs - - TNodeToIDMap::iterator n_id = nodeToNetgenID.begin(); - for ( ; n_id != nodeToNetgenID.end(); ++n_id ) + StdMeshers_ProxyMesh::Ptr proxyMesh( new StdMeshers_ProxyMesh( aMesh )); + if ( aMesh.NbQuadrangles() > 0 ) { - const SMDS_MeshNode* node = n_id->first; - - Netgen_point [ 0 ] = node->X(); - Netgen_point [ 1 ] = node->Y(); - Netgen_point [ 2 ] = node->Z(); - Ng_AddPoint(Netgen_mesh, Netgen_point); - n_id->second = ++Netgen_NbOfNodes; // set netgen ID + StdMeshers_QuadToTriaAdaptor* Adaptor = new StdMeshers_QuadToTriaAdaptor; + Adaptor->Compute(aMesh); + proxyMesh.reset( Adaptor ); } - // set triangles - list< const SMDS_MeshElement* >::iterator tria = triangles.begin(); - for ( ; tria != triangles.end(); ++tria) + // maps nodes to ng ID + typedef map< const SMDS_MeshNode*, int, TIDCompare > TNodeToIDMap; + typedef TNodeToIDMap::value_type TN2ID; + TNodeToIDMap nodeToNetgenID; + + SMDS_ElemIteratorPtr fIt = proxyMesh->GetFaces(); + while( fIt->more()) { - int i = 0; - SMDS_ElemIteratorPtr triangleNodesIt = (*tria)->nodesIterator(); - while ( triangleNodesIt->more() ) { - const SMDS_MeshNode * node = - static_cast(triangleNodesIt->next()); - if(aHelper->IsMedium(node)) - continue; - Netgen_triangle[ i ] = nodeToNetgenID[ node ]; - ++i; + // check element + const SMDS_MeshElement* elem = fIt->next(); + if ( !elem ) + return error( COMPERR_BAD_INPUT_MESH, "Null element encounters"); + if ( elem->NbCornerNodes() != 3 ) + return error( COMPERR_BAD_INPUT_MESH, "Not triangle element encounters"); + + // add three nodes of triangle + for ( int iN = 0; iN < 3; ++iN ) + { + const SMDS_MeshNode* node = elem->GetNode( iN ); + int& ngID = nodeToNetgenID.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; } Ng_AddSurfaceElement(Netgen_mesh, NG_TRIG, Netgen_triangle); } + proxyMesh.reset(); // delete tmp faces // vector of nodes in which node index == netgen ID vector< const SMDS_MeshNode* > nodeVec ( nodeToNetgenID.size() + 1 ); // insert old nodes into nodeVec - for ( n_id = nodeToNetgenID.begin(); n_id != nodeToNetgenID.end(); ++n_id ) + TNodeToIDMap::iterator n_id = nodeToNetgenID.begin(); + for ( ; n_id != nodeToNetgenID.end(); ++n_id ) nodeVec.at( n_id->second ) = n_id->first; nodeToNetgenID.clear(); diff --git a/src/NETGENPlugin/NETGENPlugin_NETGEN_3D.hxx b/src/NETGENPlugin/NETGENPlugin_NETGEN_3D.hxx index 520cb44..39cd271 100644 --- a/src/NETGENPlugin/NETGENPlugin_NETGEN_3D.hxx +++ b/src/NETGENPlugin/NETGENPlugin_NETGEN_3D.hxx @@ -26,7 +26,6 @@ // Created : lundi 27 Janvier 2003 // Author : Nadir BOUHAMOU (CEA) // Project : SALOME -// $Header$ //============================================================================= // #ifndef _NETGENPlugin_NETGEN_3D_HXX_ @@ -35,10 +34,11 @@ #include "NETGENPlugin_Defs.hxx" #include "SMESH_3D_Algo.hxx" -#include "SMESH_Mesh.hxx" -#include "StdMeshers_MaxElementVolume.hxx" #include "Utils_SALOME_Exception.hxx" +class StdMeshers_ViscousLayers; +class StdMeshers_MaxElementVolume; + class NETGENPLUGIN_EXPORT NETGENPlugin_NETGEN_3D: public SMESH_3D_Algo { public: @@ -63,6 +63,7 @@ protected: double _maxElementVolume; const StdMeshers_MaxElementVolume* _hypMaxElementVolume; + const StdMeshers_ViscousLayers* _viscousLayersHyp; }; #endif -- 2.39.2