Merge from BR_Dev_For_6_3_1 03/06/2011

diff --git a/idl/SMESH_Gen.idl b/idl/SMESH_Gen.idl
index a8c315c1927a05f4a70feb8760c6e65ae6da4dd7..9d70e066e478107c4560f28cc3bc379e1bd9e589 100644 (file)
--- a/idl/SMESH_Gen.idl
+++ b/idl/SMESH_Gen.idl
@@ -70,7 +70,8 @@ module SMESH
   const long Tag_EdgeGroups             = 12;
   const long Tag_FaceGroups             = 13;
   const long Tag_VolumeGroups           = 14;
   const long Tag_EdgeGroups             = 12;
   const long Tag_FaceGroups             = 13;
   const long Tag_VolumeGroups           = 14;

-  const long Tag_LastGroup              = 14;
+  const long Tag_0DElementsGroups       = 15;
+  const long Tag_LastGroup              = 15;

 
   /*!
    * Hypothesis definintion error
 
   /*!
    * Hypothesis definintion error

diff --git a/src/SMDS/SMDS_VolumeTool.hxx b/src/SMDS/SMDS_VolumeTool.hxx
index b31bed788cb517035270eff83f97f6db65571931..9cc0ae829eec00b44abbc09abd1966b72f3e5691 100644 (file)
--- a/src/SMDS/SMDS_VolumeTool.hxx
+++ b/src/SMDS/SMDS_VolumeTool.hxx
@@ -139,12 +139,16 @@ class SMDS_EXPORT SMDS_VolumeTool
   // To comfort link iteration, the array
   // length == NbFaceNodes( faceIndex ) + 1 and
   // the last node index == the first one.
   // To comfort link iteration, the array
   // length == NbFaceNodes( faceIndex ) + 1 and
   // the last node index == the first one.

+  // NOTE: for the quadratic volume, node indoces are in the order the nodes encounter
+  // in face boundary and not the order they are in the mesh face

 
   const SMDS_MeshNode** GetFaceNodes( int faceIndex );
   // Return the array of face nodes.
   // To comfort link iteration, the array
   // length == NbFaceNodes( faceIndex ) + 1 and
   // the last node == the first one.
 
   const SMDS_MeshNode** GetFaceNodes( int faceIndex );
   // Return the array of face nodes.
   // To comfort link iteration, the array
   // length == NbFaceNodes( faceIndex ) + 1 and
   // the last node == the first one.

+  // NOTE: for the quadratic volume, nodes are in the order they encounter in face boundary
+  // and not the order they are in the mesh face

   // WARNING: do not modify the array, some methods
   //          work basing on its contents
 
   // WARNING: do not modify the array, some methods
   //          work basing on its contents
 

diff --git a/src/SMESH/SMESH_Gen.cxx b/src/SMESH/SMESH_Gen.cxx
index 56da2ac331be6e4424520f4588ba891ccbe4ed5b..8b2840eb1e5cb6d3d7983528f02fadb9b2c0e33b 100644 (file)
--- a/src/SMESH/SMESH_Gen.cxx
+++ b/src/SMESH/SMESH_Gen.cxx
@@ -25,14 +25,18 @@
 //  Author : Paul RASCLE, EDF
 //  Module : SMESH
 //
 //  Author : Paul RASCLE, EDF
 //  Module : SMESH
 //

-#define CHRONODEF
+
+//#define CHRONODEF
+

 #include "SMESH_Gen.hxx"
 #include "SMESH_Gen.hxx"

-#include "SMESH_subMesh.hxx"
-#include "SMESH_HypoFilter.hxx"
-#include "SMESHDS_Document.hxx"
+
+#include "SMDS_Mesh.hxx"

 #include "SMDS_MeshElement.hxx"
 #include "SMDS_MeshNode.hxx"
 #include "SMDS_MeshElement.hxx"
 #include "SMDS_MeshNode.hxx"

-#include "SMDS_Mesh.hxx"
+#include "SMESHDS_Document.hxx"
+#include "SMESH_HypoFilter.hxx"
+#include "SMESH_MesherHelper.hxx"
+#include "SMESH_subMesh.hxx"

 
 #include "utilities.h"
 #include "OpUtil.hxx"
 
 #include "utilities.h"
 #include "OpUtil.hxx"


@@ -348,6 +352,15 @@ bool SMESH_Gen::Compute(SMESH_Mesh &          aMesh,
   MESSAGE("Number of cell objects " << SMDS_MeshCell::nbCells);
   //myMesh->dumpGrid();
   //aMesh.GetMeshDS()->Modified();
   MESSAGE("Number of cell objects " << SMDS_MeshCell::nbCells);
   //myMesh->dumpGrid();
   //aMesh.GetMeshDS()->Modified();

+
+  // fix quadratic mesh by bending iternal links near concave boundary
+  if ( aShape.IsSame( aMesh.GetShapeToMesh() ) &&
+       !aShapesId ) // not preview
+  {
+    SMESH_MesherHelper aHelper( aMesh );
+    if ( aHelper.IsQuadraticMesh() != SMESH_MesherHelper::LINEAR )
+      aHelper.FixQuadraticElements();
+  }

   return ret;
 }
 
   return ret;
 }
 

diff --git a/src/SMESH/SMESH_MeshEditor.cxx b/src/SMESH/SMESH_MeshEditor.cxx
index 873fb5799ec85f015c95350ef8603d2941e7b47c..db3d48077d56a5d71a5d1b9d286e4587a2cf623e 100644 (file)
--- a/src/SMESH/SMESH_MeshEditor.cxx
+++ b/src/SMESH/SMESH_MeshEditor.cxx
@@ -4019,93 +4019,111 @@ void SMESH_MeshEditor::makeWalls (TNodeOfNodeListMap &     mapNewNodes,
         for ( int iStep = 0; iStep < nbSteps; iStep++ )  {
           vTool.Set( *v );
           vTool.SetExternalNormal();
         for ( int iStep = 0; iStep < nbSteps; iStep++ )  {
           vTool.Set( *v );
           vTool.SetExternalNormal();

+          const int nextShift = vTool.IsForward() ? +1 : -1;

           list< int >::iterator ind = freeInd.begin();
           list< const SMDS_MeshElement* >::iterator srcEdge = srcEdges.begin();
           for ( ; ind != freeInd.end(); ++ind, ++srcEdge ) // loop on free faces
           {
             const SMDS_MeshNode** nodes = vTool.GetFaceNodes( *ind );
             int nbn = vTool.NbFaceNodes( *ind );
           list< int >::iterator ind = freeInd.begin();
           list< const SMDS_MeshElement* >::iterator srcEdge = srcEdges.begin();
           for ( ; ind != freeInd.end(); ++ind, ++srcEdge ) // loop on free faces
           {
             const SMDS_MeshNode** nodes = vTool.GetFaceNodes( *ind );
             int nbn = vTool.NbFaceNodes( *ind );

-            switch ( nbn ) {
-            case 3: { ///// triangle
-              const SMDS_MeshFace * f = aMesh->FindFace( nodes[ 0 ], nodes[ 1 ], nodes[ 2 ]);
-              if ( !f )
-                myLastCreatedElems.Append(aMesh->AddFace( nodes[ 0 ], nodes[ 1 ], nodes[ 2 ] ));
-              else if ( nodes[ 1 ] != f->GetNodeWrap( f->GetNodeIndex( nodes[ 0 ] ) + 1 ))
+            if ( ! (*v)->IsPoly() )
+              switch ( nbn ) {
+              case 3: { ///// triangle
+                const SMDS_MeshFace * f = aMesh->FindFace( nodes[ 0 ], nodes[ 1 ], nodes[ 2 ]);
+                if ( !f ||
+                     nodes[ 1 ] != f->GetNodeWrap( f->GetNodeIndex( nodes[ 0 ]) + nextShift ))

                 {
                 {

-                  myLastCreatedElems.Append(aMesh->AddFace( nodes[ 0 ], nodes[ 1 ], nodes[ 2 ] ));
-                  aMesh->RemoveElement(f);
+                  const SMDS_MeshNode* newOrder[3] = { nodes[ 1 - nextShift ],
+                                                       nodes[ 1 ],
+                                                       nodes[ 1 + nextShift ] };
+                  if ( f )
+                    aMesh->ChangeElementNodes( f, &newOrder[0], nbn );
+                  else
+                    myLastCreatedElems.Append(aMesh->AddFace( newOrder[ 0 ], newOrder[ 1 ],
+                                                              newOrder[ 2 ] ));

                 }
                 }

-              break;
-            }
-            case 4: { ///// quadrangle
-              const SMDS_MeshFace * f = aMesh->FindFace( nodes[ 0 ], nodes[ 1 ], nodes[ 2 ], nodes[ 3 ]);
-              if ( !f )
-                myLastCreatedElems.Append(aMesh->AddFace( nodes[ 0 ], nodes[ 1 ], nodes[ 2 ], nodes[ 3 ] ));
-              else if ( nodes[ 1 ] != f->GetNodeWrap( f->GetNodeIndex( nodes[ 0 ] ) + 1 ))
+                break;
+              }
+              case 4: { ///// quadrangle
+                const SMDS_MeshFace * f =
+                  aMesh->FindFace( nodes[ 0 ], nodes[ 1 ], nodes[ 2 ], nodes[ 3 ]);
+                if ( !f ||
+                     nodes[ 1 ] != f->GetNodeWrap( f->GetNodeIndex( nodes[ 0 ]) + nextShift ))

                 {
                 {

-                  myLastCreatedElems.Append(aMesh->AddFace( nodes[ 0 ], nodes[ 1 ], nodes[ 2 ], nodes[ 3 ] ));
-                  aMesh->RemoveElement(f);
-                }
-              break;
-            }
-            default:
-              if( (*v)->IsQuadratic() ) {
-                if(nbn==6) { /////// quadratic triangle
-                  const SMDS_MeshFace * f = aMesh->FindFace( nodes[0], nodes[2], nodes[4],
-                                                             nodes[1], nodes[3], nodes[5] );
-                  if ( !f ) {
-                    myLastCreatedElems.Append(aMesh->AddFace(nodes[0], nodes[2], nodes[4],
-                                                             nodes[1], nodes[3], nodes[5]));
-                  }
-                  else if ( nodes[ 2 ] != f->GetNodeWrap( f->GetNodeIndex( nodes[ 0 ] ) + 1 )) {
-                    const SMDS_MeshNode** tmpnodes = new const SMDS_MeshNode*[6];
-                    tmpnodes[0] = nodes[0];
-                    tmpnodes[1] = nodes[2];
-                    tmpnodes[2] = nodes[4];
-                    tmpnodes[3] = nodes[1];
-                    tmpnodes[4] = nodes[3];
-                    tmpnodes[5] = nodes[5];
-                    myLastCreatedElems.Append(aMesh->AddFace(nodes[0], nodes[2], nodes[4],
-                                                             nodes[1], nodes[3], nodes[5]));
-                    aMesh->RemoveElement(f);
-                  }
+                  const SMDS_MeshNode* newOrder[4] = { nodes[ 0 ], nodes[ 2-nextShift ],
+                                                       nodes[ 2 ], nodes[ 2+nextShift ] };
+                  if ( f )
+                    aMesh->ChangeElementNodes( f, &newOrder[0], nbn );
+                  else
+                    myLastCreatedElems.Append(aMesh->AddFace( newOrder[ 0 ], newOrder[ 1 ],
+                                                              newOrder[ 2 ], newOrder[ 3 ]));

                 }
                 }

-                else {       /////// quadratic quadrangle
-                  const SMDS_MeshFace * f = aMesh->FindFace( nodes[0], nodes[2], nodes[4], nodes[6],
-                                                             nodes[1], nodes[3], nodes[5], nodes[7] );
-                  if ( !f ) {
-                    myLastCreatedElems.Append(aMesh->AddFace(nodes[0], nodes[2], nodes[4], nodes[6],
-                                                             nodes[1], nodes[3], nodes[5], nodes[7]));
-                  }
-                  else if ( nodes[ 2 ] != f->GetNodeWrap( f->GetNodeIndex( nodes[ 0 ] ) + 1 )) {
-                    const SMDS_MeshNode** tmpnodes = new const SMDS_MeshNode*[8];
-                    tmpnodes[0] = nodes[0];
-                    tmpnodes[1] = nodes[2];
-                    tmpnodes[2] = nodes[4];
-                    tmpnodes[3] = nodes[6];
-                    tmpnodes[4] = nodes[1];
-                    tmpnodes[5] = nodes[3];
-                    tmpnodes[6] = nodes[5];
-                    tmpnodes[7] = nodes[7];
-                    myLastCreatedElems.Append(aMesh->AddFace(nodes[0], nodes[2], nodes[4], nodes[6],
-                                                             nodes[1], nodes[3], nodes[5], nodes[7]));
-                    aMesh->RemoveElement(f);
-                  }
+                break;
+              }
+              case 6: { /////// quadratic triangle
+                const SMDS_MeshFace * f = aMesh->FindFace( nodes[0], nodes[2], nodes[4],
+                                                           nodes[1], nodes[3], nodes[5] );
+                if ( !f ||
+                     nodes[2] != f->GetNodeWrap( f->GetNodeIndex( nodes[0] ) + 2*nextShift ))
+                {
+                  const SMDS_MeshNode* newOrder[6] = { nodes[2 - 2*nextShift],
+                                                       nodes[2],
+                                                       nodes[2 + 2*nextShift],
+                                                       nodes[3 - 2*nextShift],
+                                                       nodes[3],
+                                                       nodes[3 + 2*nextShift]};
+                  if ( f )
+                    aMesh->ChangeElementNodes( f, &newOrder[0], nbn );
+                  else
+                    myLastCreatedElems.Append(aMesh->AddFace( newOrder[ 0 ],
+                                                              newOrder[ 1 ],
+                                                              newOrder[ 2 ],
+                                                              newOrder[ 3 ],
+                                                              newOrder[ 4 ],
+                                                              newOrder[ 5 ] ));

                 }
                 }

+                break;

               }
               }

-              else { //////// polygon
-                vector<const SMDS_MeshNode*> polygon_nodes ( nodes, &nodes[nbn] );
-                const SMDS_MeshFace * f = aMesh->FindFace( polygon_nodes );
-                if ( !f )
+              default:       /////// quadratic quadrangle
+                const SMDS_MeshFace * f = aMesh->FindFace( nodes[0], nodes[2], nodes[4], nodes[6],
+                                                           nodes[1], nodes[3], nodes[5], nodes[7] );
+                if ( !f ||
+                     nodes[ 2 ] != f->GetNodeWrap( f->GetNodeIndex( nodes[ 0 ] ) + 2*nextShift ))
+                {
+                  const SMDS_MeshNode* newOrder[8] = { nodes[0],
+                                                       nodes[4 - 2*nextShift],
+                                                       nodes[4],
+                                                       nodes[4 + 2*nextShift],
+                                                       nodes[1],
+                                                       nodes[5 - 2*nextShift],
+                                                       nodes[5],
+                                                       nodes[5 + 2*nextShift] };
+                  if ( f )
+                    aMesh->ChangeElementNodes( f, &newOrder[0], nbn );
+                  else
+                    myLastCreatedElems.Append(aMesh->AddFace(newOrder[ 0 ], newOrder[ 1 ],
+                                                             newOrder[ 2 ], newOrder[ 3 ],
+                                                             newOrder[ 4 ], newOrder[ 5 ],
+                                                             newOrder[ 6 ], newOrder[ 7 ]));
+                }
+              } // switch ( nbn )
+
+            else { //////// polygon
+
+              vector<const SMDS_MeshNode*> polygon_nodes ( nodes, &nodes[nbn] );
+              const SMDS_MeshFace * f = aMesh->FindFace( polygon_nodes );
+              if ( !f ||
+                   nodes[ 1 ] != f->GetNodeWrap( f->GetNodeIndex( nodes[ 0 ] ) + nextShift ))
+              {
+                if ( !vTool.IsForward() )
+                  std::reverse( polygon_nodes.begin(), polygon_nodes.end());
+                if ( f )
+                  aMesh->ChangeElementNodes( f, &polygon_nodes[0], nbn );
+                else

                   myLastCreatedElems.Append(aMesh->AddPolygonalFace(polygon_nodes));
                   myLastCreatedElems.Append(aMesh->AddPolygonalFace(polygon_nodes));

-                else if ( nodes[ 1 ] != f->GetNodeWrap( f->GetNodeIndex( nodes[ 0 ] ) + 1 ))
-                  {
-                  // TODO problem ChangeElementNodes : not the same number of nodes, not the same type
-                  MESSAGE("ChangeElementNodes");
-                  aMesh->ChangeElementNodes( f, nodes, nbn );
-                  }

               }
             }
               }
             }

+

             while ( srcElements.Length() < myLastCreatedElems.Length() )
               srcElements.Append( *srcEdge );
 
             while ( srcElements.Length() < myLastCreatedElems.Length() )
               srcElements.Append( *srcEdge );
 


@@ -4114,8 +4132,9 @@ void SMESH_MeshEditor::makeWalls (TNodeOfNodeListMap &     mapNewNodes,
           // go to the next volume
           iVol = 0;
           while ( iVol++ < nbVolumesByStep ) v++;
           // go to the next volume
           iVol = 0;
           while ( iVol++ < nbVolumesByStep ) v++;

-        }
-      }
+
+        } // loop on steps
+      } // loop on volumes of one step

     } // sweep free links into faces
 
     // Make a ceiling face with a normal external to a volume
     } // sweep free links into faces
 
     // Make a ceiling face with a normal external to a volume


@@ -7368,8 +7387,11 @@ void SMESH_MeshEditor::MergeNodes (TListOfListOfNodes & theGroupsOfNodes)
       }
       curNodes[ iCur ] = n;
       bool isUnique = nodeSet.insert( n ).second;
       }
       curNodes[ iCur ] = n;
       bool isUnique = nodeSet.insert( n ).second;

-      if ( isUnique )
+      if ( isUnique ) {

         uniqueNodes[ iUnique++ ] = n;
         uniqueNodes[ iUnique++ ] = n;

+        if ( nbRepl && iRepl[ nbRepl-1 ] == iCur )
+          --nbRepl; // n do not stick to a node of the elem
+      }

       iCur++;
     }
 
       iCur++;
     }
 


@@ -7481,7 +7503,7 @@ void SMESH_MeshEditor::MergeNodes (TListOfListOfNodes & theGroupsOfNodes)
         }
 
         continue;
         }
 
         continue;

-      }
+      } // poly element

 
       // Regular elements
       // TODO not all the possible cases are solved. Find something more generic?
 
       // Regular elements
       // TODO not all the possible cases are solved. Find something more generic?


@@ -7654,8 +7676,8 @@ void SMESH_MeshEditor::MergeNodes (TListOfListOfNodes & theGroupsOfNodes)
         isOk = false;
         SMDS_VolumeTool hexa (elem);
         hexa.SetExternalNormal();
         isOk = false;
         SMDS_VolumeTool hexa (elem);
         hexa.SetExternalNormal();

-        if ( nbUniqueNodes == 4 && nbRepl == 6 ) {
-          //////////////////////// ---> tetrahedron
+        if ( nbUniqueNodes == 4 && nbRepl == 4 ) {
+          //////////////////////// HEX ---> 1 tetrahedron

           for ( int iFace = 0; iFace < 6; iFace++ ) {
             const int *ind = hexa.GetFaceNodesIndices( iFace ); // indices of face nodes
             if (curNodes[ind[ 0 ]] == curNodes[ind[ 1 ]] &&
           for ( int iFace = 0; iFace < 6; iFace++ ) {
             const int *ind = hexa.GetFaceNodesIndices( iFace ); // indices of face nodes
             if (curNodes[ind[ 0 ]] == curNodes[ind[ 1 ]] &&


@@ -7665,12 +7687,9 @@ void SMESH_MeshEditor::MergeNodes (TListOfListOfNodes & theGroupsOfNodes)
               int iOppFace = hexa.GetOppFaceIndex( iFace );
               ind = hexa.GetFaceNodesIndices( iOppFace );
               int nbStick = 0;
               int iOppFace = hexa.GetOppFaceIndex( iFace );
               ind = hexa.GetFaceNodesIndices( iOppFace );
               int nbStick = 0;

-              iUnique = 2; // reverse a tetrahedron bottom

               for ( iCur = 0; iCur < 4 && nbStick < 2; iCur++ ) {
                 if ( curNodes[ind[ iCur ]] == curNodes[ind[ iCur + 1 ]] )
                   nbStick++;
               for ( iCur = 0; iCur < 4 && nbStick < 2; iCur++ ) {
                 if ( curNodes[ind[ iCur ]] == curNodes[ind[ iCur + 1 ]] )
                   nbStick++;

-                else if ( iUnique >= 0 )
-                  uniqueNodes[ iUnique-- ] = curNodes[ind[ iCur ]];

               }
               if ( nbStick == 1 ) {
                 // ... and the opposite one - into a triangle.
               }
               if ( nbStick == 1 ) {
                 // ... and the opposite one - into a triangle.


@@ -7683,6 +7702,45 @@ void SMESH_MeshEditor::MergeNodes (TListOfListOfNodes & theGroupsOfNodes)
             }
           }
         }
             }
           }
         }

+        else if ( nbUniqueNodes == 6 && nbRepl == 2 ) {
+          //////////////////////// HEX ---> 1 prism
+          int nbTria = 0, iTria[3];
+          const int *ind; // indices of face nodes
+          // look for triangular faces
+          for ( int iFace = 0; iFace < 6 && nbTria < 3; iFace++ ) {
+            ind = hexa.GetFaceNodesIndices( iFace );
+            TIDSortedNodeSet faceNodes;
+            for ( iCur = 0; iCur < 4; iCur++ )
+              faceNodes.insert( curNodes[ind[iCur]] );
+            if ( faceNodes.size() == 3 )
+              iTria[ nbTria++ ] = iFace;
+          }
+          // check if triangles are opposite
+          if ( nbTria == 2 && iTria[0] == hexa.GetOppFaceIndex( iTria[1] ))
+          {
+            isOk = true;
+            // set nodes of the bottom triangle
+            ind = hexa.GetFaceNodesIndices( iTria[ 0 ]);
+            vector<int> indB;
+            for ( iCur = 0; iCur < 4; iCur++ )
+              if ( ind[iCur] != iRepl[0] && ind[iCur] != iRepl[1])
+                indB.push_back( ind[iCur] );
+            if ( !hexa.IsForward() )
+              std::swap( indB[0], indB[2] );
+            for ( iCur = 0; iCur < 3; iCur++ )
+              uniqueNodes[ iCur ] = curNodes[indB[iCur]];
+            // set nodes of the top triangle
+            const int *indT = hexa.GetFaceNodesIndices( iTria[ 1 ]);
+            for ( iCur = 0; iCur < 3; ++iCur )
+              for ( int j = 0; j < 4; ++j )
+                if ( hexa.IsLinked( indB[ iCur ], indT[ j ] ))
+                {
+                  uniqueNodes[ iCur + 3 ] = curNodes[ indT[ j ]];
+                  break;
+                }
+          }
+          break;
+        }

         else if (nbUniqueNodes == 5 && nbRepl == 4 ) {
           //////////////////// HEXAHEDRON ---> 2 tetrahedrons
           for ( int iFace = 0; iFace < 6; iFace++ ) {
         else if (nbUniqueNodes == 5 && nbRepl == 4 ) {
           //////////////////// HEXAHEDRON ---> 2 tetrahedrons
           for ( int iFace = 0; iFace < 6; iFace++ ) {


@@ -7820,6 +7878,10 @@ void SMESH_MeshEditor::MergeNodes (TListOfListOfNodes & theGroupsOfNodes)
             }
           }
         } // if ( nbUniqueNodes == 6 && nbRepl == 4 )
             }
           }
         } // if ( nbUniqueNodes == 6 && nbRepl == 4 )

+        else
+        {
+          MESSAGE("MergeNodes() removes hexahedron "<< elem);
+        }

         break;
       } // HEXAHEDRON
 
         break;
       } // HEXAHEDRON
 


@@ -9366,7 +9428,7 @@ void SMESH_MeshEditor::ConvertToQuadratic(const bool theForce3d)
     }
   }
 
     }
   }
 

-  if ( !theForce3d  && !getenv("NO_FixQuadraticElements"))
+  if ( !theForce3d )

   { // setenv NO_FixQuadraticElements to know if FixQuadraticElements() is guilty of bad conversion
     aHelper.SetSubShape(0); // apply FixQuadraticElements() to the whole mesh
     aHelper.FixQuadraticElements();
   { // setenv NO_FixQuadraticElements to know if FixQuadraticElements() is guilty of bad conversion
     aHelper.SetSubShape(0); // apply FixQuadraticElements() to the whole mesh
     aHelper.FixQuadraticElements();

diff --git a/src/SMESH/SMESH_MesherHelper.cxx b/src/SMESH/SMESH_MesherHelper.cxx
index 5d26bb47295a3363ba87d783c15ff36ac0545add..698446b5c53268c5c592ee118fce14bbd3a95dee 100644 (file)
--- a/src/SMESH/SMESH_MesherHelper.cxx
+++ b/src/SMESH/SMESH_MesherHelper.cxx
@@ -40,7 +40,7 @@
 #include <GeomAPI_ProjectPointOnCurve.hxx>
 #include <GeomAPI_ProjectPointOnSurf.hxx>
 #include <Geom_Curve.hxx>
 #include <GeomAPI_ProjectPointOnCurve.hxx>
 #include <GeomAPI_ProjectPointOnSurf.hxx>
 #include <Geom_Curve.hxx>

-//#include <Geom_RectangularTrimmedSurface.hxx>
+#include <Geom_RectangularTrimmedSurface.hxx>

 #include <Geom_Surface.hxx>
 #include <ShapeAnalysis.hxx>
 #include <TopExp.hxx>
 #include <Geom_Surface.hxx>
 #include <ShapeAnalysis.hxx>
 #include <TopExp.hxx>


@@ -735,7 +735,14 @@ gp_XY SMESH_MesherHelper::GetMiddleUV(const Handle(Geom_Surface)& surface,
                                       const gp_XY&                p1,
                                       const gp_XY&                p2)
 {
                                       const gp_XY&                p1,
                                       const gp_XY&                p2)
 {

-  return applyIn2D( surface, p1, p2, & AverageUV );
+  // NOTE:
+  // the proper place of getting basic surface seems to be in applyIn2D()
+  // but we put it here to decrease a risk of regressions just before releasing a version
+  Handle(Geom_Surface) surf = surface;
+  while ( !surf.IsNull() && surf->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface )))
+    surf = Handle(Geom_RectangularTrimmedSurface)::DownCast( surf )->BasisSurface();
+
+  return applyIn2D( surf, p1, p2, & AverageUV );

 }
 
 //=======================================================================
 }
 
 //=======================================================================


@@ -1830,6 +1837,57 @@ double SMESH_MesherHelper::GetOtherParam(const double param) const
   return fabs(param-myPar1[i]) < fabs(param-myPar2[i]) ? myPar2[i] : myPar1[i];
 }
 
   return fabs(param-myPar1[i]) < fabs(param-myPar2[i]) ? myPar2[i] : myPar1[i];
 }
 

+namespace {
+
+  //=======================================================================
+  /*!
+   * \brief Iterator on ancestors of the given type
+   */
+  //=======================================================================
+
+  struct TAncestorsIterator : public SMDS_Iterator<const TopoDS_Shape*>
+  {
+    TopTools_ListIteratorOfListOfShape _ancIter;
+    TopAbs_ShapeEnum                   _type;
+    TopTools_MapOfShape                _encountered;
+    TAncestorsIterator( const TopTools_ListOfShape& ancestors, TopAbs_ShapeEnum type)
+      : _ancIter( ancestors ), _type( type )
+    {
+      if ( _ancIter.More() ) {
+        if ( _ancIter.Value().ShapeType() != _type ) next();
+        else _encountered.Add( _ancIter.Value() );
+      }
+    }
+    virtual bool more()
+    {
+      return _ancIter.More();
+    }
+    virtual const TopoDS_Shape* next()
+    {
+      const TopoDS_Shape* s = _ancIter.More() ? & _ancIter.Value() : 0;
+      if ( _ancIter.More() )
+        for ( _ancIter.Next();  _ancIter.More(); _ancIter.Next())
+          if ( _ancIter.Value().ShapeType() == _type && _encountered.Add( _ancIter.Value() ))
+            break;
+      return s;
+    }
+  };
+
+} // namespace
+
+//=======================================================================
+/*!
+ * \brief Return iterator on ancestors of the given type
+ */
+//=======================================================================
+
+PShapeIteratorPtr SMESH_MesherHelper::GetAncestors(const TopoDS_Shape& shape,
+                                                   const SMESH_Mesh&   mesh,
+                                                   TopAbs_ShapeEnum    ancestorType)
+{
+  return PShapeIteratorPtr( new TAncestorsIterator( mesh.GetAncestors(shape), ancestorType));
+}
+

 //#include <Perf_Meter.hxx>
 
 //=======================================================================
 //#include <Perf_Meter.hxx>
 
 //=======================================================================


