X-Git-Url: http://git.salome-platform.org/gitweb/?p=modules%2Fsmesh.git;a=blobdiff_plain;f=src%2FStdMeshers%2FStdMeshers_RadialQuadrangle_1D2D.cxx;h=b861f3701c87d909fc01924f1445e02bc75b9e38;hp=cb6aa8f7fa83e2b978a51e5d1ddb72f8dedda9ea;hb=59627b07d70f4caa4c768be6805334d2610fa54c;hpb=a3ee66b0371b361bf97d6261c61a35995b4595bb

diff --git a/src/StdMeshers/StdMeshers_RadialQuadrangle_1D2D.cxx b/src/StdMeshers/StdMeshers_RadialQuadrangle_1D2D.cxx
index cb6aa8f7f..b861f3701 100644
--- a/src/StdMeshers/StdMeshers_RadialQuadrangle_1D2D.cxx
+++ b/src/StdMeshers/StdMeshers_RadialQuadrangle_1D2D.cxx
@@ -1,4 +1,4 @@
-// Copyright (C) 2007-2014  CEA/DEN, EDF R&D, OPEN CASCADE
+// Copyright (C) 2007-2015  CEA/DEN, EDF R&D, OPEN CASCADE
 //
 // This library is free software; you can redistribute it and/or
 // modify it under the terms of the GNU Lesser General Public
@@ -75,10 +75,12 @@ StdMeshers_RadialQuadrangle_1D2D::StdMeshers_RadialQuadrangle_1D2D(int hypId,
 
   _compatibleHypothesis.push_back("LayerDistribution2D");
   _compatibleHypothesis.push_back("NumberOfLayers2D");
-  myNbLayerHypo = 0;
-  myDistributionHypo = 0;
   _requireDiscreteBoundary = false;
   _supportSubmeshes = true;
+  _neededLowerHyps[ 1 ] = true;  // suppress warning on hiding a global 1D algo
+
+  myNbLayerHypo      = 0;
+  myDistributionHypo = 0;
 }
 
 
@@ -278,6 +280,38 @@ namespace
     }
     return nbe;
   }
+  //================================================================================
+  /*!
+   * \brief Checks if the common vertex between LinEdge's lies inside the circle
+   *  and not outside
+   *  \param [in] CircEdge - 
+   *  \param [in] LinEdge1 - 
+   *  \param [in] LinEdge2 - 
+   *  \return bool - false if there are 3 EDGEs and the corner is outside
+   */
+  //================================================================================
+
+  bool isCornerInsideCircle(const TopoDS_Edge& CircEdge,
+                            const TopoDS_Edge& LinEdge1,
+                            const TopoDS_Edge& LinEdge2)
+  {
+    if ( !CircEdge.IsNull() &&
+         !LinEdge1.IsNull() &&
+         !LinEdge2.IsNull() )
+    {
+      Handle(Geom_Circle) aCirc = Handle(Geom_Circle)::DownCast( getCurve( CircEdge ));
+      TopoDS_Vertex aCommonV;
+      if ( !aCirc.IsNull() &&
+           TopExp::CommonVertex( LinEdge1, LinEdge2, aCommonV ))
+      {
+        gp_Pnt aCommonP = BRep_Tool::Pnt( aCommonV );
+        gp_Pnt  aCenter = aCirc->Location();
+        double     dist = aCenter.Distance( aCommonP );
+        return dist < 0.1 * aCirc->Radius();
+      }
+    }
+    return true;
+  }
 
 //================================================================================
 //================================================================================
