X-Git-Url: http://git.salome-platform.org/gitweb/?a=blobdiff_plain;f=src%2FNETGENPlugin%2FNETGENPlugin_NETGEN_2D_ONLY.cxx;h=acb11d6dae5ee78158b471a0625586332c891041;hb=8b4cb29237b43e0050f2bd4702830889cbb9f048;hp=fb48a06edfffe95d0392c53fb8d58b1ca9c38b11;hpb=5e1cfe70db0966728aff50c2d2796bf13bb3614f;p=plugins%2Fnetgenplugin.git diff --git a/src/NETGENPlugin/NETGENPlugin_NETGEN_2D_ONLY.cxx b/src/NETGENPlugin/NETGENPlugin_NETGEN_2D_ONLY.cxx index fb48a06..acb11d6 100644 --- a/src/NETGENPlugin/NETGENPlugin_NETGEN_2D_ONLY.cxx +++ b/src/NETGENPlugin/NETGENPlugin_NETGEN_2D_ONLY.cxx @@ -1,24 +1,22 @@ -// Copyright (C) 2007-2008 CEA/DEN, EDF R&D, OPEN CASCADE +// Copyright (C) 2007-2012 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 +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. // -// This library is free software; you can redistribute it and/or -// modify it under the terms of the GNU Lesser General Public -// License as published by the Free Software Foundation; either -// version 2.1 of the License. +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. // -// This library is distributed in the hope that it will be useful, -// but WITHOUT ANY WARRANTY; without even the implied warranty of -// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -// Lesser General Public License for more details. +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA // -// You should have received a copy of the GNU Lesser General Public -// License along with this library; if not, write to the Free Software -// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA -// -// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com // + // File : NETGENPlugin_NETGEN_2D_ONLY.cxx // Author : Edward AGAPOV (OCC) // Project : SALOME @@ -26,6 +24,7 @@ #include "NETGENPlugin_NETGEN_2D_ONLY.hxx" #include "NETGENPlugin_Mesher.hxx" +#include "NETGENPlugin_Hypothesis_2D.hxx" #include "SMDS_MeshElement.hxx" #include "SMDS_MeshNode.hxx" @@ -39,6 +38,7 @@ #include "StdMeshers_LengthFromEdges.hxx" #include "StdMeshers_QuadranglePreference.hxx" +#include #include #include @@ -46,6 +46,7 @@ #include #include +#include /* Netgen include files @@ -53,19 +54,23 @@ namespace nglib { #include } +#ifndef OCCGEOMETRY #define OCCGEOMETRY +#endif #include #include //#include namespace netgen { extern int OCCGenerateMesh (OCCGeometry&, Mesh*&, int, int, char*); - /*extern*/ MeshingParameters mparam; + extern MeshingParameters mparam; } using namespace std; using namespace netgen; using namespace nglib; +//#define DUMP_SEGMENTS + //============================================================================= /*! * @@ -78,16 +83,18 @@ NETGENPlugin_NETGEN_2D_ONLY::NETGENPlugin_NETGEN_2D_ONLY(int hypId, int studyId, { MESSAGE("NETGENPlugin_NETGEN_2D_ONLY::NETGENPlugin_NETGEN_2D_ONLY"); _name = "NETGEN_2D_ONLY"; - + _shapeType = (1 << TopAbs_FACE);// 1 bit /shape type _compatibleHypothesis.push_back("MaxElementArea"); _compatibleHypothesis.push_back("LengthFromEdges"); _compatibleHypothesis.push_back("QuadranglePreference"); + _compatibleHypothesis.push_back("NETGEN_Parameters_2D"); _hypMaxElementArea = 0; _hypLengthFromEdges = 0; _hypQuadranglePreference = 0; + _hypParameters = 0; } //============================================================================= @@ -138,21 +145,21 @@ bool NETGENPlugin_NETGEN_2D_ONLY::CheckHypothesis (SMESH_Mesh& aMesh, _hypLengthFromEdges = static_cast (hyp); else if ( hypName == "QuadranglePreference" ) _hypQuadranglePreference = static_cast(hyp); + else if ( hypName == "NETGEN_Parameters_2D" ) + _hypParameters = static_cast(hyp); else { aStatus = HYP_INCOMPATIBLE; return false; } } - if ( _hypMaxElementArea && _hypLengthFromEdges ) { + int nbHyps = bool(_hypMaxElementArea) + bool(_hypLengthFromEdges) + bool(_hypParameters ); + if ( nbHyps > 1 ) aStatus = HYP_CONCURENT; - return false; - } - - if ( _hypMaxElementArea || _hypLengthFromEdges ) + else if ( nbHyps == 1) aStatus = HYP_OK; - return aStatus == HYP_OK; + return ( aStatus == HYP_OK ); } //================================================================================ @@ -162,7 +169,7 @@ bool NETGENPlugin_NETGEN_2D_ONLY::CheckHypothesis (SMESH_Mesh& aMesh, */ //================================================================================ -static TError AddSegmentsToMesh(netgen::Mesh& ngMesh, +static TError addSegmentsToMesh(netgen::Mesh& ngMesh, OCCGeometry& geom, const TSideVector& wires, SMESH_MesherHelper& helper, @@ -172,20 +179,27 @@ static TError AddSegmentsToMesh(netgen::Mesh& ngMesh, // Check wires and count nodes // ---------------------------- int nbNodes = 0; + double totalLength = 0; for ( int iW = 0; iW < wires.size(); ++iW ) { StdMeshers_FaceSidePtr wire = wires[ iW ]; if ( wire->MissVertexNode() ) - return TError - (new SMESH_ComputeError(COMPERR_BAD_INPUT_MESH, "Missing nodes on vertices")); - + { + // Commented for issue 0020960. It worked for the case, let's wait for case where it doesn't. + // It seems that there is no reason for this limitation +// return TError +// (new SMESH_ComputeError(COMPERR_BAD_INPUT_MESH, "Missing nodes on vertices")); + if (getenv("USER") && string("eap")==getenv("USER")) + cout << "Warning: NETGENPlugin_NETGEN_2D_ONLY : try to work with missing nodes on vertices"<& uvPtVec = wire->GetUVPtStruct(); if ( uvPtVec.size() != wire->NbPoints() ) return TError (new SMESH_ComputeError(COMPERR_BAD_INPUT_MESH, SMESH_Comment("Unexpected nb of points on wire ") << iW << ": " << uvPtVec.size()<<" != "<NbPoints())); - nbNodes += wire->NbSegments(); + nbNodes += wire->NbPoints(); + totalLength += wire->Length(); } nodeVec.reserve( nbNodes ); @@ -197,37 +211,64 @@ static TError AddSegmentsToMesh(netgen::Mesh& ngMesh, // bb.Increase (bb.Diam()/10); // ngMesh.SetLocalH (bb.PMin(), bb.PMax(), 0.5); // set grading + // map for nodes on vertices since they can be shared between wires + // ( issue 0020676, face_int_box.brep) + map node2ngID; + const int faceID = 1, solidID = 0; - ngMesh.AddFaceDescriptor (FaceDescriptor(faceID, solidID, solidID, 0)); + if ( ngMesh.GetNFD() < 1 ) + ngMesh.AddFaceDescriptor (FaceDescriptor(faceID, solidID, solidID, 0)); for ( int iW = 0; iW < wires.size(); ++iW ) { StdMeshers_FaceSidePtr wire = wires[ iW ]; const vector& uvPtVec = wire->GetUVPtStruct(); + const int nbSegments = wire->NbPoints() - 1; + + // compute length of every segment + vector segLen( nbSegments ); + for ( int i = 0; i < nbSegments; ++i ) + segLen[i] = SMESH_TNodeXYZ( uvPtVec[ i ].node ).Distance( uvPtVec[ i+1 ].node ); - int firstPointID = ngMesh.GetNP() + 1; int edgeID = 1, posID = -2; - for ( int i = 0; i < wire->NbSegments(); ++i ) // loop on segments + bool isInternalWire = false; + for ( int i = 0; i < nbSegments; ++i ) // loop on segments { // Add the first point of a segment const SMDS_MeshNode * n = uvPtVec[ i ].node; - const int posShapeID = n->GetPosition()->GetShapeId(); + const int posShapeID = n->getshapeId(); + bool onVertex = ( n->GetPosition()->GetTypeOfPosition() == SMDS_TOP_VERTEX ); // skip nodes on degenerated edges if ( helper.IsDegenShape( posShapeID ) && - helper.IsDegenShape( uvPtVec[ i+1 ].node->GetPosition()->GetShapeId() )) + helper.IsDegenShape( uvPtVec[ i+1 ].node->getshapeId() )) continue; - nodeVec.push_back( n ); - - MeshPoint mp( Point<3> (n->X(), n->Y(), n->Z()) ); - ngMesh.AddPoint ( mp, 1, EDGEPOINT ); + int ngID1 = ngMesh.GetNP() + 1, ngID2 = ngID1+1; + if ( onVertex ) + ngID1 = node2ngID.insert( make_pair( n, ngID1 )).first->second; + if ( ngID1 > ngMesh.GetNP() ) + { + MeshPoint mp( Point<3> (n->X(), n->Y(), n->Z()) ); + ngMesh.AddPoint ( mp, 1, EDGEPOINT ); + nodeVec.push_back( n ); + } + else + { + ngID2 = ngMesh.GetNP() + 1; + if ( i > 0 ) // prev segment belongs to same wire + { + Segment& prevSeg = ngMesh.LineSegment( ngMesh.GetNSeg() ); + prevSeg[1] = ngID1; + } + } // Add the segment Segment seg; - seg.p1 = ngMesh.GetNP(); // ng node id - seg.p2 = seg.p1 + 1; // ng node id + seg[0] = ngID1; // ng node id + seg[1] = ngID2; // ng node id + seg.edgenr = ngMesh.GetNSeg() + 1;// segment id seg.si = faceID; // = geom.fmap.FindIndex (face); @@ -240,7 +281,7 @@ static TError AddSegmentsToMesh(netgen::Mesh& ngMesh, seg.epgeominfo[ iEnd ].v = pnt.v; // find out edge id and node parameter on edge - bool onVertex = ( pnt.node->GetPosition()->GetTypeOfPosition() == SMDS_TOP_VERTEX ); + onVertex = ( pnt.node->GetPosition()->GetTypeOfPosition() == SMDS_TOP_VERTEX ); if ( onVertex || posShapeID != posID ) { // get edge id @@ -250,29 +291,74 @@ static TError AddSegmentsToMesh(netgen::Mesh& ngMesh, const TopoDS_Edge& edge = wire->Edge( wire->EdgeIndex( normParam )); edgeID = geom.emap.FindIndex( edge ); posID = posShapeID; - if ( onVertex ) // param on curve is different on each of two edges - seg.epgeominfo[ iEnd ].dist = helper.GetNodeU( edge, pnt.node ); + isInternalWire = ( edge.Orientation() == TopAbs_INTERNAL ); + // if ( onVertex ) // param on curve is different on each of two edges + // seg.epgeominfo[ iEnd ].dist = helper.GetNodeU( edge, pnt.node ); } seg.epgeominfo[ iEnd ].edgenr = edgeID; // = geom.emap.FindIndex(edge); } ngMesh.AddSegment (seg); + { + // restrict size of elements near the segment + SMESH_TNodeXYZ np1( n ), np2( uvPtVec[ i+1 ].node ); + // get an average size of adjacent segments to avoid sharp change of + // element size (regression on issue 0020452, note 0010898) + int iPrev = SMESH_MesherHelper::WrapIndex( i-1, nbSegments ); + int iNext = SMESH_MesherHelper::WrapIndex( i+1, nbSegments ); + double avgH = ( segLen[ iPrev ] + segLen[ i ] + segLen[ iNext ]) / 3; + + NETGENPlugin_Mesher::RestrictLocalSize( ngMesh, 0.5*(np1+np2), avgH ); + } +#ifdef DUMP_SEGMENTS + cout << "Segment: " << seg.edgenr << endl + << "\tp1: " << seg[0] << endl + << "\tp2: " << seg[1] << endl + << "\tp0 param: " << seg.epgeominfo[ 0 ].dist << endl + << "\tp0 uv: " << seg.epgeominfo[ 0 ].u <<", "<< seg.epgeominfo[ 0 ].v << endl + << "\tp0 edge: " << seg.epgeominfo[ 0 ].edgenr << endl + << "\tp1 param: " << seg.epgeominfo[ 1 ].dist << endl + << "\tp1 uv: " << seg.epgeominfo[ 1 ].u <<", "<< seg.epgeominfo[ 1 ].v << endl + << "\tp1 edge: " << seg.epgeominfo[ 1 ].edgenr << endl; +#endif + if ( isInternalWire ) + { + swap (seg[0], seg[1]); + swap( seg.epgeominfo[0], seg.epgeominfo[1] ); + seg.edgenr = ngMesh.GetNSeg() + 1; // segment id + ngMesh.AddSegment (seg); +#ifdef DUMP_SEGMENTS + cout << "Segment: " << seg.edgenr << endl << "\tis REVRESE of the previous one" << endl; +#endif + } + } // loop on segments on a wire -// cout << "Segment: " << seg.edgenr << endl -// << "\tp1: " << seg.p1 << endl -// << "\tp2: " << seg.p2 << endl -// << "\tp0 param: " << seg.epgeominfo[ 0 ].dist << endl -// << "\tp0 uv: " << seg.epgeominfo[ 0 ].u <<", "<< seg.epgeominfo[ 0 ].v << endl -// << "\tp0 edge: " << seg.epgeominfo[ 0 ].edgenr << endl -// << "\tp1 param: " << seg.epgeominfo[ 1 ].dist << endl -// << "\tp1 uv: " << seg.epgeominfo[ 1 ].u <<", "<< seg.epgeominfo[ 1 ].v << endl -// << "\tp1 edge: " << seg.epgeominfo[ 1 ].edgenr << endl; + // close chain of segments + if ( nbSegments > 0 ) + { + Segment& lastSeg = ngMesh.LineSegment( ngMesh.GetNSeg() - int( isInternalWire)); + const SMDS_MeshNode * lastNode = uvPtVec.back().node; + lastSeg[1] = node2ngID.insert( make_pair( lastNode, lastSeg[1] )).first->second; + if ( lastSeg[1] > ngMesh.GetNP() ) + { + MeshPoint mp( Point<3> (lastNode->X(), lastNode->Y(), lastNode->Z()) ); + ngMesh.AddPoint ( mp, 1, EDGEPOINT ); + nodeVec.push_back( lastNode ); + } + if ( isInternalWire ) + { + Segment& realLastSeg = ngMesh.LineSegment( ngMesh.GetNSeg() ); + realLastSeg[0] = lastSeg[1]; + } } - Segment& seg = ngMesh.LineSegment( ngMesh.GetNSeg() ); - seg.p2 = firstPointID; - } - ngMesh.CalcSurfacesOfNode(); + } // loop on wires of a face + + // add a segment instead of internal vertex + NETGENPlugin_Internals intShapes( *helper.GetMesh(), helper.GetSubShape(), /*is3D=*/false ); + NETGENPlugin_Mesher::addIntVerticesInFaces( geom, ngMesh, nodeVec, intShapes ); + + ngMesh.CalcSurfacesOfNode(); return TError(); } @@ -286,6 +372,9 @@ static TError AddSegmentsToMesh(netgen::Mesh& ngMesh, bool NETGENPlugin_NETGEN_2D_ONLY::Compute(SMESH_Mesh& aMesh, const TopoDS_Shape& aShape) { +#ifdef WITH_SMESH_CANCEL_COMPUTE + netgen::multithread.terminate = 0; +#endif MESSAGE("NETGENPlugin_NETGEN_2D_ONLY::Compute()"); SMESHDS_Mesh* meshDS = aMesh.GetMeshDS(); @@ -311,60 +400,73 @@ bool NETGENPlugin_NETGEN_2D_ONLY::Compute(SMESH_Mesh& aMesh, return error(COMPERR_BAD_INPUT_MESH, SMESH_Comment("Too few segments: ")<NbSegments()); - // ------------------------- - // Make input netgen mesh - // ------------------------- - - Ng_Init(); - netgen::Mesh * ngMesh = new netgen::Mesh (); + // -------------------- + // compute edge length + // -------------------- + NETGENPlugin_Mesher aMesher( &aMesh, aShape, /*isVolume=*/false); netgen::OCCGeometry occgeo; - NETGENPlugin_Mesher::PrepareOCCgeometry( occgeo, F, aMesh ); + aMesher.PrepareOCCgeometry( occgeo, F, aMesh ); occgeo.fmap.Clear(); // face can be reversed, which is wrong in this case (issue 19978) occgeo.fmap.Add( F ); - vector< const SMDS_MeshNode* > nodeVec; - problem = AddSegmentsToMesh( *ngMesh, occgeo, wires, helper, nodeVec ); - if ( problem && !problem->IsOK() ) { - delete ngMesh; Ng_Exit(); - return error( problem ); + if ( _hypParameters ) + { + aMesher.SetParameters(_hypParameters); } - - // -------------------- - // compute edge length - // -------------------- - - double edgeLength = 0; - if (_hypLengthFromEdges || !_hypLengthFromEdges && !_hypMaxElementArea) + else { - int nbSegments = 0; - for ( int iW = 0; iW < nbWires; ++iW ) + double edgeLength = 0; + if (_hypLengthFromEdges || (!_hypLengthFromEdges && !_hypMaxElementArea)) { - edgeLength += wires[ iW ]->Length(); - nbSegments += wires[ iW ]->NbSegments(); + int nbSegments = 0; + for ( int iW = 0; iW < nbWires; ++iW ) + { + edgeLength += wires[ iW ]->Length(); + nbSegments += wires[ iW ]->NbSegments(); + } + if ( nbSegments ) + edgeLength /= nbSegments; } - if ( nbSegments ) - edgeLength /= nbSegments; - } - if ( _hypMaxElementArea ) - { - double maxArea = _hypMaxElementArea->GetMaxArea(); - edgeLength = sqrt(2. * maxArea/sqrt(3.0)); + if ( _hypMaxElementArea ) + { + double maxArea = _hypMaxElementArea->GetMaxArea(); + edgeLength = sqrt(2. * maxArea/sqrt(3.0)); + } + if ( edgeLength < DBL_MIN ) + edgeLength = occgeo.GetBoundingBox().Diam(); + + netgen::mparam.maxh = edgeLength; + netgen::mparam.minh = aMesher.GetDefaultMinSize( aShape, netgen::mparam.maxh ); + netgen::mparam.quad = _hypQuadranglePreference ? 1 : 0; + netgen::mparam.grading = 0.7; // very coarse mesh by default } - if ( edgeLength < DBL_MIN ) - edgeLength = occgeo.GetBoundingBox().Diam(); +#ifdef NETGEN_NEW + occgeo.face_maxh = netgen::mparam.maxh; +#endif - //cout << " edgeLength = " << edgeLength << endl; + // ------------------------- + // Make input netgen mesh + // ------------------------- + + NETGENPlugin_NetgenLibWrapper ngLib; + netgen::Mesh * ngMesh = (netgen::Mesh*) ngLib._ngMesh; - netgen::mparam.maxh = edgeLength; - netgen::mparam.quad = _hypQuadranglePreference ? 1 : 0; - //ngMesh->SetGlobalH ( edgeLength ); + Box<3> bb = occgeo.GetBoundingBox(); + bb.Increase (bb.Diam()/10); + ngMesh->SetLocalH (bb.PMin(), bb.PMax(), netgen::mparam.grading); + ngMesh->SetGlobalH (netgen::mparam.maxh); + + vector< const SMDS_MeshNode* > nodeVec; + problem = addSegmentsToMesh( *ngMesh, occgeo, wires, helper, nodeVec ); + if ( problem && !problem->IsOK() ) + return error( problem ); // ------------------------- // Generate surface mesh // ------------------------- - char *optstr; + char *optstr = 0; int startWith = MESHCONST_MESHSURFACE; int endWith = MESHCONST_OPTSURFACE; int err = 1; @@ -374,20 +476,26 @@ bool NETGENPlugin_NETGEN_2D_ONLY::Compute(SMESH_Mesh& aMesh, OCC_CATCH_SIGNALS; #endif err = netgen::OCCGenerateMesh(occgeo, ngMesh, startWith, endWith, optstr); +#ifdef WITH_SMESH_CANCEL_COMPUTE + if(netgen::multithread.terminate) + return false; +#endif + if ( err ) + error(SMESH_Comment("Error in netgen::OCCGenerateMesh() at ") << netgen::multithread.task); } - catch (Standard_Failure& ex) { - string comment = ex.DynamicType()->Name(); - if ( ex.GetMessageString() && strlen( ex.GetMessageString() )) { - comment += ": "; - comment += ex.GetMessageString(); - } - error(COMPERR_OCC_EXCEPTION, comment); - } - catch (NgException exc) { - error( SMESH_Comment("NgException: ") << exc.What() ); + catch (Standard_Failure& ex) + { + SMESH_Comment str("Exception in netgen::OCCGenerateMesh()"); + str << " at " << netgen::multithread.task + << ": " << ex.DynamicType()->Name(); + if ( ex.GetMessageString() && strlen( ex.GetMessageString() )) + str << ": " << ex.GetMessageString(); + error(str); } catch (...) { - error(COMPERR_EXCEPTION,"Exception in netgen::OCCGenerateMesh()"); + SMESH_Comment str("Exception in netgen::OCCGenerateMesh()"); + str << " at " << netgen::multithread.task; + error(str); } // ---------------------------------------------------- @@ -404,19 +512,26 @@ bool NETGENPlugin_NETGEN_2D_ONLY::Compute(SMESH_Mesh& aMesh, for ( int i = nbInputNodes + 1; i <= nbNodes; ++i ) { const MeshPoint& ngPoint = ngMesh->Point(i); +#ifdef NETGEN_NEW + SMDS_MeshNode * node = meshDS->AddNode(ngPoint(0), ngPoint(1), ngPoint(2)); +#else SMDS_MeshNode * node = meshDS->AddNode(ngPoint.X(), ngPoint.Y(), ngPoint.Z()); +#endif nodeVec[ i-1 ] = node; } // create faces bool reverse = ( aShape.Orientation() == TopAbs_REVERSED ); - for ( int i = 1; i <= nbFaces ; ++i ) + int i,j; + for ( i = 1; i <= nbFaces ; ++i ) { const Element2d& elem = ngMesh->SurfaceElement(i); vector nodes( elem.GetNP() ); - for (int j=1; j <= elem.GetNP(); ++j) + for (j=1; j <= elem.GetNP(); ++j) { int pind = elem.PNum(j); + if ( pind-1 < 0 ) + break; const SMDS_MeshNode* node = nodeVec.at(pind-1); if ( reverse ) nodes[ nodes.size()-j ] = node; @@ -428,17 +543,108 @@ bool NETGENPlugin_NETGEN_2D_ONLY::Compute(SMESH_Mesh& aMesh, meshDS->SetNodeOnFace((SMDS_MeshNode*)node, faceID, pgi.u, pgi.v); } } - SMDS_MeshFace* face = 0; - if ( elem.GetType() == TRIG ) - face = helper.AddFace(nodes[0],nodes[1],nodes[2]); - else - face = helper.AddFace(nodes[0],nodes[1],nodes[2],nodes[3]); + if ( j > elem.GetNP() ) + { + SMDS_MeshFace* face = 0; + if ( elem.GetType() == TRIG ) + face = helper.AddFace(nodes[0],nodes[1],nodes[2]); + else + face = helper.AddFace(nodes[0],nodes[1],nodes[2],nodes[3]); + } + } + + return !err; +} + +#ifdef WITH_SMESH_CANCEL_COMPUTE +void NETGENPlugin_NETGEN_2D_ONLY::CancelCompute() +{ + SMESH_Algo::CancelCompute(); + netgen::multithread.terminate = 1; +} +#endif + +//============================================================================= +/*! + * + */ +//============================================================================= + +bool NETGENPlugin_NETGEN_2D_ONLY::Evaluate(SMESH_Mesh& aMesh, + const TopoDS_Shape& aShape, + MapShapeNbElems& aResMap) +{ + TopoDS_Face F = TopoDS::Face(aShape); + if(F.IsNull()) + return false; + + // collect info from edges + int nb0d = 0, nb1d = 0; + bool IsQuadratic = false; + bool IsFirst = true; + double fullLen = 0.0; + TopTools_MapOfShape tmpMap; + for (TopExp_Explorer exp(F, TopAbs_EDGE); exp.More(); exp.Next()) { + TopoDS_Edge E = TopoDS::Edge(exp.Current()); + if( tmpMap.Contains(E) ) + continue; + tmpMap.Add(E); + SMESH_subMesh *aSubMesh = aMesh.GetSubMesh(exp.Current()); + MapShapeNbElemsItr anIt = aResMap.find(aSubMesh); + if( anIt==aResMap.end() ) { + SMESH_subMesh *sm = aMesh.GetSubMesh(F); + SMESH_ComputeErrorPtr& smError = sm->GetComputeError(); + smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,"Submesh can not be evaluated",this)); + return false; + } + std::vector aVec = (*anIt).second; + nb0d += aVec[SMDSEntity_Node]; + nb1d += Max(aVec[SMDSEntity_Edge],aVec[SMDSEntity_Quad_Edge]); + double aLen = SMESH_Algo::EdgeLength(E); + fullLen += aLen; + if(IsFirst) { + IsQuadratic = (aVec[SMDSEntity_Quad_Edge] > aVec[SMDSEntity_Edge]); + IsFirst = false; + } } + tmpMap.Clear(); - Ng_DeleteMesh((nglib::Ng_Mesh*)ngMesh); - Ng_Exit(); + // compute edge length + double ELen = 0; + if (_hypLengthFromEdges || !_hypLengthFromEdges && !_hypMaxElementArea) { + if ( nb1d > 0 ) + ELen = fullLen / nb1d; + } + if ( _hypMaxElementArea ) { + double maxArea = _hypMaxElementArea->GetMaxArea(); + ELen = sqrt(2. * maxArea/sqrt(3.0)); + } + GProp_GProps G; + BRepGProp::SurfaceProperties(F,G); + double anArea = G.Mass(); - NETGENPlugin_Mesher::RemoveTmpFiles(); + const int hugeNb = numeric_limits::max()/10; + if ( anArea / hugeNb > ELen*ELen ) + { + SMESH_subMesh *sm = aMesh.GetSubMesh(F); + SMESH_ComputeErrorPtr& smError = sm->GetComputeError(); + smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,"Submesh can not be evaluated.\nToo small element length",this)); + return false; + } + int nbFaces = (int) ( anArea / ( ELen*ELen*sqrt(3.) / 4 ) ); + int nbNodes = (int) ( ( nbFaces*3 - (nb1d-1)*2 ) / 6 + 1 ); + std::vector aVec(SMDSEntity_Last); + for(int i=SMDSEntity_Node; i