@@ -1882,7 +1940,7 @@ namespace { // Structures used by FixQuadraticElements()
     void Move(const gp_Vec& move, bool sum=false) const
     { _nodeMove += move; _nbMoves += sum ? (_nbMoves==0) : 1; }
     gp_XYZ Move() const { return _nodeMove.XYZ() / _nbMoves; }
     void Move(const gp_Vec& move, bool sum=false) const
     { _nodeMove += move; _nbMoves += sum ? (_nbMoves==0) : 1; }
     gp_XYZ Move() const { return _nodeMove.XYZ() / _nbMoves; }

-    bool IsMoved() const { return (_nbMoves > 0 && !IsStraight()); }
+    bool IsMoved() const { return (_nbMoves > 0 /*&& !IsStraight()*/); }

     bool IsStraight() const
     { return isStraightLink( (XYZ(node1())-XYZ(node2())).SquareModulus(),
                              _nodeMove.SquareMagnitude());
     bool IsStraight() const
     { return isStraightLink( (XYZ(node1())-XYZ(node2())).SquareModulus(),
                              _nodeMove.SquareMagnitude());


@@ -1928,6 +1986,8 @@ namespace { // Structures used by FixQuadraticElements()
     const QLink* operator->() const { return _qlink; }
 
     gp_Vec Normal() const;
     const QLink* operator->() const { return _qlink; }
 
     gp_Vec Normal() const;

+
+    bool IsStraight() const;

   };
   // --------------------------------------------------------------------
   typedef list< TChainLink > TChain;
   };
   // --------------------------------------------------------------------
   typedef list< TChainLink > TChain;


@@ -2105,9 +2165,10 @@ namespace { // Structures used by FixQuadraticElements()
 //               chLink = chain.insert( chain.begin(), TChainLink(face->_sides[i]));
             // add a face to a chained link and put a continues face in the queue
             chLink->SetFace( face );
 //               chLink = chain.insert( chain.begin(), TChainLink(face->_sides[i]));
             // add a face to a chained link and put a continues face in the queue
             chLink->SetFace( face );

-            if ( face->_sides[i]->MediumPos() >= pos )
+            if ( face->_sides[i]->MediumPos() == pos )

               if ( const QFace* contFace = face->_sides[i]->GetContinuesFace( face ))
               if ( const QFace* contFace = face->_sides[i]->GetContinuesFace( face ))

-                faces.push_back( contFace );
+                if ( contFace->_sides.size() == 3 )
+                  faces.push_back( contFace );

           }
         }
         faces.pop_front();
           }
         }
         faces.pop_front();