@@ -436,7 +470,7 @@ bool StdMeshers_RadialQuadrangle_1D2D::Compute(SMESH_Mesh&         aMesh,
   TopoDS_Edge CircEdge, LinEdge1, LinEdge2;
   int nbe = analyseFace( aShape, CircEdge, LinEdge1, LinEdge2 );
   Handle(Geom_Circle) aCirc = Handle(Geom_Circle)::DownCast( getCurve( CircEdge ));
-  if( nbe>3 || nbe < 1 || aCirc.IsNull() )
+  if( nbe > 3 || nbe < 1 || aCirc.IsNull() )
     return error("The face must be a full circle or a part of circle (i.e. the number "
                  "of edges is less or equal to 3 and one of them is a circle curve)");
 
@@ -503,8 +537,8 @@ bool StdMeshers_RadialQuadrangle_1D2D::Compute(SMESH_Mesh&         aMesh,
     gp_Vec2d aVec2d(PC,p2dV);
     Nodes1.resize( myLayerPositions.size()+1 );
     Nodes2.resize( myLayerPositions.size()+1 );
-    int i = 0;
-    for(; i<myLayerPositions.size(); i++) {
+    size_t i = 0;
+    for ( ; i < myLayerPositions.size(); i++ ) {
       gp_Pnt P( P0.X() + aVec.X()*myLayerPositions[i],
                 P0.Y() + aVec.Y()*myLayerPositions[i],
                 P0.Z() + aVec.Z()*myLayerPositions[i] );
@@ -614,7 +648,7 @@ bool StdMeshers_RadialQuadrangle_1D2D::Compute(SMESH_Mesh&         aMesh,
       meshDS->SetNodeOnEdge(NC, edgeID, cp);
       Nodes1.resize( myLayerPositions.size()+1 );
       Nodes2.resize( myLayerPositions.size()+1 );
-      int i = 0;
+      size_t i = 0;
       for(; i<myLayerPositions.size(); i++) {
         gp_Pnt P( P0.X() + aVec.X()*myLayerPositions[i],
                   P0.Y() + aVec.Y()*myLayerPositions[i],
@@ -669,15 +703,17 @@ bool StdMeshers_RadialQuadrangle_1D2D::Compute(SMESH_Mesh&         aMesh,
     // segments of line
     double fp, lp;
     Handle(Geom_Circle) aCirc = Handle(Geom_Circle)::DownCast( getCurve( CircEdge ));
-    Handle(Geom_Line)  aLine1 = Handle(Geom_Line)::DownCast( getCurve( LinEdge1 ));
-    Handle(Geom_Line)  aLine2 = Handle(Geom_Line)::DownCast( getCurve( LinEdge2 ));
-    if( aCirc.IsNull() || aLine1.IsNull() || aLine2.IsNull() )
+    Handle(Geom_Line)  aLine1 = Handle(Geom_Line  )::DownCast( getCurve( LinEdge1 ));
+    Handle(Geom_Line)  aLine2 = Handle(Geom_Line  )::DownCast( getCurve( LinEdge2 ));
+    if ( aCirc.IsNull() || aLine1.IsNull() || aLine2.IsNull() )
+      return error(COMPERR_BAD_SHAPE);
+    if ( !isCornerInsideCircle( CircEdge, LinEdge1, LinEdge2 ))
       return error(COMPERR_BAD_SHAPE);
 
     if ( !algo1d->ComputeCircularEdge( aMesh, CircEdge ))
       return error( algo1d->GetComputeError() );
     map< double, const SMDS_MeshNode* > theNodes;
-    if ( !GetSortedNodesOnEdge(aMesh.GetMeshDS(),CircEdge,true,theNodes))
+    if ( !GetSortedNodesOnEdge( aMesh.GetMeshDS(), CircEdge, true, theNodes ))
       return error("Circular edge is incorrectly meshed");
 
     myHelper->IsQuadraticSubMesh( aShape );
@@ -738,8 +774,8 @@ bool StdMeshers_RadialQuadrangle_1D2D::Compute(SMESH_Mesh&         aMesh,
       bool nodesFromP0ToP1 = ( theNodes.rbegin()->second == NF );
       NC = const_cast<SMDS_MeshNode*>
         ( nodesFromP0ToP1 ? theNodes.begin()->second : theNodes.rbegin()->second );
-      int i = 0, ir = Nodes1.size()-1;
-      int * pi = nodesFromP0ToP1 ? &i : &ir;
+      size_t i = 0, ir = Nodes1.size()-1;
+      size_t * pi = nodesFromP0ToP1 ? &i : &ir;
       itn = theNodes.begin();
       if ( nodesFromP0ToP1 ) ++itn;
       for ( ; i < Nodes1.size(); ++i, --ir, ++itn )
@@ -781,8 +817,8 @@ bool StdMeshers_RadialQuadrangle_1D2D::Compute(SMESH_Mesh&         aMesh,
         meshDS->SetNodeOnVertex(NC, vertID);
       }
       double dp = lp-fp;
-      int i = 0;
-      for(; i<myLayerPositions.size(); i++) {
+      size_t i = 0;
+      for ( ; i < myLayerPositions.size(); i++ ) {
         gp_Pnt P( P0.X() + aVec.X()*myLayerPositions[i],
                   P0.Y() + aVec.Y()*myLayerPositions[i],
                   P0.Z() + aVec.Z()*myLayerPositions[i] );
@@ -790,7 +826,7 @@ bool StdMeshers_RadialQuadrangle_1D2D::Compute(SMESH_Mesh&         aMesh,
         SMDS_MeshNode * node = meshDS->AddNode(P.X(), P.Y(), P.Z());
         Nodes1[i] = node;
         double param;
-        if(!ori)
+        if ( !ori )
           param = fp + dp*(1-myLayerPositions[i]);
         else
           param = fp + dp*myLayerPositions[i];
@@ -804,11 +840,11 @@ bool StdMeshers_RadialQuadrangle_1D2D::Compute(SMESH_Mesh&         aMesh,
       // create 1D elements on edge
       SMDS_MeshEdge* ME = myHelper->AddEdge( NC, Nodes1[0] );
       if(ME) meshDS->SetMeshElementOnShape(ME, edgeID);
-      for(i=1; i<Nodes1.size(); i++) {
+      for ( i = 1; i < Nodes1.size(); i++ ) {
         ME = myHelper->AddEdge( Nodes1[i-1], Nodes1[i] );
         if(ME) meshDS->SetMeshElementOnShape(ME, edgeID);
       }
-      if (nbe==2 && LinEdge1.Orientation() == TopAbs_INTERNAL )
+      if ( nbe == 2 && LinEdge1.Orientation() == TopAbs_INTERNAL )
         Nodes2 = Nodes1;
     }
     markEdgeAsComputedByMe( LinEdge1, aMesh.GetSubMesh( F ));
@@ -820,8 +856,8 @@ bool StdMeshers_RadialQuadrangle_1D2D::Compute(SMESH_Mesh&         aMesh,
         return error("Invalid mesh on a straight edge");
 
       bool nodesFromP0ToP2 = ( theNodes.rbegin()->second == NL );
-      int i = 0, ir = Nodes1.size()-1;
-      int * pi = nodesFromP0ToP2 ? &i : &ir;
+      size_t i = 0, ir = Nodes1.size()-1;
+      size_t * pi = nodesFromP0ToP2 ? &i : &ir;
       itn = theNodes.begin();
       if ( nodesFromP0ToP2 ) ++itn;
       for ( ; i < Nodes2.size(); ++i, --ir, ++itn )
@@ -850,7 +886,7 @@ bool StdMeshers_RadialQuadrangle_1D2D::Compute(SMESH_Mesh&         aMesh,
         PC = PL;
       }
       double dp = lp-fp;
-      for(int i=0; i<myLayerPositions.size(); i++) {
+      for ( size_t i = 0; i < myLayerPositions.size(); i++ ) {
         gp_Pnt P( P0.X() + aVec.X()*myLayerPositions[i],
                   P0.Y() + aVec.Y()*myLayerPositions[i],
                   P0.Z() + aVec.Z()*myLayerPositions[i] );
@@ -869,8 +905,8 @@ bool StdMeshers_RadialQuadrangle_1D2D::Compute(SMESH_Mesh&         aMesh,
       Nodes2[ myLayerPositions.size() ] = NL;
       // create 1D elements on edge
       SMDS_MeshEdge* ME = myHelper->AddEdge( NC, Nodes2[0] );
-      if(ME) meshDS->SetMeshElementOnShape(ME, edgeID);
-      for(int i=1; i<Nodes2.size(); i++) {
+      if ( ME ) meshDS->SetMeshElementOnShape(ME, edgeID);
+      for ( size_t i = 1; i < Nodes2.size(); i++ ) {
         ME = myHelper->AddEdge( Nodes2[i-1], Nodes2[i] );
         if(ME) meshDS->SetMeshElementOnShape(ME, edgeID);
       }
@@ -895,7 +931,6 @@ bool StdMeshers_RadialQuadrangle_1D2D::Compute(SMESH_Mesh&         aMesh,
   //cout<<"Nodes1.size() = "<<Nodes1.size()<<"   Pnts2d1.Length() = "<<Pnts2d1.Length()<<endl;
   for(; i<Angles.Length(); i++) {
     vector< const SMDS_MeshNode* > tmpNodes;
-    tmpNodes.reserve(Nodes1.size());
     gp_Trsf aTrsf;
     gp_Ax1 theAxis(P0,gp_Dir(Axis));
     aTrsf.SetRotation( theAxis, Angles.Value(i) );
@@ -913,12 +948,12 @@ bool StdMeshers_RadialQuadrangle_1D2D::Compute(SMESH_Mesh&         aMesh,
       aTrsf2d.Transforms( cx, cy );
       // set node on face
       meshDS->SetNodeOnFace( node, faceID, cx, cy );
-      tmpNodes[j-1] = node;
+      tmpNodes.push_back(node);
     }
     // create faces
-    tmpNodes[Points.Length()] = CNodes[i];
+    tmpNodes.push_back( CNodes[i] );
     // quad
-    for(j=0; j<Nodes1.size()-1; j++) {
+    for ( j = 0; j < (int)Nodes1.size() - 1; j++ ) {
       SMDS_MeshFace* MF;
       if(IsForward)
         MF = myHelper->AddFace( tmpNodes[j], Nodes1[j],
@@ -934,14 +969,14 @@ bool StdMeshers_RadialQuadrangle_1D2D::Compute(SMESH_Mesh&         aMesh,
       MF = myHelper->AddFace( NC, Nodes1[0], tmpNodes[0] );
     else
       MF = myHelper->AddFace( NC, tmpNodes[0], Nodes1[0] );
-    if(MF) meshDS->SetMeshElementOnShape(MF, faceID);
-    for(j=0; j<Nodes1.size(); j++) {
+    if ( MF ) meshDS->SetMeshElementOnShape(MF, faceID);
+    for ( j = 0; j < (int) Nodes1.size(); j++ ) {
       Nodes1[j] = tmpNodes[j];
     }
   }
   // create last faces
   // quad
-  for(i=0; i<Nodes1.size()-1; i++) {
+  for ( i = 0; i < (int)Nodes1.size()-1; i++ ) {
     SMDS_MeshFace* MF;
     if(IsForward)
       MF = myHelper->AddFace( Nodes2[i], Nodes1[i],
@@ -1071,7 +1106,7 @@ bool StdMeshers_RadialQuadrangle_1D2D::computeLayerPositions(const gp_Pnt&
       }
       else {
         
-        if ( myLayerPositions.size() != nbNodes )
+        if ((int) myLayerPositions.size() != nbNodes )
           return error("Radial edge is meshed by other algorithm");
       }
     }
@@ -1161,7 +1196,7 @@ bool StdMeshers_RadialQuadrangle_1D2D::Evaluate(SMESH_Mesh& aMesh,
     ok = !aResMap.count( aMesh.GetSubMesh(LinEdge1) );
     if ( !ok ) {
       const vector<int>& aVec = aResMap[ aMesh.GetSubMesh(LinEdge1) ];
-      ok = ( aVec[SMDSEntity_Node] == myLayerPositions.size() );
+      ok = ( aVec[SMDSEntity_Node] == (int) myLayerPositions.size() );
     }
     if(ok) {
       ok = algo1d->EvaluateCircularEdge( aMesh, CircEdge, aResMap );
@@ -1208,7 +1243,7 @@ bool StdMeshers_RadialQuadrangle_1D2D::Evaluate(SMESH_Mesh& aMesh,
       // other curve not line
       return error(COMPERR_BAD_SHAPE);
     }
-    int nbLayers = myLayerPositions.size();
+    size_t nbLayers = myLayerPositions.size();
     computeLayerPositions( P0, P1, LinEdge2 );
     if ( nbLayers != myLayerPositions.size() )
       return error("Different hypotheses apply to radial edges");
@@ -1219,7 +1254,7 @@ bool StdMeshers_RadialQuadrangle_1D2D::Evaluate(SMESH_Mesh& aMesh,
         ok = true; // override other 1d hyps
       else {
         const vector<int>& aVec = aResMap[ aMesh.GetSubMesh(LinEdge1) ];
-        ok = ( aVec[SMDSEntity_Node] == myLayerPositions.size() );
+        ok = ( aVec[SMDSEntity_Node] == (int) myLayerPositions.size() );
       }
     }
     if( ok && aResMap.count( aMesh.GetSubMesh(LinEdge2) )) {
@@ -1227,7 +1262,7 @@ bool StdMeshers_RadialQuadrangle_1D2D::Evaluate(SMESH_Mesh& aMesh,
         ok = true; // override other 1d hyps
       else {
         const vector<int>& aVec = aResMap[ aMesh.GetSubMesh(LinEdge2) ];
-        ok = ( aVec[SMDSEntity_Node] == myLayerPositions.size() );
+        ok = ( aVec[SMDSEntity_Node] == (int) myLayerPositions.size() );
       }
     }
     if(ok) {
@@ -1290,18 +1325,20 @@ bool StdMeshers_RadialQuadrangle_1D2D::Evaluate(SMESH_Mesh& aMesh,
 
 //================================================================================
 /*!
- * \brief Return true if applied compute mesh on this shape
+ * \brief Return true if the algorithm can compute mesh on this shape
  */
 //================================================================================
 
 bool StdMeshers_RadialQuadrangle_1D2D::IsApplicable( const TopoDS_Shape & aShape, bool toCheckAll )
 {
   int nbFoundFaces = 0;
-  for (TopExp_Explorer exp( aShape, TopAbs_FACE ); exp.More(); exp.Next(), ++nbFoundFaces ){
+  for (TopExp_Explorer exp( aShape, TopAbs_FACE ); exp.More(); exp.Next(), ++nbFoundFaces )
+  {
     TopoDS_Edge CircEdge, LinEdge1, LinEdge2;
     int nbe = analyseFace( exp.Current(), CircEdge, LinEdge1, LinEdge2 );
     Handle(Geom_Circle) aCirc = Handle(Geom_Circle)::DownCast( getCurve( CircEdge ));
-    bool ok = ( nbe <= 3 && nbe >= 1 && !aCirc.IsNull() );
+    bool ok = ( nbe <= 3 && nbe >= 1 && !aCirc.IsNull() &&
+                isCornerInsideCircle( CircEdge, LinEdge1, LinEdge2 ));
     if( toCheckAll  && !ok ) return false;
     if( !toCheckAll && ok  ) return true;
   }