@@ -2130,10 +2191,11 @@ namespace { // Structures used by FixQuadraticElements()
     // propagate from quadrangle to neighbour faces
     if ( link->MediumPos() >= pos ) {
       int nbLinkFaces = link->_faces.size();
     // propagate from quadrangle to neighbour faces
     if ( link->MediumPos() >= pos ) {
       int nbLinkFaces = link->_faces.size();

-      if ( nbLinkFaces == 4 || (nbLinkFaces < 4 && link->OnBoundary())) {
+      if ( nbLinkFaces == 4 || (/*nbLinkFaces < 4 && */link->OnBoundary())) {

         // hexahedral mesh or boundary quadrangles - goto a continous face
         if ( const QFace* f = link->GetContinuesFace( this ))
         // hexahedral mesh or boundary quadrangles - goto a continous face
         if ( const QFace* f = link->GetContinuesFace( this ))

-          return f->GetLinkChain( *chLink, chain, pos, error );
+          if ( f->_sides.size() == 4 )
+            return f->GetLinkChain( *chLink, chain, pos, error );

       }
       else {
         TChainLink chLink(link); // side face of prismatic mesh - visit all faces of iSide
       }
       else {
         TChainLink chLink(link); // side face of prismatic mesh - visit all faces of iSide


@@ -2327,7 +2389,7 @@ namespace { // Structures used by FixQuadraticElements()
     gp_Vec linkDir2(0,0,0);
     try {
       OCC_CATCH_SIGNALS;
     gp_Vec linkDir2(0,0,0);
     try {
       OCC_CATCH_SIGNALS;

-      if ( f1 )
+      if ( f1 && theLink->MediumPos() <= (*link1)->MediumPos() )

         len1 = f1->MoveByBoundary
           ( *link1, theRefVec, theLinks, theFaceHelper, len1, theStep-1, &linkDir1, theSign);
       else
         len1 = f1->MoveByBoundary
           ( *link1, theRefVec, theLinks, theFaceHelper, len1, theStep-1, &linkDir1, theSign);
       else


@@ -2338,7 +2400,7 @@ namespace { // Structures used by FixQuadraticElements()
     }
     try {
       OCC_CATCH_SIGNALS;
     }
     try {
       OCC_CATCH_SIGNALS;

-      if ( f2 )
+      if ( f2 && theLink->MediumPos() <= (*link2)->MediumPos() )

         len2 = f2->MoveByBoundary
           ( *link2, theRefVec, theLinks, theFaceHelper, len2, theStep-1, &linkDir2, theSign);
       else
         len2 = f2->MoveByBoundary
           ( *link2, theRefVec, theLinks, theFaceHelper, len2, theStep-1, &linkDir2, theSign);
       else


@@ -2417,7 +2479,9 @@ namespace { // Structures used by FixQuadraticElements()
 
     if ( _faces.empty() )
       return;
 
     if ( _faces.empty() )
       return;

-    int iFaceCont = -1;
+    int iFaceCont = -1, nbBoundary = 0, iBoundary[2]={-1,-1};
+    if ( _faces[0]->IsBoundary() )
+      iBoundary[ nbBoundary++ ] = 0;

     for ( int iF = 1; iFaceCont < 0 && iF < _faces.size(); ++iF )
     {
       // look for a face bounding none of volumes bound by _faces[0]
     for ( int iF = 1; iFaceCont < 0 && iF < _faces.size(); ++iF )
     {
       // look for a face bounding none of volumes bound by _faces[0]


@@ -2428,8 +2492,21 @@ namespace { // Structures used by FixQuadraticElements()
                     _faces[iF]->_volumes[iV] == _faces[0]->_volumes[1]);
       if ( !sameVol )
         iFaceCont = iF;
                     _faces[iF]->_volumes[iV] == _faces[0]->_volumes[1]);
       if ( !sameVol )
         iFaceCont = iF;

+      if ( _faces[iF]->IsBoundary() )
+        iBoundary[ nbBoundary++ ] = iF;
+    }
+    // Set continues faces: arrange _faces to have
+    // _faces[0] continues to _faces[1]
+    // _faces[2] continues to _faces[3]
+    if ( nbBoundary == 2 ) // bnd faces are continues
+    {
+      if (( iBoundary[0] < 2 ) != ( iBoundary[1] < 2 ))
+      {
+        int iNear0 = iBoundary[0] < 2 ? 1-iBoundary[0] : 5-iBoundary[0];
+        std::swap( _faces[ iBoundary[1] ], _faces[iNear0] );
+      }

     }
     }

-    if ( iFaceCont > 0 ) // continues faces found, set one by the other
+    else if ( iFaceCont > 0 ) // continues faces found

     {
       if ( iFaceCont != 1 )
         std::swap( _faces[1], _faces[iFaceCont] );
     {
       if ( iFaceCont != 1 )
         std::swap( _faces[1], _faces[iFaceCont] );


@@ -2479,6 +2556,27 @@ namespace { // Structures used by FixQuadraticElements()
     if (_qfaces[1]) norm += _qfaces[1]->_normal;
     return norm;
   }
     if (_qfaces[1]) norm += _qfaces[1]->_normal;
     return norm;
   }

+  //================================================================================
+  /*!
+   * \brief Test link curvature taking into account size of faces
+   */
+  //================================================================================
+
+  bool TChainLink::IsStraight() const
+  {
+    bool isStraight = _qlink->IsStraight();
+    if ( isStraight && _qfaces[0] && !_qfaces[1] )
+    {
+      int i = _qfaces[0]->LinkIndex( _qlink );
+      int iOpp = ( i + 2 ) % _qfaces[0]->_sides.size();
+      gp_XYZ mid1 = _qlink->MiddlePnt();
+      gp_XYZ mid2 = _qfaces[0]->_sides[ iOpp ]->MiddlePnt();
+      double faceSize2 = (mid1-mid2).SquareModulus();
+      isStraight = _qlink->_nodeMove.SquareMagnitude() < 1/3./3. * faceSize2;
+    }
+    return isStraight;
+  }
+  

   //================================================================================
   /*!
    * \brief Move medium nodes of vertical links of pentahedrons adjacent by side faces
   //================================================================================
   /*!
    * \brief Move medium nodes of vertical links of pentahedrons adjacent by side faces


@@ -2497,7 +2595,7 @@ namespace { // Structures used by FixQuadraticElements()
         bndLinks1.insert( lnk->_qlink );
       else
         interLinks.insert( lnk->_qlink );
         bndLinks1.insert( lnk->_qlink );
       else
         interLinks.insert( lnk->_qlink );

-      isCurved = isCurved || !(*lnk)->IsStraight();
+      isCurved = isCurved || !lnk->IsStraight();

     }
     if ( !isCurved )
       return; // no need to move
     }
     if ( !isCurved )
       return; // no need to move


@@ -2546,7 +2644,7 @@ namespace { // Structures used by FixQuadraticElements()
 
     for ( linkIt = linkSet.begin(); linkIt != linksEnd; ++linkIt)
     {
 
     for ( linkIt = linkSet.begin(); linkIt != linksEnd; ++linkIt)
     {

-      if ( linkIt->IsBoundary() && !(*linkIt)->IsStraight() && linkIt->_qfaces[0])
+      if ( linkIt->IsBoundary() && !linkIt->IsStraight() && linkIt->_qfaces[0])

       {
         // move iff a boundary link is bent towards inside of a face (issue 0021084)
         const QFace* face = linkIt->_qfaces[0];
       {
         // move iff a boundary link is bent towards inside of a face (issue 0021084)
         const QFace* face = linkIt->_qfaces[0];


@@ -2757,6 +2855,10 @@ namespace { // Structures used by FixQuadraticElements()
 
 void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly)
 {
 
 void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly)
 {

+  // setenv NO_FixQuadraticElements to know if FixQuadraticElements() is guilty of bad conversion
+  if ( getenv("NO_FixQuadraticElements") )
+    return;
+

   // 0. Apply algorithm to solids or geom faces
   // ----------------------------------------------
   if ( myShape.IsNull() ) {
   // 0. Apply algorithm to solids or geom faces
   // ----------------------------------------------
   if ( myShape.IsNull() ) {


@@ -2789,7 +2891,7 @@ void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly)
     }
     // fix nodes on geom faces
 #ifdef _DEBUG_
     }
     // fix nodes on geom faces
 #ifdef _DEBUG_

-    //int nbfaces = faces.Extent();
+    int nbfaces = faces.Extent(); /*avoid "unused varianbles": */ nbfaces++, nbfaces--; 

 #endif
     for ( TopTools_MapIteratorOfMapOfShape fIt( faces ); fIt.More(); fIt.Next() ) {
       MSG("FIX FACE " << nbfaces-- << " #" << GetMeshDS()->ShapeToIndex(fIt.Key()));
 #endif
     for ( TopTools_MapIteratorOfMapOfShape fIt( faces ); fIt.More(); fIt.Next() ) {
       MSG("FIX FACE " << nbfaces-- << " #" << GetMeshDS()->ShapeToIndex(fIt.Key()));


@@ -2831,15 +2933,18 @@ void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly)
   bool isCurved = false;
   //bool hasRectFaces = false;
   //set<int> nbElemNodeSet;
   bool isCurved = false;
   //bool hasRectFaces = false;
   //set<int> nbElemNodeSet;

+  SMDS_VolumeTool volTool;
+
+  TIDSortedNodeSet apexOfPyramid;
+  const int apexIndex = 4;

 
   if ( elemType == SMDSAbs_Volume )
   {
 
   if ( elemType == SMDSAbs_Volume )
   {

-    SMDS_VolumeTool volTool;

     while ( elemIt->more() ) // loop on volumes
     {
       const SMDS_MeshElement* vol = elemIt->next();
       if ( !vol->IsQuadratic() || !volTool.Set( vol ))
     while ( elemIt->more() ) // loop on volumes
     {
       const SMDS_MeshElement* vol = elemIt->next();
       if ( !vol->IsQuadratic() || !volTool.Set( vol ))

-        return; //continue;
+        return;

       for ( int iF = 0; iF < volTool.NbFaces(); ++iF ) // loop on faces of volume
       {
         int nbN = volTool.NbFaceNodes( iF );
       for ( int iF = 0; iF < volTool.NbFaces(); ++iF ) // loop on faces of volume
       {
         int nbN = volTool.NbFaceNodes( iF );


@@ -2873,6 +2978,9 @@ void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly)
                                                faceNodes[4],faceNodes[6] );
 #endif
       }
                                                faceNodes[4],faceNodes[6] );
 #endif
       }

+      // collect pyramid apexes for further correction
+      if ( vol->NbCornerNodes() == 5 )
+        apexOfPyramid.insert( vol->GetNode( apexIndex ));

     }
     set< QLink >::iterator pLink = links.begin();
     for ( ; pLink != links.end(); ++pLink )
     }
     set< QLink >::iterator pLink = links.begin();
     for ( ; pLink != links.end(); ++pLink )


@@ -2907,9 +3015,9 @@ void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly)
     return; // no curved edges of faces
 
   // 3. Compute displacement of medium nodes
     return; // no curved edges of faces
 
   // 3. Compute displacement of medium nodes

-  // -------------------------------------
+  // ---------------------------------------

 
 

-  // two loops on faces: the first is to treat boundary links, the second is for internal ones
+  // two loops on QFaces: the first is to treat boundary links, the second is for internal ones

   TopLoc_Location loc;
   // not treat boundary of volumic submesh
   int isInside = ( elemType == SMDSAbs_Volume && volumeOnly ) ? 1 : 0;
   TopLoc_Location loc;
   // not treat boundary of volumic submesh
   int isInside = ( elemType == SMDSAbs_Volume && volumeOnly ) ? 1 : 0;


@@ -2921,7 +3029,7 @@ void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly)
     for ( pFace = faces.begin(); pFace != faces.end(); ++pFace ) {
       if ( bool(isInside) == pFace->IsBoundary() )
         continue;
     for ( pFace = faces.begin(); pFace != faces.end(); ++pFace ) {
       if ( bool(isInside) == pFace->IsBoundary() )
         continue;

-      for ( int dir = 0; dir < 2; ++dir ) // 2 directions of propagation from quadrangle
+      for ( int dir = 0; dir < 2; ++dir ) // 2 directions of propagation from the quadrangle

       {
         MSG( "CHAIN");
         // make chain of links connected via continues faces
       {
         MSG( "CHAIN");
         // make chain of links connected via continues faces


@@ -2954,12 +3062,12 @@ void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly)
         {
           TChain& chain = chains[iC];
           if ( chain.empty() ) continue;
         {
           TChain& chain = chains[iC];
           if ( chain.empty() ) continue;

-          if ( chain.front()->IsStraight() && chain.back()->IsStraight() ) {
+          if ( chain.front().IsStraight() && chain.back().IsStraight() ) {

             MSG("3D straight - ignore");
             continue;
           }
           if ( chain.front()->MediumPos() > bndPos ||
             MSG("3D straight - ignore");
             continue;
           }
           if ( chain.front()->MediumPos() > bndPos ||

-               chain.back()->MediumPos() > bndPos ) {
+               chain.back() ->MediumPos() > bndPos ) {

             MSG("Internal chain - ignore");
             continue;
           }
             MSG("Internal chain - ignore");
             continue;
           }


@@ -2989,9 +3097,11 @@ void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly)
 
           TopoDS_Face face;
           bool checkUV = true;
 
           TopoDS_Face face;
           bool checkUV = true;

-          if ( !isInside ) {
-            // compute node displacement of end links in parametric space of face
-            const SMDS_MeshNode* nodeOnFace = (*(++chain.begin()))->_mediumNode;
+          if ( !isInside )
+          {
+            // compute node displacement of end links of chain in parametric space of face
+            TChainLink& linkOnFace = *(++chain.begin());
+            const SMDS_MeshNode* nodeOnFace = linkOnFace->_mediumNode;

             TopoDS_Shape f = GetSubShapeByNode( nodeOnFace, GetMeshDS() );
             if ( !f.IsNull() && f.ShapeType() == TopAbs_FACE )
             {
             TopoDS_Shape f = GetSubShapeByNode( nodeOnFace, GetMeshDS() );
             if ( !f.IsNull() && f.ShapeType() == TopAbs_FACE )
             {


@@ -3010,14 +3120,24 @@ void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly)
                 ( is1 ? move1 : move0 ).SetCoord( uvMove.X(), uvMove.Y(), 0 );
                 if ( !is1 ) // correct nodeOnFace for move1 (issue 0020919)
                   nodeOnFace = (*(++chain.rbegin()))->_mediumNode;
                 ( is1 ? move1 : move0 ).SetCoord( uvMove.X(), uvMove.Y(), 0 );
                 if ( !is1 ) // correct nodeOnFace for move1 (issue 0020919)
                   nodeOnFace = (*(++chain.rbegin()))->_mediumNode;

-                isStraight[is1] = isStraightLink( (uv2-uv1).SquareModulus(),uvMove.SquareModulus());
+                isStraight[is1] = isStraightLink( (uv2-uv1).SquareModulus(),
+                                                  10 * uvMove.SquareModulus());

               }
               }

-//               if ( move0.SquareMagnitude() < straightTol2 &&
-//                    move1.SquareMagnitude() < straightTol2 ) {

               if ( isStraight[0] && isStraight[1] ) {
                 MSG("2D straight - ignore");
                 continue; // straight - no need to move nodes of internal links
               }
               if ( isStraight[0] && isStraight[1] ) {
                 MSG("2D straight - ignore");
                 continue; // straight - no need to move nodes of internal links
               }

+
+              // check if a chain is already fixed
+              gp_XY uvm = GetNodeUV( face, linkOnFace->_mediumNode, 0, &checkUV);
+              gp_XY uv1 = GetNodeUV( face, linkOnFace->node1(), nodeOnFace, &checkUV);
+              gp_XY uv2 = GetNodeUV( face, linkOnFace->node2(), nodeOnFace, &checkUV);
+              gp_XY uv12 = GetMiddleUV( surf, uv1, uv2);
+              if (( uvm - uv12 ).SquareModulus() > 1e-10 )
+              {
+                MSG("Already fixed - ignore");
+                continue;
+              }

             }
           }
           gp_Trsf trsf;
             }
           }
           gp_Trsf trsf;


@@ -3087,60 +3207,106 @@ void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly)
   }
 
   // 4. Move nodes
   }
 
   // 4. Move nodes

-  // -----------
+  // -------------

 
 

+//   vector<const SMDS_MeshElement*> vols( 100 );
+//   vector<double>                  volSize( 100 );
+//   int nbVols;
+//   bool ok;

   for ( pLink = links.begin(); pLink != links.end(); ++pLink ) {
     if ( pLink->IsMoved() ) {
   for ( pLink = links.begin(); pLink != links.end(); ++pLink ) {
     if ( pLink->IsMoved() ) {

-      //gp_Pnt p = pLink->MediumPnt() + pLink->Move();

       gp_Pnt p = pLink->MiddlePnt() + pLink->Move();
       GetMeshDS()->MoveNode( pLink->_mediumNode, p.X(), p.Y(), p.Z());
       gp_Pnt p = pLink->MiddlePnt() + pLink->Move();
       GetMeshDS()->MoveNode( pLink->_mediumNode, p.X(), p.Y(), p.Z());

-    }
-  }
-}
+      //
+//       gp_Pnt pNew = pLink->MiddlePnt() + pLink->Move();
+//       if ( pLink->MediumPos() != SMDS_TOP_3DSPACE )
+//       {
+//         // avoid making distorted volumes near boundary
+//         SMDS_ElemIteratorPtr volIt =
+//           (*pLink)._mediumNode->GetInverseElementIterator( SMDSAbs_Volume );
+//         for ( nbVols = 0; volIt->more() && volTool.Set( volIt->next() ); ++nbVols )
+//         {
+//           vols   [ nbVols ] = volTool.Element();
+//           volSize[ nbVols ] = volTool.GetSize();
+//         }
+//         gp_Pnt pOld = pLink->MediumPnt();
+//         const_cast<SMDS_MeshNode*>( pLink->_mediumNode )->setXYZ( pNew.X(), pNew.Y(), pNew.Z() );
+//         ok = true;
+//         while ( nbVols-- && ok )
+//         {
+//           volTool.Set( vols[ nbVols ]);
+//           ok = ( volSize[ nbVols ] * volTool.GetSize() > 1e-20 ); 
+//         }
+//         if ( !ok )
+//         {
+//           const_cast<SMDS_MeshNode*>( pLink->_mediumNode )->setXYZ( pOld.X(), pOld.Y(), pOld.Z() );
+//           MSG( "Do NOT move \t" << pLink->_mediumNode->GetID()
+//                << " because of distortion of volume " << vols[ nbVols+1 ]->GetID());
+//           continue;
+//         }
+//       }
+//       GetMeshDS()->MoveNode( pLink->_mediumNode, pNew.X(), pNew.Y(), pNew.Z() );
+    }
+  }
+
+  //return;
+
+  // issue 0020982
+  // Move the apex of pyramid together with the most curved link
+
+  TIDSortedNodeSet::iterator apexIt = apexOfPyramid.begin();
+  for ( ; apexIt != apexOfPyramid.end(); ++apexIt )
+  {
+    SMESH_TNodeXYZ apex = *apexIt;

 
 

-//=======================================================================
-/*!
- * \brief Iterator on ancestors of the given type
- */
-//=======================================================================
+    gp_Vec maxMove( 0,0,0 );
+    double maxMoveSize2 = 0;

 
 

-struct TAncestorsIterator : public SMDS_Iterator<const TopoDS_Shape*>
-{
-  TopTools_ListIteratorOfListOfShape _ancIter;
-  TopAbs_ShapeEnum                   _type;
-  TopTools_MapOfShape                _encountered;
-  TAncestorsIterator( const TopTools_ListOfShape& ancestors, TopAbs_ShapeEnum type)
-    : _ancIter( ancestors ), _type( type )
-  {
-    if ( _ancIter.More() ) {
-      if ( _ancIter.Value().ShapeType() != _type ) next();
-      else _encountered.Add( _ancIter.Value() );
+    // shift of node index to get medium nodes between the base nodes
+    const int base2MediumShift = 5;
+
+    // find maximal movement of medium node
+    SMDS_ElemIteratorPtr volIt = apex._node->GetInverseElementIterator( SMDSAbs_Volume );
+    vector< const SMDS_MeshElement* > pyramids;
+    while ( volIt->more() )
+    {
+      const SMDS_MeshElement* pyram = volIt->next();
+      if ( pyram->GetEntityType() != SMDSEntity_Quad_Pyramid ) continue;
+      pyramids.push_back( pyram );
+
+      for ( int iBase = 0; iBase < apexIndex; ++iBase )
+      {
+        SMESH_TNodeXYZ medium = pyram->GetNode( iBase + base2MediumShift );
+        if ( medium._node->GetPosition()->GetTypeOfPosition() != SMDS_TOP_3DSPACE )
+        {
+          SMESH_TNodeXYZ n1 = pyram->GetNode( iBase );
+          SMESH_TNodeXYZ n2 = pyram->GetNode( ( iBase+1 ) % 4 );
+          gp_Pnt middle = 0.5 * ( n1 + n2 );
+          gp_Vec move( middle, medium );
+          double moveSize2 = move.SquareMagnitude();
+          if ( moveSize2 > maxMoveSize2 )
+            maxMove = move, maxMoveSize2 = moveSize2;
+        }
+      }

     }
     }

-  }
-  virtual bool more()
-  {
-    return _ancIter.More();
-  }
-  virtual const TopoDS_Shape* next()
-  {
-    const TopoDS_Shape* s = _ancIter.More() ? & _ancIter.Value() : 0;
-    if ( _ancIter.More() )
-      for ( _ancIter.Next();  _ancIter.More(); _ancIter.Next())
-        if ( _ancIter.Value().ShapeType() == _type && _encountered.Add( _ancIter.Value() ))
-          break;
-    return s;
-  }
-};

 
 

-//=======================================================================
-/*!
- * \brief Return iterator on ancestors of the given type
- */
-//=======================================================================
+    // move the apex
+    if ( maxMoveSize2 > 1e-20 )
+    {
+      apex += maxMove.XYZ();
+      GetMeshDS()->MoveNode( apex._node, apex.X(), apex.Y(), apex.Z());

 
 

-PShapeIteratorPtr SMESH_MesherHelper::GetAncestors(const TopoDS_Shape& shape,
-                                                   const SMESH_Mesh&   mesh,
-                                                   TopAbs_ShapeEnum    ancestorType)
-{
-  return PShapeIteratorPtr( new TAncestorsIterator( mesh.GetAncestors(shape), ancestorType));
+      // move medium nodes neighboring the apex to the middle
+      const int base2MediumShift_2 = 9;
+      for ( unsigned i = 0; i < pyramids.size(); ++i )
+        for ( int iBase = 0; iBase < apexIndex; ++iBase )
+        {
+          SMESH_TNodeXYZ         base = pyramids[i]->GetNode( iBase );
+          const SMDS_MeshNode* medium = pyramids[i]->GetNode( iBase + base2MediumShift_2 );
+          gp_XYZ middle = 0.5 * ( apex + base );
+          GetMeshDS()->MoveNode( medium, middle.X(), middle.Y(), middle.Z());
+        }
+    }
+  }

 }
 }

+

diff --git a/src/SMESH_I/SMESH_Gen_i.cxx b/src/SMESH_I/SMESH_Gen_i.cxx
index 19c78bf8fbab9dd92b2622f8e7c9e9bd03def8db..60a6136fe072ff0ea69410820041b807ebed00fc 100644 (file)
--- a/src/SMESH_I/SMESH_Gen_i.cxx
+++ b/src/SMESH_I/SMESH_Gen_i.cxx
@@ -2444,6 +2444,8 @@ SMESH::SMESH_Mesh_ptr SMESH_Gen_i::CopyMesh(SMESH::SMESH_IDSource_ptr meshPart,
     }
   }
 
     }
   }
 

+  newMeshDS->Modified();
+

   *pyDump << newMesh << " = " << this
           << ".CopyMesh( " << meshPart << ", "
           << "'" << meshName << "', "
   *pyDump << newMesh << " = " << this
           << ".CopyMesh( " << meshPart << ", "
           << "'" << meshName << "', "


@@ -3138,7 +3140,7 @@ SALOMEDS::TMPFile* SMESH_Gen_i::Save( SALOMEDS::SComponent_ptr theComponent,
 
             // groups root sub-branch
             SALOMEDS::SObject_var myGroupsBranch;
 
             // groups root sub-branch
             SALOMEDS::SObject_var myGroupsBranch;

-            for ( int i = GetNodeGroupsTag(); i <= GetVolumeGroupsTag(); i++ ) {
+            for ( int i = GetNodeGroupsTag(); i <= Get0DElementsGroupsTag(); i++ ) {

               found = gotBranch->FindSubObject( i, myGroupsBranch );
               if ( found ) {
                 char name_group[ 30 ];
               found = gotBranch->FindSubObject( i, myGroupsBranch );
               if ( found ) {
                 char name_group[ 30 ];


@@ -3150,6 +3152,8 @@ SALOMEDS::TMPFile* SMESH_Gen_i::Save( SALOMEDS::SComponent_ptr theComponent,
                   strcpy( name_group, "Groups of Faces" );
                 else if ( i == GetVolumeGroupsTag() )
                   strcpy( name_group, "Groups of Volumes" );
                   strcpy( name_group, "Groups of Faces" );
                 else if ( i == GetVolumeGroupsTag() )
                   strcpy( name_group, "Groups of Volumes" );

+                else if ( i == Get0DElementsGroupsTag() )
+                  strcpy( name_group, "Groups of 0D Elements" );

 
                 aGroup = new HDFgroup( name_group, aTopGroup );
                 aGroup->CreateOnDisk();
 
                 aGroup = new HDFgroup( name_group, aTopGroup );
                 aGroup->CreateOnDisk();


@@ -3979,7 +3983,8 @@ bool SMESH_Gen_i::Load( SALOMEDS::SComponent_ptr theComponent,
         }
       }
 
         }
       }
 

-      // try to get applied algorithms
+      // Try to get applied ALGORITHMS (mesh is not cleared by algo addition because
+      // nodes and elements are not yet put into sub-meshes)

       if ( aTopGroup->ExistInternalObject( "Applied Algorithms" ) ) {
         aGroup = new HDFgroup( "Applied Algorithms", aTopGroup );
         aGroup->OpenOnDisk();
       if ( aTopGroup->ExistInternalObject( "Applied Algorithms" ) ) {
         aGroup = new HDFgroup( "Applied Algorithms", aTopGroup );
         aGroup->OpenOnDisk();


@@ -4127,21 +4132,6 @@ bool SMESH_Gen_i::Load( SALOMEDS::SComponent_ptr theComponent,
               if ( aSubMesh->_is_nil() )
                 continue;
 
               if ( aSubMesh->_is_nil() )
                 continue;
 

-              // VSR: Get submesh data from MED convertor
-              //                  int anInternalSubmeshId = aSubMesh->GetId(); // this is not a persistent ID, it's an internal one computed from sub-shape
-              //                  if (myNewMeshImpl->_mapSubMesh.find(anInternalSubmeshId) != myNewMeshImpl->_mapSubMesh.end()) {
-              //                    if(MYDEBUG) MESSAGE("VSR - SMESH_Gen_i::Load(): loading from MED file submesh with ID = " <<
-              //                            subid << " for subshape # " << anInternalSubmeshId);
-              //                    SMESHDS_SubMesh* aSubMeshDS =
-              //                      myNewMeshImpl->_mapSubMesh[anInternalSubmeshId]->CreateSubMeshDS();
-              //                    if ( !aSubMeshDS ) {
-              //                      if(MYDEBUG) MESSAGE("VSR - SMESH_Gen_i::Load(): FAILED to create a submesh for subshape # " <<
-              //                              anInternalSubmeshId << " in current mesh!");
-              //                    }
-              //                    else
-              //                      myReader.GetSubMesh( aSubMeshDS, subid );
-              //                  }
-

               // try to get applied algorithms
               if ( aSubGroup->ExistInternalObject( "Applied Algorithms" ) ) {
                 // open "applied algorithms" HDF group
               // try to get applied algorithms
               if ( aSubGroup->ExistInternalObject( "Applied Algorithms" ) ) {
                 // open "applied algorithms" HDF group


@@ -4161,8 +4151,6 @@ bool SMESH_Gen_i::Load( SALOMEDS::SComponent_ptr theComponent,
                     aDataset->ReadFromDisk( refFromFile );
                     aDataset->CloseOnDisk();
 
                     aDataset->ReadFromDisk( refFromFile );
                     aDataset->CloseOnDisk();
 

-                    //SALOMEDS::SObject_var hypSO = myCurrentStudy->FindObjectID( refFromFile );
-                    //CORBA::Object_var hypObject = SObjectToObject( hypSO );

                     int id = atoi( refFromFile );
                     string anIOR = myStudyContext->getIORbyOldId( id );
                     if ( !anIOR.empty() ) {
                     int id = atoi( refFromFile );
                     string anIOR = myStudyContext->getIORbyOldId( id );
                     if ( !anIOR.empty() ) {


@@ -4198,8 +4186,6 @@ bool SMESH_Gen_i::Load( SALOMEDS::SComponent_ptr theComponent,
                     aDataset->ReadFromDisk( refFromFile );
                     aDataset->CloseOnDisk();
 
                     aDataset->ReadFromDisk( refFromFile );
                     aDataset->CloseOnDisk();
 

-                    //SALOMEDS::SObject_var hypSO = myCurrentStudy->FindObjectID( refFromFile );
-                    //CORBA::Object_var hypObject = SObjectToObject( hypSO );

                     int id = atoi( refFromFile );
                     string anIOR = myStudyContext->getIORbyOldId( id );
                     if ( !anIOR.empty() ) {
                     int id = atoi( refFromFile );
                     string anIOR = myStudyContext->getIORbyOldId( id );
                     if ( !anIOR.empty() ) {


@@ -4227,11 +4213,11 @@ bool SMESH_Gen_i::Load( SALOMEDS::SComponent_ptr theComponent,
 
       if(hasData) {
 
 
       if(hasData) {
 

-        // Read sub-meshes from MED
-        // -------------------------
+        // Read sub-meshes
+        // ----------------

         if(MYDEBUG) MESSAGE("Create all sub-meshes");
         bool submeshesInFamilies = ( ! aTopGroup->ExistInternalObject( "Submeshes" ));
         if(MYDEBUG) MESSAGE("Create all sub-meshes");
         bool submeshesInFamilies = ( ! aTopGroup->ExistInternalObject( "Submeshes" ));

-        if ( submeshesInFamilies )
+        if ( submeshesInFamilies ) // from MED

         {
           // old way working before fix of PAL 12992
           myReader.CreateAllSubMeshes();
         {
           // old way working before fix of PAL 12992
           myReader.CreateAllSubMeshes();


@@ -4430,7 +4416,7 @@ bool SMESH_Gen_i::Load( SALOMEDS::SComponent_ptr theComponent,
       } // if ( hasData )
 
       // try to get groups
       } // if ( hasData )
 
       // try to get groups

-      for ( int ii = GetNodeGroupsTag(); ii <= GetVolumeGroupsTag(); ii++ ) {
+      for ( int ii = GetNodeGroupsTag(); ii <= Get0DElementsGroupsTag(); ii++ ) {

         char name_group[ 30 ];
         if ( ii == GetNodeGroupsTag() )
           strcpy( name_group, "Groups of Nodes" );
         char name_group[ 30 ];
         if ( ii == GetNodeGroupsTag() )
           strcpy( name_group, "Groups of Nodes" );


@@ -4440,6 +4426,8 @@ bool SMESH_Gen_i::Load( SALOMEDS::SComponent_ptr theComponent,
           strcpy( name_group, "Groups of Faces" );
         else if ( ii == GetVolumeGroupsTag() )
           strcpy( name_group, "Groups of Volumes" );
           strcpy( name_group, "Groups of Faces" );
         else if ( ii == GetVolumeGroupsTag() )
           strcpy( name_group, "Groups of Volumes" );

+        else if ( ii == Get0DElementsGroupsTag() )
+          strcpy( name_group, "Groups of 0D Elements" );

 
         if ( aTopGroup->ExistInternalObject( name_group ) ) {
           aGroup = new HDFgroup( name_group, aTopGroup );
 
         if ( aTopGroup->ExistInternalObject( name_group ) ) {
           aGroup = new HDFgroup( name_group, aTopGroup );

diff --git a/src/SMESH_I/SMESH_Gen_i.hxx b/src/SMESH_I/SMESH_Gen_i.hxx
index bc999e3db6990738e35c234784f7a18398692f64..1087dc02368690e409fe8c594afce12cfdda37bb 100644 (file)
--- a/src/SMESH_I/SMESH_Gen_i.hxx
+++ b/src/SMESH_I/SMESH_Gen_i.hxx
@@ -461,6 +461,7 @@ public:
   static long GetEdgeGroupsTag();
   static long GetFaceGroupsTag();
   static long GetVolumeGroupsTag();
   static long GetEdgeGroupsTag();
   static long GetFaceGroupsTag();
   static long GetVolumeGroupsTag();

+  static long Get0DElementsGroupsTag();

 
   // publishing methods
   SALOMEDS::SComponent_ptr PublishComponent(SALOMEDS::Study_ptr theStudy);
 
   // publishing methods
   SALOMEDS::SComponent_ptr PublishComponent(SALOMEDS::Study_ptr theStudy);

diff --git a/src/SMESH_I/SMESH_Gen_i_1.cxx b/src/SMESH_I/SMESH_Gen_i_1.cxx
index b19214c452620f96e3d3217d437022ebcbacbd5e..41095815cd0fb06b190b9cd8645a3f86ddcfcf0a 100644 (file)
--- a/src/SMESH_I/SMESH_Gen_i_1.cxx
+++ b/src/SMESH_I/SMESH_Gen_i_1.cxx
@@ -136,6 +136,11 @@ long SMESH_Gen_i::GetVolumeGroupsTag()
   return SMESH::Tag_VolumeGroups;
 }
 
   return SMESH::Tag_VolumeGroups;
 }
 

+long SMESH_Gen_i::Get0DElementsGroupsTag()
+{
+  return SMESH::Tag_0DElementsGroups;
+}
+

 //=============================================================================
 /*!
  *  SMESH_Gen_i::CanPublishInStudy
 //=============================================================================
 /*!
  *  SMESH_Gen_i::CanPublishInStudy

diff --git a/src/SMESH_I/SMESH_Group_i.cxx b/src/SMESH_I/SMESH_Group_i.cxx
index 6c495af8ff8fe96461d9715600496c7a466dfe7f..6312b7389cfa1b6882f06e896baf8e00c5cab54c 100644 (file)
--- a/src/SMESH_I/SMESH_Group_i.cxx
+++ b/src/SMESH_I/SMESH_Group_i.cxx
@@ -508,7 +508,7 @@ void SMESH_GroupBase_i::SetColorNumber(CORBA::Long color)
   if (aGroupDS)
   {
     aGroupDS->SetColorGroup(color);
   if (aGroupDS)
   {
     aGroupDS->SetColorGroup(color);

-    TPythonDump()<<_this()<<".SetColorGroup( "<<color<<" )";
+    TPythonDump()<<_this()<<".SetColorNumber( "<<color<<" )";

   }
   MESSAGE("set color number of a group");
   return ;
   }
   MESSAGE("set color number of a group");
   return ;


@@ -535,6 +535,21 @@ SMESH::long_array* SMESH_GroupBase_i::GetMeshInfo()
     aRes[ SMESH::Entity_Node ] = aGrpDS->Extent();
   else
     SMESH_Mesh_i::CollectMeshInfo( aGrpDS->GetElements(), aRes);
     aRes[ SMESH::Entity_Node ] = aGrpDS->Extent();
   else
     SMESH_Mesh_i::CollectMeshInfo( aGrpDS->GetElements(), aRes);

+
+//   SMDS_ElemIteratorPtr it = aGrpDS->GetElements();
+//   if ( it->more() )
+//   {
+//     cout << "START" << endl;
+//     set< const SMDS_MeshElement* > nodes;
+//     const SMDS_MeshElement* e = it->next();
+//     for ( int i = 0; i < 1000000; ++i)
+//     {
+//       SMDS_ElemIteratorPtr it = e->nodesIterator();
+//       nodes.insert( e + i );
+//     }
+//     cout << "END "<< nodes.size() << endl;
+//   }
+ 

   return aRes._retn();
 }
 
   return aRes._retn();
 }
 

diff --git a/src/StdMeshers/StdMeshers_Projection_2D.cxx b/src/StdMeshers/StdMeshers_Projection_2D.cxx
index 5ed4e15dfacc41ccaba6693d38e0acd66b65ffc0..152aa72bdb4965738e2077b70eea76e4ee6d4746 100644 (file)
--- a/src/StdMeshers/StdMeshers_Projection_2D.cxx
+++ b/src/StdMeshers/StdMeshers_Projection_2D.cxx
@@ -531,6 +531,11 @@ namespace {
 
     const SMDS_MeshNode* nullNode = 0;
 
 
     const SMDS_MeshNode* nullNode = 0;
 

+    // indices of nodes to create properly oriented faces
+    int tri1 = 1, tri2 = 2, quad1 = 1, quad3 = 3;
+    if ( trsf.Form() != gp_Identity )
+      std::swap( tri1, tri2 ), std::swap( quad1, quad3 );
+

     SMESHDS_SubMesh* srcSubDS = srcMesh->GetMeshDS()->MeshElements( srcFace );
     SMDS_ElemIteratorPtr elemIt = srcSubDS->GetElements();
     vector< const SMDS_MeshNode* > tgtNodes;
     SMESHDS_SubMesh* srcSubDS = srcMesh->GetMeshDS()->MeshElements( srcFace );
     SMDS_ElemIteratorPtr elemIt = srcSubDS->GetElements();
     vector< const SMDS_MeshNode* > tgtNodes;


@@ -556,11 +561,11 @@ namespace {
         }
         tgtNodes[i] = srcN_tgtN->second;
       }
         }
         tgtNodes[i] = srcN_tgtN->second;
       }

-      // create a new face (with reversed orientation)
+      // create a new face

       switch ( nbN )
       {
       switch ( nbN )
       {

-      case 3: helper.AddFace(tgtNodes[0], tgtNodes[2], tgtNodes[1]); break;
-      case 4: helper.AddFace(tgtNodes[0], tgtNodes[3], tgtNodes[2], tgtNodes[1]); break;
+      case 3: helper.AddFace(tgtNodes[0], tgtNodes[tri1], tgtNodes[tri2]); break;
+      case 4: helper.AddFace(tgtNodes[0], tgtNodes[quad1], tgtNodes[2], tgtNodes[quad3]); break;

       }
     }
     return true;
       }
     }
     return true;

diff --git a/src/StdMeshers/StdMeshers_QuadToTriaAdaptor.cxx b/src/StdMeshers/StdMeshers_QuadToTriaAdaptor.cxx
index 2a34a0022e0dc66185046a096c7d8e96cbe8e6a2..db32a5f1a5f92274ef99ed98d721e89459bd9711 100644 (file)
--- a/src/StdMeshers/StdMeshers_QuadToTriaAdaptor.cxx
+++ b/src/StdMeshers/StdMeshers_QuadToTriaAdaptor.cxx
@@ -188,7 +188,7 @@ namespace
     typedef SMDS_StdIterator< const SMDS_MeshElement*, SMDS_ElemIteratorPtr > TStdElemIterator;
     TStdElemIterator itEnd;
 
     typedef SMDS_StdIterator< const SMDS_MeshElement*, SMDS_ElemIteratorPtr > TStdElemIterator;
     TStdElemIterator itEnd;
 

-    // shift of node index to get medium nodes corresponding to the 4 base nodes
+    // shift of node index to get medium nodes between the 4 base nodes and the apex

     const int base2MediumShift = 9;
 
     set<const SMDS_MeshNode*>::const_iterator nIt = commonApex.begin();
     const int base2MediumShift = 9;
 
     set<const SMDS_MeshNode*>::const_iterator nIt = commonApex.begin();