Salome HOME
Merge from BR_imps_2013 14/01/2014
authorvsr <vsr@opencascade.com>
Wed, 15 Jan 2014 09:41:17 +0000 (09:41 +0000)
committervsr <vsr@opencascade.com>
Wed, 15 Jan 2014 09:41:17 +0000 (09:41 +0000)
82 files changed:
doc/salome/examples/defining_hypotheses_ex01.py
doc/salome/examples/defining_hypotheses_ex08.py
doc/salome/gui/SMESH/images/a-geometric1d.png [new file with mode: 0644]
doc/salome/gui/SMESH/images/viscous_layers_hyp.png
doc/salome/gui/SMESH/input/1d_meshing_hypo.doc
doc/salome/gui/SMESH/input/about_hypo.doc
doc/salome/gui/SMESH/input/additional_hypo.doc
doc/salome/gui/SMESH/input/tui_defining_hypotheses.doc
idl/SMESH_BasicHypothesis.idl
idl/SMESH_Mesh.idl
resources/StdMeshers.xml.in
src/SMESH/SMESH_Algo.cxx
src/SMESH/SMESH_MesherHelper.cxx
src/SMESH/SMESH_Pattern.cxx
src/SMESH/SMESH_Pattern.hxx
src/SMESHGUI/SMESHGUI_Hypotheses.cxx
src/SMESHGUI/SMESHGUI_HypothesesUtils.cxx
src/SMESHGUI/SMESHGUI_HypothesesUtils.h
src/SMESHGUI/SMESHGUI_MeshDlg.cxx
src/SMESHGUI/SMESHGUI_MeshDlg.h
src/SMESHGUI/SMESHGUI_MeshOp.cxx
src/SMESHGUI/SMESHGUI_MeshOp.h
src/SMESHGUI/SMESH_msg_en.ts
src/SMESHGUI/SMESH_msg_fr.ts
src/SMESHUtils/SMESH_Block.cxx
src/SMESHUtils/SMESH_Block.hxx
src/SMESHUtils/SMESH_MeshAlgos.cxx
src/SMESHUtils/SMESH_MeshAlgos.hxx
src/SMESHUtils/SMESH_TypeDefs.hxx
src/SMESH_I/SMESH_2smeshpy.cxx
src/SMESH_I/SMESH_DumpPython.cxx
src/SMESH_I/SMESH_Mesh_i.cxx
src/SMESH_I/SMESH_Mesh_i.hxx
src/SMESH_I/SMESH_PythonDump.hxx
src/SMESH_SWIG/StdMeshersBuilder.py
src/SMESH_SWIG/smeshBuilder.py
src/SMESH_SWIG/smesh_algorithm.py
src/StdMeshers/CMakeLists.txt
src/StdMeshers/StdMeshers_CartesianParameters3D.cxx
src/StdMeshers/StdMeshers_CartesianParameters3D.hxx
src/StdMeshers/StdMeshers_Cartesian_3D.cxx
src/StdMeshers/StdMeshers_FaceSide.cxx
src/StdMeshers/StdMeshers_FaceSide.hxx
src/StdMeshers/StdMeshers_Geometric1D.cxx [new file with mode: 0644]
src/StdMeshers/StdMeshers_Geometric1D.hxx [new file with mode: 0644]
src/StdMeshers/StdMeshers_Hexa_3D.cxx
src/StdMeshers/StdMeshers_Prism_3D.cxx
src/StdMeshers/StdMeshers_Projection_2D.cxx
src/StdMeshers/StdMeshers_Propagation.cxx
src/StdMeshers/StdMeshers_Propagation.hxx
src/StdMeshers/StdMeshers_QuadrangleParams.cxx
src/StdMeshers/StdMeshers_QuadrangleParams.hxx
src/StdMeshers/StdMeshers_Quadrangle_2D.cxx
src/StdMeshers/StdMeshers_Quadrangle_2D.hxx
src/StdMeshers/StdMeshers_Regular_1D.cxx
src/StdMeshers/StdMeshers_Regular_1D.hxx
src/StdMeshers/StdMeshers_Reversible1D.cxx [new file with mode: 0644]
src/StdMeshers/StdMeshers_Reversible1D.hxx [new file with mode: 0644]
src/StdMeshers/StdMeshers_UseExisting_1D2D.cxx
src/StdMeshers/StdMeshers_ViscousLayers.cxx
src/StdMeshersGUI/StdMeshersGUI_CartesianParamCreator.cxx
src/StdMeshersGUI/StdMeshersGUI_CartesianParamCreator.h
src/StdMeshersGUI/StdMeshersGUI_StdHypothesisCreator.cxx
src/StdMeshersGUI/StdMeshers_images.ts
src/StdMeshersGUI/StdMeshers_msg_en.ts
src/StdMeshersGUI/StdMeshers_msg_fr.ts
src/StdMeshers_I/CMakeLists.txt
src/StdMeshers_I/StdMeshers_CartesianParameters3D_i.cxx
src/StdMeshers_I/StdMeshers_CartesianParameters3D_i.hxx
src/StdMeshers_I/StdMeshers_Geometric1D_i.cxx [new file with mode: 0644]
src/StdMeshers_I/StdMeshers_Geometric1D_i.hxx [new file with mode: 0644]
src/StdMeshers_I/StdMeshers_ObjRefUlils.cxx
src/StdMeshers_I/StdMeshers_ObjRefUlils.hxx
src/StdMeshers_I/StdMeshers_Propagation_i.cxx
src/StdMeshers_I/StdMeshers_Propagation_i.hxx
src/StdMeshers_I/StdMeshers_QuadrangleParams_i.cxx
src/StdMeshers_I/StdMeshers_QuadrangleParams_i.hxx
src/StdMeshers_I/StdMeshers_Reversible1D_i.cxx [new file with mode: 0644]
src/StdMeshers_I/StdMeshers_Reversible1D_i.hxx [new file with mode: 0644]
src/StdMeshers_I/StdMeshers_ViscousLayers_i.cxx
src/StdMeshers_I/StdMeshers_ViscousLayers_i.hxx
src/StdMeshers_I/StdMeshers_i.cxx

index 4bb75c5..1d5d281 100644 (file)
@@ -1,12 +1,11 @@
-# Arithmetic 1D
+# Arithmetic 1D and Geometric Progression
 
 import salome
 salome.salome_init()
-import GEOM
+
 from salome.geom import geomBuilder
 geompy = geomBuilder.New(salome.myStudy)
 
-import SMESH, SALOMEDS
 from salome.smesh import smeshBuilder
 smesh =  smeshBuilder.New(salome.myStudy)
 
@@ -21,12 +20,20 @@ hexa = smesh.Mesh(box, "Box : hexahedrical mesh")
 algo1D = hexa.Segment()
 
 # optionally reverse node distribution on certain edges
-allEdges = geompy.SubShapeAllSortedIDs( box, geompy.ShapeType["EDGE"])
+allEdges = geompy.SubShapeAllSorted( box, geompy.ShapeType["EDGE"])
 reversedEdges = [ allEdges[0], allEdges[4] ]
 
 # define "Arithmetic1D" hypothesis to cut all edges in several segments with increasing arithmetic length 
 algo1D.Arithmetic1D(1, 4, reversedEdges)
 
+# define "Geometric Progression" hypothesis on one edge to cut this edge in segments with length increasing by 20% starting from 1
+gpAlgo = hexa.Segment( allEdges[1] )
+gpAlgo.GeometricProgression( 1, 1.2 )
+
+# propagate distribution of nodes computed using "Geometric Progression" to parallel edges
+gpAlgo.PropagationOfDistribution() 
+
+
 # create a quadrangle 2D algorithm for faces
 hexa.Quadrangle()
 
index 34fcbc3..3c877f3 100644 (file)
@@ -11,7 +11,8 @@ from salome.smesh import smeshBuilder
 smesh =  smeshBuilder.New(salome.myStudy)
 
 # create a box
-box = geompy.MakeBoxDXDYDZ(10., 10., 10.)
+base = geompy.MakeSketcher("Sketcher:F 0 0:TT 10 0:TT 20 10:TT 0 10:WF", theName="F")
+box  = geompy.MakePrismDXDYDZ( base, 0,0,10 )
 geompy.addToStudy(box, "Box")
 
 # get one edge of the box to put local hypothesis on
@@ -20,7 +21,7 @@ EdgeX = geompy.GetEdgeNearPoint(box, p5)
 geompy.addToStudyInFather(box, EdgeX, "Edge [0,0,0 - 10,0,0]")
 
 # create a hexahedral mesh on the box
-hexa = smesh.Mesh(box, "Box : hexahedrical mesh")
+hexa = smesh.Mesh(box, "Propagation of hypothesis")
 
 # set global algorithms and hypotheses
 algo1D = hexa.Segment()
@@ -28,15 +29,37 @@ hexa.Quadrangle()
 hexa.Hexahedron()
 algo1D.NumberOfSegments(4)
 
-# create a sub-mesh with local 1D hypothesis and propagation
+# create a sub-mesh with local 1D hypothesis and "Propagation of 1D Hypothesis"
 algo_local = hexa.Segment(EdgeX)
 
 # define "Arithmetic1D" hypothesis to cut an edge in several segments with increasing length
 algo_local.Arithmetic1D(1, 4)
 
-# define "Propagation" hypothesis that propagates all other 1D hypotheses
-# from all edges on the opposite side of a face in case of quadrangular faces
+# define "Propagation" hypothesis that propagates "Arithmetic1D" hypothesis
+# from 'EdgeX' on opposite sides of all quadilateral faces
 algo_local.Propagation()
 
-# compute the mesh
+# compute the mesh which contains prisms
 hexa.Compute()
+
+
+# create another mesh on the box
+mesh = smesh.Mesh(box, "Propagation of distribution of nodes")
+
+# set global algorithms and hypotheses
+algo1D = mesh.Segment()
+mesh.Quadrangle()
+mesh.Hexahedron()
+algo1D.NumberOfSegments(4)
+
+# create a sub-mesh with local 1D hypothesis and "Propagation of Node Distribution"
+algo_local = mesh.Segment(EdgeX)
+algo_local.Arithmetic1D(1, 4)
+
+# define "Propagation Of Distribution" hypothesis that propagates
+# distribution of nodes generated by "Arithmetic1D" hypothesis
+# from 'EdgeX' on opposite sides of all quadilateral faces
+algo_local.PropagationOfDistribution()
+
+# compute the mesh which contains hexahedra only
+mesh.Compute()
diff --git a/doc/salome/gui/SMESH/images/a-geometric1d.png b/doc/salome/gui/SMESH/images/a-geometric1d.png
new file mode 100644 (file)
index 0000000..a60be94
Binary files /dev/null and b/doc/salome/gui/SMESH/images/a-geometric1d.png differ
index 9c99303..7178191 100644 (file)
Binary files a/doc/salome/gui/SMESH/images/viscous_layers_hyp.png and b/doc/salome/gui/SMESH/images/viscous_layers_hyp.png differ
index 470fe0c..b992bc8 100644 (file)
@@ -6,6 +6,7 @@
 <ul>
 <li>\ref adaptive_1d_anchor "Adaptive"</li>
 <li>\ref arithmetic_1d_anchor "Arithmetic 1D"</li>
+<li>\ref geometric_1d_anchor "Geometric Progression"</li>
 <li>\ref average_length_anchor "Local Length"</li>
 <li>\ref max_length_anchor "Max Size"</li>
 <li>\ref deflection_1d_anchor "Deflection 1D"</li>
@@ -55,7 +56,30 @@ picking them in the 3D viewer or by selecting the edges or groups of edges in th
 \image html b-ithmetic1d.png "Arithmetic 1D hypothesis - the size of mesh elements gradually increases"
 
 <b>See Also</b> a sample TUI Script of a 
-\ref tui_1d_arithmetic "Defining Arithmetic 1D hypothesis" operation.  
+\ref tui_1d_arithmetic "Defining Arithmetic 1D and Geometric Progression hypothesis" operation.  
+
+<br>
+\anchor geometric_1d_anchor
+<h2>Geometric Progression hypothesis</h2>
+
+<b>Geometric Progression</b> hypothesis allows to split edges into
+segments with a length that changes in geometric progression (Lk =
+Lk-1 * d) beginning from a given starting length and with a given
+common ratio.
+
+The direction of the splitting is defined by the orientation of the
+underlying geometrical edge. <b>"Reverse Edges"</b> list box allows to
+specify the edges for which the splitting should be made in the
+direction opposing to their orientation. This list box is enabled only
+if the geometry object is selected for the meshing. In this case the
+user can select edges to be reversed either directly picking them in
+the 3D viewer or by selecting the edges or groups of edges in the
+Object Browser. 
+
+\image html a-geometric1d.png
+
+<b>See Also</b> a sample TUI Script of a 
+\ref tui_1d_arithmetic "Defining Arithmetic 1D and Geometric Progression hypothesis" operation.  
 
 <br>
 \anchor deflection_1d_anchor
index 09768cc..65d5c7a 100644 (file)
@@ -20,6 +20,7 @@ In \b MESH there are the following Basic Hypotheses:
 <li>\ref max_length_anchor "Max Size"</li>
 <li>\ref adaptive_1d_anchor "Adaptive"</li>
 <li>\ref arithmetic_1d_anchor "Arithmetic 1D"</li>
+<li>\ref geometric_1d_anchor "Geometric 1D"</li>
 <li>\ref start_and_end_length_anchor "Start and end length"</li>
 <li>\ref deflection_1d_anchor "Deflection 1D"</li>
 <li>\ref automatic_length_anchor "Automatic Length"</li>
@@ -41,6 +42,7 @@ There also exist
 with other hypotheses:
 <ul>
 <li>\ref propagation_anchor "Propagation of 1D Hypothesis on opposite edges"</li>
+<li>\ref propagofdistribution_anchor "Propagation of Node Distribution on Opposite Edges"</li>
 <li>\ref viscous_layers_anchor "Viscous layers"</li>
 <li>\ref quadratic_mesh_anchor "Quadratic mesh"</li>
 <li>\ref non_conform_allowed_anchor "Non conform mesh allowed"</li>
index 938c787..adccc49 100644 (file)
@@ -38,6 +38,19 @@ has been locally defined on the opposite edge.
 <br><b>See Also</b> a sample TUI Script of a 
 \ref tui_propagation "Propagation hypothesis" operation
 
+\anchor propagofdistribution_anchor
+<h2>Propagation of Node Distribution on Opposite Edges</h2>
+
+<b>Propagation of Node Distribution on Opposite Edges</b> allows to propagate
+distribution of nodes onto an opposite edge. If a local hypothesis and
+propagation are defined on an edge of a quadrangular face, the
+opposite edge will have the same number of nodes and the same
+relations between segment lengths, unless another hypothesis
+has been locally defined on the opposite edge.
+<br><b>See Also</b> a sample TUI Script of a 
+\ref tui_propagation "Propagation hypothesis" operation
+
 \anchor quadrangle_preference_anchor
 <h2>Quadrangle Preference</h2>
 
@@ -68,29 +81,29 @@ computations.
 <li><b>Number of layers</b> - defines the number of element layers.</li>
 <li><b>Stretch factor</b> - defines the growth factor of element height
   from the mesh boundary inwards.</li>
-<li><b>Specified Edges are</b> - defines how the shapes specified by
+<li><b>Specified Faces/Edges are</b> - defines how the shapes specified by
   the next parameter are used.
-<li><b>Faces without layers</b> and <b>Edges with/without layers</b> - 
-  in the 3D case it defines geometrical faces on which element layers
-  should not be constructed; in the 2D case it defines geometrical edges
-  on which element layers either should be or should not be
-  constructed, depending on the value of the previous parameter
-  (<b>Specified Edges are</b>). 
+<li><b> Faces/Edges with/without layers</b> - 
+  defines geometrical faces or edges on which element layers
+  either should be or should not be constructed, depending on the
+  value of the previous parameter (<b>Specified Faces/Edges are</b>). 
+  Faces (or edges) can be selected either in the Object Browser or in
+  the VTK Viewer.
   \note A mesh shown in the 3D Viewer can prevent selection of faces
   and edges, just hide the mesh to avoid this. To avoid a long wait when a
   geometry with many faces (or edges) is displayed, the number of faces
   (edges) shown at a time is limited by the value of "Sub-shapes
   preview chunk size" preference (in Preferences/Mesh/General tab).<br>
 
-  Whatever shapes are specified by this
-  parameter, the element layers are not constructed on geometrical
-  faces shared by several solids in 3D case and edges shared by
-  several faces in 2D case. In other words the element layers can be
-  constructed on boundary faces and edges, and are not constructed on
-  internal faces and edges. There is an exception to this rule in 2D
-  case: if "Viscous Layers 2D" hypothesis is assigned to a sub-mesh,
-  the element layers can be constructed on boundary edges of the shape
-  of this sub-mesh.
+  If faces/edges without layers are specified, the element layers are
+  not constructed on geometrical faces shared by several solids in 3D
+  case and edges shared by several faces in 2D case. In other words,
+  in this mode the element layers can be constructed on boundary faces
+  and edges only, and are not constructed on internal faces and
+  edges. There is an exception to this rule: if a hypothesis is
+  assigned to a sub-mesh, the element layers can be constructed on
+  boundary faces/edges of the shape of this sub-mesh, at same time
+  possibly being internal faces/edges within the whole model.
   \image html viscous_layers_on_submesh.png 2D viscous layers constructed on boundary edges of a sub-mesh on a disk face.
 
 </li>
@@ -101,5 +114,4 @@ computations.
 <br><b>See also</b> a sample TUI script of a \ref tui_viscous_layers
 "Viscous layers construction".
 
-
 */
index 333017d..863bf76 100644 (file)
@@ -9,6 +9,7 @@ This page provides example codes of \ref tui_defining_meshing_algos
   <ul>
     <li>\ref tui_1d_adaptive "Adaptive 1D" hypothesis</li>
     <li>\ref tui_1d_arithmetic "Arithmetic 1D" hypothesis</li>
+    <li>\ref tui_1d_arithmetic "Geometric Progression" hypothesis</li>
     <li>\ref tui_deflection_1d "Deflection 1D and Number of Segments" hypotheses</li>
     <li>\ref tui_start_and_end_length "Start and End Length" hypotheses</li>
     <li>\ref tui_average_length "Local Length"</li>
@@ -44,7 +45,7 @@ This page provides example codes of \ref tui_defining_meshing_algos
 
 <br>
 \anchor tui_1d_arithmetic
-<h3>Arithmetic 1D</h3>
+<h3>Arithmetic 1D and Geometric Progression</h3>
 \tui_script{defining_hypotheses_ex01.py}
 
 <br>
index 798475b..186b066 100644 (file)
@@ -128,9 +128,35 @@ module StdMeshers
   };
 
   /*!
+   * Common inteface of 1D hypotheses that can be reversed
+   */
+  interface Reversible1D
+  {
+    /*!
+     * Set list of edges to reverse
+     */
+    void SetReversedEdges( in SMESH::long_array list );
+    
+    /*!
+     * Returns list of edges to reverse
+     */
+    SMESH::long_array GetReversedEdges();
+    
+    /*!
+     * Set entry of the main object
+     */
+    void SetObjectEntry( in string entry );
+    
+    /*!
+     * Get the entry of the main object
+     */
+    string GetObjectEntry();
+  };
+
+  /*!
    * StdMeshers_NumberOfSegments: interface of "Nb. Segments" hypothesis
    */
-  interface StdMeshers_NumberOfSegments : SMESH::SMESH_Hypothesis
+  interface StdMeshers_NumberOfSegments : SMESH::SMESH_Hypothesis, Reversible1D
   {
     /*!
      * Builds and returns point distribution according to passed density function
@@ -209,32 +235,12 @@ module StdMeshers
      */
     long ConversionMode()
       raises (SALOME::SALOME_Exception);
-
-    /*!
-     * Set list of edges to reverse
-     */
-    void SetReversedEdges( in SMESH::long_array list );
-    
-    /*!
-     * Returns list of edges to reverse
-     */
-    SMESH::long_array GetReversedEdges();
-    
-    /*!
-     * Set entry of the main object
-     */
-    void SetObjectEntry( in string entry );
-    
-    /*!
-     * Get the entry of the main object
-     */
-    string GetObjectEntry();
   };
 
   /*!
    * StdMeshers_Arithmetic1D: interface of "Arithmetic 1D" hypothesis
    */
-  interface StdMeshers_Arithmetic1D : SMESH::SMESH_Hypothesis
+  interface StdMeshers_Arithmetic1D : SMESH::SMESH_Hypothesis, Reversible1D
   {
     /*!
      * Sets <start segment length> or <end segment length> parameter value
@@ -260,26 +266,36 @@ module StdMeshers
      * Returns <start segment length> or <end segment length> parameter value
      */
     double GetLength(in boolean isStartLength);
-    
+
+  };
+
+  /*!
+   * StdMeshers_Arithmetic1D: interface of "Geometric 1D" hypothesis
+   */
+  interface StdMeshers_Geometric1D : SMESH::SMESH_Hypothesis, Reversible1D
+  {
     /*!
-     * Set list of edges to reverse
+     * Sets length of the first segment
      */
-    void SetReversedEdges( in SMESH::long_array list );
-    
+    void SetStartLength(in double length) 
+      raises (SALOME::SALOME_Exception);
+
     /*!
-     * Returns list of edges to reverse
+     * Sets value of Common Ratio
      */
-    SMESH::long_array GetReversedEdges();
-    
+    void SetCommonRatio(in double factor)
+      raises (SALOME::SALOME_Exception);
+
     /*!
-     * Set entry of the main object
+     * Returns length of the first segment
      */
-    void SetObjectEntry( in string entry );
-    
+    double GetStartLength();
+
     /*!
-     * Get the entry of the main object
+     * Returns value of Common Ratio
      */
-    string GetObjectEntry();
+    double GetCommonRatio();
+
   };
 
   /*!
@@ -319,7 +335,7 @@ module StdMeshers
   /*!
    * StdMeshers_StartEndLength: interface of "Start and End Length" hypothesis
    */
-  interface StdMeshers_StartEndLength : SMESH::SMESH_Hypothesis
+  interface StdMeshers_StartEndLength : SMESH::SMESH_Hypothesis, Reversible1D
   {
     /*!
      * Sets <start segment length> or <end segment length> parameter value
@@ -346,25 +362,6 @@ module StdMeshers
      */
     double GetLength(in boolean isStartLength);
 
-    /*!
-     * Set list of edges to reverse
-     */
-    void SetReversedEdges( in SMESH::long_array list );
-    
-    /*!
-     * Returns list of edges to reverse
-     */
-    SMESH::long_array GetReversedEdges();
-    
-    /*!
-     * Set entry of the main object
-     */
-    void SetObjectEntry( in string entry );
-    
-    /*!
-     * Get the entry of the main object
-     */
-    string GetObjectEntry();
   };
 
 
@@ -388,7 +385,7 @@ module StdMeshers
   /*!
    * StdMeshers_FixedPoints1D: interface of "Fixed points 1D" hypothesis
    */
-  interface StdMeshers_FixedPoints1D : SMESH::SMESH_Hypothesis
+  interface StdMeshers_FixedPoints1D : SMESH::SMESH_Hypothesis, Reversible1D
   {
     /*!
      * Sets some points on edge using parameter on curve from 0 to 1
@@ -410,26 +407,7 @@ module StdMeshers
      * Returns list of numbers of segments
      */
     SMESH::long_array GetNbSegments();
-    
-    /*!
-     * Set list of edges to reverse
-     */
-    void SetReversedEdges( in SMESH::long_array list );
-    
-    /*!
-     * Returns list of edges to reverse
-     */
-    SMESH::long_array GetReversedEdges();
-    
-    /*!
-     * Set entry of the main object
-     */
-    void SetObjectEntry( in string entry );
-    
-    /*!
-     * Get the entry of the main object
-     */
-    string GetObjectEntry();
+
   };
 
   /*!
@@ -483,7 +461,8 @@ module StdMeshers
   };
 
   /*!
-   * StdMeshers_Propagation: interface of "Propagation" hypothesis.
+   * StdMeshers_Propagation: interface of "Propagation of 1D Hyp. on
+   * Opposite Edges" hypothesis.
    * Presence of this hypothesis on any edge propagates any other 1D
    * hypothesis from this edge on all edges, opposite to it.
    * It concerns only edges of quadrangle faces.
@@ -493,6 +472,17 @@ module StdMeshers
   };
 
   /*!
+   * StdMeshers_Propagation: interface of "Propagation of Node
+   * Distribution on Opposite Edges" hypothesis.
+   * Presence of this hypothesis on any edge propagates distribution of nodes
+   * from this edge on all edges, opposite to it.
+   * It concerns only edges of quadrangle faces.
+   */
+  interface StdMeshers_PropagOfDistribution : SMESH::SMESH_Hypothesis
+  {
+  };
+
+  /*!
    * StdMeshers_QuadranglePreference: interface of "QuadranglePreference" hypothesis.
    * This hypothesis is used by StdMeshers_Quadrangle_2D algorithm.
    * Presence of this hypothesis forces construction of quadrangles if the number
@@ -807,6 +797,22 @@ module StdMeshers
      * Get the type of quadrangulation
      */
     QuadType GetQuadType();
+
+    /*!
+     * Set positions of enforced nodes
+     */
+    void SetEnforcedNodes(in GEOM::ListOfGO vertices, in SMESH::nodes_array points)
+      raises (SALOME::SALOME_Exception);
+
+    /*!
+     * Returns positions of enforced nodes
+     */
+    void GetEnforcedNodes(out GEOM::ListOfGO vertices, out SMESH::nodes_array points);
+
+    /*!
+     * Returns entries of shapes defining enforced nodes
+     */
+    SMESH::string_array GetEnfVertices();
   };
 
   /*!
@@ -866,6 +872,14 @@ module StdMeshers
     SMESH::long_array GetIgnoreFaces();
 
     /*!
+     * Set faces either to exclude from treatment or to make the Viscous Layers on.
+     */
+    void SetFaces(in SMESH::long_array faceIDs,
+                  in boolean           toIgnore) raises (SALOME::SALOME_Exception);
+    SMESH::long_array GetFaces();
+    boolean           GetIsToIgnoreFaces();
+
+    /*!
      * Set total thickness of layers of prisms
      */
     void SetTotalThickness(in double thickness) raises (SALOME::SALOME_Exception);
@@ -936,7 +950,7 @@ module StdMeshers
     /*!
      * Set size threshold. A polyhedral cell got by cutting an initial
      * hexahedron by geometry boundary is considered small and is removed if
-     * it's size is \athreshold times less than the size of the initial hexahedron. 
+     * it's size is \a threshold times less than the size of the initial hexahedron. 
      * threshold must be > 1.0
      */
     void SetSizeThreshold(in double threshold) raises (SALOME::SALOME_Exception);
@@ -971,6 +985,13 @@ module StdMeshers
     void GetGridSpacing(out SMESH::string_array spaceFunctions,
                         out SMESH::double_array internalPoints,
                         in short                axis) raises (SALOME::SALOME_Exception);
+    /*!
+     * Enables implementation of geometrical edges into the mesh. If this feature
+     * is disabled, sharp edges of the shape are lost ("smoothed") in the mesh if
+     * they don't coincide with the grid lines
+     */
+    void SetToAddEdges(in boolean toAdd);
+    boolean GetToAddEdges();
 
     /*!
      * \brief Computes node coordinates by spacing functions
@@ -978,13 +999,15 @@ module StdMeshers
      *  \param x1 - upper coordinate
      *  \param spaceFuns - space functions
      *  \param points - internal points
-     *  \param coords - the computed coordinates
+     *  \param axisName - e.g. "X"
+     *  \return the computed coordinates
      */
     SMESH::double_array ComputeCoordinates(in double              x0,
                                            in double              x1,
                                            in SMESH::string_array spaceFuns,
                                            in SMESH::double_array points,
-                                           in string              axisName ) raises (SALOME::SALOME_Exception);
+                                           in string              axisName ) 
+      raises (SALOME::SALOME_Exception);
   };
 
   /*!
index fb85b15..ff60e68 100644 (file)
@@ -926,6 +926,11 @@ module SMESH
     long_array GetElemFaceNodes(in long elemId, in short faceIndex);
 
     /*!
+     * Returns three components of normal of given mesh face (or an empty array in KO case)
+     */
+    double_array GetFaceNormal(in long faceId);
+
+    /*!
      * Returns an element based on all given nodes.
      */
     long FindElementByNodes(in long_array nodes);
index d80b4e7..d6f6527 100644 (file)
                 icon-id  ="mesh_hypo_length.png"
                 dim      ="1"/>
 
+    <hypothesis type     ="GeometricProgression"
+                label-id ="Geometric Progression"
+                icon-id  ="mesh_hypo_length.png"
+                dim      ="1"/>
+
     <hypothesis type     ="FixedPoints1D"
                 label-id ="Fixed Points 1D"
                 icon-id  ="mesh_hypo_length.png"
                 dim      ="1"
                 auxiliary="true"/>
 
+    <hypothesis type     ="PropagOfDistribution"
+                label-id ="Propagation of Node Distribution on Opposite Edges"
+                icon-id  ="mesh_hypo_length.png"
+                dim      ="1"
+                auxiliary="true"/>
+
     <hypothesis type     ="AutomaticLength"
                 label-id ="Automatic Length"
                 icon-id  ="mesh_hypo_length.png"
     <algorithm type     ="Regular_1D"
               label-id ="Wire Discretisation"
               icon-id  ="mesh_algo_regular.png"
-               hypos    ="Adaptive1D,LocalLength,MaxLength,Arithmetic1D,StartEndLength,NumberOfSegments,Deflection1D,AutomaticLength,FixedPoints1D"
-               opt-hypos="Propagation,QuadraticMesh"
+               hypos    ="Adaptive1D,LocalLength,MaxLength,Arithmetic1D,GeometricProgression,StartEndLength,NumberOfSegments,Deflection1D,AutomaticLength,FixedPoints1D"
+               opt-hypos="Propagation,PropagOfDistribution,QuadraticMesh"
                input    ="VERTEX"
                output   ="EDGE"
                dim      ="1">
         <hypo>LocalLength=LocalLength(SetLength(1),,SetPrecision(1))</hypo>
         <hypo>MaxLength=MaxSize(SetLength(1))</hypo>
         <hypo>Arithmetic1D=Arithmetic1D(SetStartLength(),SetEndLength(),SetReversedEdges())</hypo>
+        <hypo>GeometricProgression=GeometricProgression(SetStartLength(),SetCommonRatio(),SetReversedEdges())</hypo>
         <hypo>StartEndLength=StartEndLength(SetStartLength(),SetEndLength(),SetReversedEdges())</hypo>
         <hypo>Deflection1D=Deflection1D(SetDeflection())</hypo>
         <hypo>Adaptive1D=Adaptive(SetMinSize(),SetMaxSize(),SetDeflection())</hypo>
         <hypo>AutomaticLength=AutomaticLength(SetFineness())</hypo>
         <hypo>FixedPoints1D=FixedPoints1D(SetPoints(),SetNbSegments(),SetReversedEdges())</hypo>
         <hypo>Propagation=Propagation()</hypo>
+        <hypo>PropagOfDistribution=PropagationOfDistribution()</hypo>
         <hypo>QuadraticMesh=QuadraticMesh()</hypo>
       </python-wrap>
     </algorithm>
     <algorithm type     ="CompositeSegment_1D"
               label-id ="Composite Side Discretisation"
               icon-id  ="mesh_algo_regular.png"
-               hypos    ="Adaptive1D,LocalLength,MaxLength,Arithmetic1D,StartEndLength,NumberOfSegments,Deflection1D,AutomaticLength,FixedPoints1D"
-               opt-hypos="Propagation,QuadraticMesh"
+               hypos    ="Adaptive1D,LocalLength,MaxLength,Arithmetic1D,GeometricProgression,StartEndLength,NumberOfSegments,Deflection1D,AutomaticLength,FixedPoints1D"
+               opt-hypos="Propagation,PropagOfDistribution,QuadraticMesh"
                input    ="VERTEX"
                output   ="EDGE"
                dim      ="1">
         <hypo>LocalLength=LocalLength(SetLength(), ,SetPrecision())</hypo>
         <hypo>MaxLength=MaxSize(SetLength())</hypo>
         <hypo>Arithmetic1D=Arithmetic1D(SetStartLength(),SetEndLength(),SetReversedEdges())</hypo>
+        <hypo>GeometricProgression=GeometricProgression(SetStartLength(),SetCommonRatio(),SetReversedEdges())</hypo>
         <hypo>StartEndLength=StartEndLength(SetStartLength(),SetEndLength(),SetReversedEdges())</hypo>
         <hypo>Deflection1D=Deflection1D(SetDeflection())</hypo>
         <hypo>Adaptive1D=Adaptive(SetMinSize(),SetMaxSize(),SetDeflection())</hypo>
         <hypo>AutomaticLength=AutomaticLength(SetFineness())</hypo>
         <hypo>FixedPoints1D=FixedPoints1D(SetPoints(),SetNbSegments(),SetReversedEdges())</hypo>
         <hypo>Propagation=Propagation()</hypo>
+        <hypo>PropagOfDistribution=PropagationOfDistribution()</hypo>
         <hypo>QuadraticMesh=QuadraticMesh()</hypo>
       </python-wrap>
     </algorithm>
                label-id ="Hexahedron (i,j,k)"
                icon-id  ="mesh_algo_hexa.png"
                input    ="QUAD"
+               output   ="HEXA,PENTA"
               need-geom="false"
                opt-hypos="ViscousLayers"
                dim      ="3">
                label-id="3D Extrusion"
                icon-id ="mesh_algo_hexa.png"
                input   ="QUAD,TRIA"
+               output  ="HEXA,PENTA,OCTA,POLYHEDRON"
                dim     ="3">
       <python-wrap>
         <algo>Prism_3D=Prism()</algo>
                icon-id ="mesh_algo_hexa.png"
                hypos   ="NumberOfLayers, LayerDistribution"
                input   ="QUAD,TRIA"
+               output  ="HEXA,PENTA,OCTA,POLYHEDRON"
                dim     ="3">
       <python-wrap>
         <algo>RadialPrism_3D=Prism('RadialPrism_3D')</algo>
                icon-id ="mesh_algo_quad.png"
                hypos   ="NumberOfLayers2D, LayerDistribution2D"
                input   ="EDGE"
-               output  ="QUAD,TRIA"
+               output  ="QUAD"
                dim     ="2">
       <python-wrap>
         <algo>RadialQuadrangle_1D2D=Quadrangle(algo=smeshBuilder.RADIAL_QUAD)</algo>
                icon-id          ="mesh_algo_hexa.png"
                hypos            ="CartesianParameters3D"
                support-submeshes="false"
+               output           ="HEXA"
                dim              ="3">
       <python-wrap>
         <algo>Cartesian_3D=BodyFitted()</algo>
index 19a1e4b..8beda03 100644 (file)
@@ -410,7 +410,7 @@ bool SMESH_Algo::GetSortedNodesOnEdge(const SMESHDS_Mesh*                   theM
     return false;
 
   SMESHDS_SubMesh * eSubMesh = theMesh->MeshElements( theEdge );
-  if ( !eSubMesh || !eSubMesh->GetElements()->more() )
+  if ( !eSubMesh || ( eSubMesh->NbElements()==0 && eSubMesh->NbNodes() == 0))
     return false; // edge is not meshed
 
   int nbNodes = 0;
index b7d8759..06f2470 100644 (file)
@@ -896,7 +896,7 @@ bool SMESH_MesherHelper::CheckNodeU(const TopoDS_Edge&   E,
                                     const bool           force,
                                     double               distXYZ[4]) const
 {
-  int shapeID = n->getshapeId();
+  int  shapeID = n->getshapeId();
   bool infinit = Precision::IsInfinite( u );
   if ( force || toCheckPosOnShape( shapeID ) || infinit )
   {
index 63e6f36..5d0f5e5 100644 (file)
@@ -2641,10 +2641,10 @@ bool SMESH_Pattern::Apply (const SMDS_MeshFace* theFace,
 
   list< const SMDS_MeshNode* > nodes;
   list< const SMDS_MeshNode* >::iterator n = nodes.end();
-  SMDS_ElemIteratorPtr noIt = theFace->nodesIterator();
+  SMDS_NodeIteratorPtr noIt = theFace->nodeIterator();
   int iSub = 0;
   while ( noIt->more() && iSub < nbFaceNodes ) {
-    const SMDS_MeshNode* node = smdsNode( noIt->next() );
+    const SMDS_MeshNode* node = noIt->next();
     nodes.push_back( node );
     if ( iSub++ == theNodeIndexOnKeyPoint1 )
       n = --nodes.end();
@@ -2661,7 +2661,7 @@ bool SMESH_Pattern::Apply (const SMDS_MeshFace* theFace,
   list< gp_XYZ > xyzList;
   myOrderedNodes.resize( nbFaceNodes );
   for ( iSub = 0, n = nodes.begin(); n != nodes.end(); ++n ) {
-    xyzList.push_back( gp_XYZ( (*n)->X(), (*n)->Y(), (*n)->Z() ));
+    xyzList.push_back( SMESH_TNodeXYZ( *n ));
     myOrderedNodes[ iSub++] = *n;
   }
 
@@ -2963,11 +2963,6 @@ bool SMESH_Pattern::Apply (SMESH_Mesh*                     theMesh,
   myXYZ.resize( myPoints.size() * theFaces.size(), undefinedXYZ() );
   myElements.reserve( theFaces.size() );
 
-  // to find point index
-  map< TPoint*, int > pointIndex;
-  for ( int i = 0; i < myPoints.size(); i++ )
-    pointIndex.insert( make_pair( & myPoints[ i ], i ));
-
   int ind1 = 0; // lowest point index for a face
 
   // meshed geometry
@@ -3019,7 +3014,7 @@ bool SMESH_Pattern::Apply (SMESH_Mesh*                     theMesh,
     {
       list< TPoint* > & linkPoints = getShapePoints( eID++ );
       const SMDS_MeshNode* n1 = myOrderedNodes[ i ];
-      const SMDS_MeshNode* n2 = myOrderedNodes[ i + 1 == nbNodes ? 0 : i + 1 ];
+      const SMDS_MeshNode* n2 = myOrderedNodes[( i+1 ) % nbNodes ];
       // make a link and a node set
       TNodeSet linkSet, node1Set;
       linkSet.insert( n1 );
@@ -3028,7 +3023,7 @@ bool SMESH_Pattern::Apply (SMESH_Mesh*                     theMesh,
       list< TPoint* >::iterator p = linkPoints.begin();
       {
         // map the first link point to n1
-        int nId = pointIndex[ *p ] + ind1;
+        int nId = ( *p - &myPoints[0] ) + ind1;
         myXYZIdToNodeMap[ nId ] = n1;
         list< list< int > >& groups = myIdsOnBoundary[ node1Set ];
         groups.push_back(list< int > ());
@@ -3040,7 +3035,7 @@ bool SMESH_Pattern::Apply (SMESH_Mesh*                     theMesh,
       list< int >& indList = groups.back();
       // add points to the map excluding the end points
       for ( p++; *p != linkPoints.back(); p++ )
-        indList.push_back( pointIndex[ *p ] + ind1 );
+        indList.push_back( ( *p - &myPoints[0] ) + ind1 );
     }
     ind1 += myPoints.size();
   }
@@ -3443,7 +3438,7 @@ void SMESH_Pattern::mergePoints (const bool uniteGroups)
       Bnd_Box box;
       TNodeSet::const_iterator n = nodes.begin();
       for ( ; n != nodes.end(); ++n )
-        box.Add( gp_Pnt( (*n)->X(), (*n)->Y(), (*n)->Z() ));
+        box.Add( gp_Pnt( SMESH_TNodeXYZ( *n )));
       double x, y, z, X, Y, Z;
       box.Get( x, y, z, X, Y, Z );
       gp_Pnt p( x, y, z ), P( X, Y, Z );
@@ -3454,7 +3449,7 @@ void SMESH_Pattern::mergePoints (const bool uniteGroups)
     bool unite = ( uniteGroups && nodes.size() == 2 );
     map< double, int > distIndMap;
     const SMDS_MeshNode* node = *nodes.begin();
-    gp_Pnt P( node->X(), node->Y(), node->Z() );
+    gp_Pnt P = SMESH_TNodeXYZ( node );
 
     // compare points, replace indices
 
@@ -3928,32 +3923,142 @@ bool SMESH_Pattern::MakeMesh(SMESH_Mesh* theMesh,
                                              myXYZ[ i ].Y(),
                                              myXYZ[ i ].Z());
     }
-  }
+    if ( theMesh->HasShapeToMesh() )
+    {
+      // set nodes on EDGEs (IMP 22368)
+      SMESH_MesherHelper helper( *theMesh );
+      helper.ToFixNodeParameters( true );
+      map< TNodeSet, list< list< int > > >::iterator idListIt = myIdsOnBoundary.begin();
+      for ( ; idListIt != myIdsOnBoundary.end(); idListIt++ )
+      {
+        list<list< int > >& groups = idListIt->second;
+        const TNodeSet&      nodes = idListIt->first;
+        if ( nodes.size() != 2 )
+          continue; // not a link
+        const SMDS_MeshNode* n1 = *nodes.begin();
+        const SMDS_MeshNode* n2 = *nodes.rbegin();
+        TopoDS_Shape S1 = helper.GetSubShapeByNode( n1, aMeshDS );
+        TopoDS_Shape S2 = helper.GetSubShapeByNode( n2, aMeshDS );
+        if ( S1.IsNull() || S1.ShapeType() < TopAbs_EDGE ||
+             S2.IsNull() || S2.ShapeType() < TopAbs_EDGE )
+          continue;
+        TopoDS_Shape S;
+        if ( S1.ShapeType() == TopAbs_EDGE )
+        {
+          if ( S1 == S2 || helper.IsSubShape( S2, S1 ))
+            S = S1;
+        }
+        else if ( S2.ShapeType() == TopAbs_EDGE )
+        {
+          if ( helper.IsSubShape( S1, S2 ))
+            S = S2;
+        }
+        else
+        {
+          S = helper.GetCommonAncestor( S1, S2, *theMesh, TopAbs_EDGE );
+        } 
+        if ( S.IsNull() )
+          continue;
+        const TopoDS_Edge & E = TopoDS::Edge( S );
+        helper.SetSubShape( E );
+        list<list< int > >::iterator g = groups.begin();
+        for ( ; g != groups.end(); ++g )
+        {
+          list< int >& ids = *g;
+          list< int >::iterator id = ids.begin();
+          for ( ; id != ids.end(); ++id )
+            if ( nodesVector[ *id ] && nodesVector[ *id ]->getshapeId() < 1 )
+            {
+              double u = 1e100;
+              aMeshDS->SetNodeOnEdge( nodesVector[ *id ], E, u );
+              helper.CheckNodeU( E, nodesVector[ *id ], u, 1e-7, true );
+            }
+        }
+      }
+    }
+  }  // if ( onMeshElements )
+
   else
   {
     nodesVector.resize( myPoints.size(), 0 );
 
-    // to find point index
-    map< TPoint*, int > pointIndex;
-    for ( int i = 0; i < myPoints.size(); i++ )
-      pointIndex.insert( make_pair( & myPoints[ i ], i ));
+    // find existing nodes on EDGEs and VERTEXes (IMP 22368)
+    map< int, list< TPoint* > >::iterator idPointIt = myShapeIDToPointsMap.begin();
+    if ( !myShapeIDMap.IsEmpty() && aMeshDS->NbNodes() > 0 )
+
+      for ( ; idPointIt != myShapeIDToPointsMap.end(); idPointIt++ )
+      {
+        const TopoDS_Shape&    S = myShapeIDMap( idPointIt->first );
+        list< TPoint* > & points = idPointIt->second;
+        if ( points.empty() )
+          continue;
+
+        switch ( S.ShapeType() )
+        {
+        case TopAbs_VERTEX:
+        {
+          int pIndex = points.back() - &myPoints[0];
+          if ( !nodesVector[ pIndex ] )
+            nodesVector[ pIndex ] = SMESH_Algo::VertexNode( TopoDS::Vertex( S ), aMeshDS );
+          break;
+        }
+        case TopAbs_EDGE:
+        {
+          const TopoDS_Edge& edge = TopoDS::Edge( S );
+          map< double, const SMDS_MeshNode* > paramsOfNodes;
+          if ( !SMESH_Algo::GetSortedNodesOnEdge( aMeshDS, edge,
+                                                  /*ignoreMediumNodes=*/false,
+                                                  paramsOfNodes )
+               || paramsOfNodes.size() < 3 )
+            break;
+          // points on VERTEXes are included with wrong myU
+          list< TPoint* >::reverse_iterator pItR = ++points.rbegin();
+          list< TPoint* >::iterator         pItF = ++points.begin();
+          const bool isForward = ( (*pItF)->myU < (*pItR)->myU );
+          map< double, const SMDS_MeshNode* >::iterator u2n    = ++paramsOfNodes.begin();
+          map< double, const SMDS_MeshNode* >::iterator u2nEnd = --paramsOfNodes.end();
+          TPoint* p;
+          while ( u2n != u2nEnd && pItF != points.end() )
+          {
+            const double         u = u2n->first;
+            const SMDS_MeshNode* n = u2n->second;
+            const double       tol = ( (++u2n)->first - u ) / 20;
+            do
+            {
+              p = ( isForward ? *pItF : *pItR );
+              if ( Abs( u - p->myU ) < tol )
+              {
+                int pIndex = p - &myPoints[0];
+                if ( !nodesVector [ pIndex ] )
+                  nodesVector [ pIndex ] = n;
+                ++pItF;
+                ++pItR;
+                break;
+              }
+            }
+            while ( p->myU < u && ( ++pItF, ++pItR != points.rend() ));
+          }
+          break;
+        }
+        default:;
+        }
+      } // end of "find existing nodes on EDGEs and VERTEXes"
 
     // loop on sub-shapes of myShape: create nodes
-    map< int, list< TPoint* > >::iterator idPointIt = myShapeIDToPointsMap.begin();
+    idPointIt = myShapeIDToPointsMap.begin();
     for ( ; idPointIt != myShapeIDToPointsMap.end(); idPointIt++ )
     {
       TopoDS_Shape S;
-      //SMESHDS_SubMesh * subMeshDS = 0;
       if ( !myShapeIDMap.IsEmpty() ) {
         S = myShapeIDMap( idPointIt->first );
-        //subMeshDS = aMeshDS->MeshElements( S );
       }
       list< TPoint* > & points = idPointIt->second;
       list< TPoint* >::iterator pIt = points.begin();
       for ( ; pIt != points.end(); pIt++ )
       {
         TPoint* point = *pIt;
-        int pIndex = pointIndex[ point ];
+        //int pIndex = pointIndex[ point ];
+        int pIndex = point - &myPoints[0];
         if ( nodesVector [ pIndex ] )
           continue;
         SMDS_MeshNode* node = aMeshDS->AddNode (point->myXYZ.X(),
@@ -4148,8 +4253,9 @@ void SMESH_Pattern::createElements(SMESH_Mesh*                            theMes
         SMDS_ElemIteratorPtr noIt = elem->nodesIterator();
         while ( noIt->more() ) {
           SMDS_MeshNode* node = const_cast<SMDS_MeshNode*>(smdsNode( noIt->next() ));
-          if (!node->getshapeId() &&
-              shellNodes.find( node ) == shellNodes.end() ) {
+          if ( node->getshapeId() < 1 &&
+               shellNodes.find( node ) == shellNodes.end() )
+          {
             if ( S.ShapeType() == TopAbs_FACE )
               aMeshDS->SetNodeOnFace( node, shapeID,
                                       Precision::Infinite(),// <- it's a sign that UV is not set
index b2606e4..b57624c 100644 (file)
@@ -356,7 +356,7 @@ private:
   // all functions assure that shapes are indexed so that first go
   // ordered vertices, then ordered edge, then faces and maybe a shell
   TopTools_IndexedMapOfOrientedShape   myShapeIDMap;
-  std::map< int, std::list< TPoint* > >     myShapeIDToPointsMap;
+  std::map< int, std::list< TPoint*> > myShapeIDToPointsMap;
 
   // for the 2d case:
   // nb of key-points in each of pattern boundaries
index fa9add8..e8ebeea 100644 (file)
@@ -551,6 +551,8 @@ QString SMESHGUI_GenericHypothesisCreator::helpPage() const
     aHelpFileName = "a1d_meshing_hypo_page.html#max_length_anchor";
   else if ( aHypType == "Arithmetic1D")
     aHelpFileName = "a1d_meshing_hypo_page.html#arithmetic_1d_anchor";
+  else if ( aHypType == "GeometricProgression")
+    aHelpFileName = "a1d_meshing_hypo_page.html#geometric_1d_anchor";
   else if ( aHypType == "FixedPoints1D")
     aHelpFileName = "a1d_meshing_hypo_page.html#fixed_points_1d_anchor";
   else if ( aHypType == "MaxElementArea")
index 2d45966..5800395 100644 (file)
@@ -296,24 +296,45 @@ namespace SMESH
   }
 
 
-  QStringList GetHypothesesSets(int maxDim)
+  QStringList GetHypothesesSets(int maxDim, const QString& MeshType)
   {
     QStringList aSetNameList;
 
     // Init list of available hypotheses, if needed
     InitAvailableHypotheses();
-
     QList<HypothesesSet*>::iterator hypoSet;
-    for ( hypoSet  = myListOfHypothesesSets.begin(); 
+    for ( hypoSet  = myListOfHypothesesSets.begin();
           hypoSet != myListOfHypothesesSets.end();
           ++hypoSet ) {
       HypothesesSet* aSet = *hypoSet;
-      if ( aSet &&
-           ( aSet->count( true ) || aSet->count( false )) &&
-           aSet->maxDim() <= maxDim)
+      bool isAvailable = false;
+      if ( !MeshType.isEmpty() )
+      {
+        if ( aSet->maxDim() != maxDim)
+          continue;
+        aSet->init( true );
+        while ( aSet->next(), aSet->more() )
+        {
+          if ( HypothesisData* hypData = SMESH::GetHypothesisData( aSet->current() ) )
+          {
+            QStringList::const_iterator inElemType = hypData->OutputTypes.begin();
+            for ( ; inElemType != hypData->OutputTypes.end(); inElemType++ )
+            {
+              if ( *inElemType == MeshType ){
+                isAvailable = true;
+                break;
+              }
+            }
+          }
+          if ( isAvailable ) break;
+        }
+      }
+      else if ( aSet && ( aSet->count( true ) || aSet->count( false )) &&
+                aSet->maxDim() <= maxDim)
       {
-        aSetNameList.append( mangledHypoSetName( aSet ));
+        isAvailable = true;
       }
+      if ( isAvailable ) aSetNameList.append( mangledHypoSetName( aSet ));
     }
     aSetNameList.removeDuplicates();
     aSetNameList.sort();
index a2c12f5..333418b 100644 (file)
@@ -71,7 +71,7 @@ namespace SMESH
                                       const bool = false,
                                       const bool = true);
   SMESHGUI_EXPORT
-  QStringList GetHypothesesSets( int maxDim );
+  QStringList GetHypothesesSets( int, const QString& );
 
   SMESHGUI_EXPORT
   HypothesesSet* GetHypothesesSet( const QString& );
index 4b030bb..94d9ce6 100644 (file)
@@ -364,6 +364,9 @@ SMESHGUI_MeshDlg::SMESHGUI_MeshDlg( const bool theToCreate, const bool theIsMesh
   // geometry
   createObject( tr( "GEOMETRY" ), mainFrame(), Geom );
   myGeomPopup = 0;
+  // mesh type
+  QLabel* anMeshTypeLbl = new QLabel( tr( "MESH_TYPE" ), this );
+  myMeshType = new QComboBox( this );
   
   // Create tab widget
   
@@ -398,10 +401,16 @@ SMESHGUI_MeshDlg::SMESHGUI_MeshDlg( const bool theToCreate, const bool theIsMesh
   aLay->addWidget( objectWg( Geom, Label ),   2, 0 );
   aLay->addWidget( objectWg( Geom, Btn ),     2, 1 );
   aLay->addWidget( objectWg( Geom, Control ), 2, 2 );
-  aLay->addWidget( myTabWg,                   4, 0, 1, 3 );
-  aLay->addWidget( myHypoSetButton,           5, 0, 1, 3 );
+  aLay->addWidget( anMeshTypeLbl,             3, 0 );
+  aLay->addWidget( myMeshType,                3, 2 );
+  aLay->addWidget( myTabWg,                   5, 0, 1, 3 );
+  aLay->addWidget( myHypoSetButton,           6, 0, 1, 3 );
   aLay->setRowMinimumHeight( 3, 20 );
 
+  myMeshType->clear();
+
+  // Connect signals and slots
+  connect( myMeshType, SIGNAL( activated( int ) ), SLOT( onChangedMeshType( int ) ) );
   // Disable controls if necessary
   setObjectShown( Mesh, false );
   if ( theToCreate )
@@ -615,3 +624,36 @@ int SMESHGUI_MeshDlg::getActiveObject()
       return i;
   return -1;
 }
+//================================================================================
+/*!
+ * \brief Sets available types of mesh
+ * \param theTypeMesh - list of available types of mesh
+ */
+//================================================================================
+void SMESHGUI_MeshDlg::setAvailableMeshType( const QStringList& theTypeMesh )
+{
+  myMeshType->clear();
+  myMeshType->addItems(theTypeMesh);
+}
+//================================================================================
+/*!
+ * \brief Emits selectMeshType( const int, const int ) signal
+ *
+ * SLOT is called when a combo box "mesh type" is selected.
+ */
+//================================================================================
+void SMESHGUI_MeshDlg::onChangedMeshType( const int isIndex )
+{
+  emit selectMeshType( Dim3D - myTabWg->currentIndex(), isIndex );
+}
+//================================================================================
+/*!
+ * \brief Get current index types of mesh
+ */
+//================================================================================
+int SMESHGUI_MeshDlg::currentMeshType( )
+{
+  return myMeshType->currentIndex( );
+}
+
+
index 6ba8e6c..08f11ff 100644 (file)
@@ -74,21 +74,26 @@ public:
   void                         enableTab(const int);
   bool                         isTabEnabled(const int) const;
   int                          getActiveObject();
+  void                         setAvailableMeshType(const QStringList& );
+  int                          currentMeshType();
 
 signals:
   void                         hypoSet( const QString& );
   void                         geomSelectionByMesh( bool );
+  void                         selectMeshType( const int, const int );
 
 private slots:  
   void                         onHypoSetPopup( QAction* );
   void                         onGeomPopup( QAction* );
   void                         onGeomSelectionButton( bool );
+  void                         onChangedMeshType( const int );
 
 private:
   QMap<int, SMESHGUI_MeshTab*> myTabs;
   QTabWidget*                  myTabWg;
   QToolButton*                 myHypoSetButton;
   QMenu*                       myGeomPopup;
+  QComboBox*                   myMeshType;
 };
 
 /*!
index 62a5c4b..eb0a30e 100644 (file)
@@ -95,6 +95,7 @@ SMESHGUI_MeshOp::SMESHGUI_MeshOp( const bool theToCreate, const bool theIsMesh )
   if ( GeometryGUI::GetGeomGen()->_is_nil() )// check that GEOM_Gen exists
     GeometryGUI::InitGeomGen();
   myIsOnGeometry = true;
+  myMaxShapeDim = -1;
 }
 
 //================================================================================
@@ -211,14 +212,13 @@ void SMESHGUI_MeshOp::startOperation()
     }
     connect( myDlg, SIGNAL( hypoSet( const QString& )), SLOT( onHypoSet( const QString& )));
     connect( myDlg, SIGNAL( geomSelectionByMesh( bool )), SLOT( onGeomSelectionByMesh( bool )));
-
+    connect( myDlg, SIGNAL( selectMeshType( const int, const int ) ), SLOT( onAlgoSetByMeshType( const int, const int)));
     if ( myToCreate )
       if ( myIsMesh ) myHelpFileName = "constructing_meshes_page.html";
       else myHelpFileName = "constructing_submeshes_page.html";
     else myHelpFileName = "editing_meshes_page.html";
   }
   SMESHGUI_SelectionOp::startOperation();
-
   // iterate through dimensions and get available algoritms, set them to the dialog
   _PTR(SComponent) aFather = SMESH::GetActiveStudyDocument()->FindComponent( "SMESH" );
   for ( int i = SMESH::DIM_0D; i <= SMESH::DIM_3D; i++ )
@@ -245,6 +245,11 @@ void SMESHGUI_MeshOp::startOperation()
     myDlg->activateObject( SMESHGUI_MeshDlg::Obj );
 
   myDlg->setCurrentTab( SMESH::DIM_3D );
+
+  QStringList TypeMeshList;
+  createMeshTypeList( TypeMeshList );
+  setAvailableMeshType( TypeMeshList );
+
   myDlg->show();
   myDlg->setGeomPopupEnabled(false);
   selectionDone();
@@ -582,7 +587,8 @@ void SMESHGUI_MeshOp::selectionDone()
         onAlgoSelected(-1, i);
       }
       myDlg->setMaxHypoDim( shapeDim );
-      myDlg->setHypoSets( SMESH::GetHypothesesSets( shapeDim ));
+      myMaxShapeDim = shapeDim;
+      myDlg->setHypoSets( SMESH::GetHypothesesSets( shapeDim, "" ));
 
       if (!myToCreate) // edition: read hypotheses
       {
@@ -669,12 +675,16 @@ void SMESHGUI_MeshOp::selectionDone()
       for (int i = SMESH::DIM_0D;i < SMESH::DIM_3D; ++i) {
         myDlg->disableTab(i);
       }
+      myMaxShapeDim = -1;
       //Hide labels and fields (Mesh ang Geometry)
       myDlg->setObjectShown( SMESHGUI_MeshDlg::Mesh, false );
       myDlg->setObjectShown( SMESHGUI_MeshDlg::Geom, false );
       myDlg->adjustSize();
       readMesh();
     }
+    QStringList TypeMeshList;
+    createMeshTypeList( TypeMeshList );
+    setAvailableMeshType( TypeMeshList );
   }
   catch ( const SALOME::SALOME_Exception& S_ex )
   {
@@ -1380,7 +1390,19 @@ void SMESHGUI_MeshOp::onAlgoSelected( const int theIndex,
     availableHyps( aDim, Algo, anAvailable, myAvailableHypData[ aDim ][ Algo ]);
     myDlg->tab( aDim )->setAvailableHyps( Algo, anAvailable );
   }
-
+  // check that tab enable, if algorithm building needed algo is one less than dimension
+  if ( algoData && myIsOnGeometry && !algoData->InputTypes.isEmpty() &&
+     ( aDim > SMESH::DIM_0D ) && !isAccessibleDim( aDim - 1 ) ){
+    myDlg->enableTab( aDim - 1 );
+  }
+  if ( (myDlg->currentMeshType() != MT_ANY) &&
+          (( !algoData && ( aDim > SMESH::DIM_0D ) && isAccessibleDim( aDim - 1 )) ||
+       ( algoData && myIsOnGeometry && algoData->InputTypes.isEmpty() &&
+       ( aDim > SMESH::DIM_0D ) && isAccessibleDim( aDim - 1 ) ) ) ){
+    for (int i = aDim - 1; i >= SMESH::DIM_0D; i--){
+      if ( isAccessibleDim( i ) ) myDlg->disableTab( i );
+    }
+  }
   // check that algorithms of other dimentions are compatible with
   // the selected one
 
@@ -2343,3 +2365,193 @@ void SMESHGUI_MeshOp::selectObject( _PTR(SObject) theSObj ) const
     sm->setSelectedObjects( anIOList, false );
   }
 }
+//================================================================================
+/*!
+ * \brief Create available list types of mesh
+  * \param theTypeMesh - Output list of available types of mesh
+ */
+//================================================================================
+void SMESHGUI_MeshOp::createMeshTypeList( QStringList& theTypeMesh)
+{
+  theTypeMesh.clear();
+  theTypeMesh.append( tr( "MT_ANY" ) );
+  if ( myIsOnGeometry &&  ( myMaxShapeDim >= 2 || myMaxShapeDim == -1 ) )
+  {
+    theTypeMesh.append( tr( "MT_TRIANGULAR" ) );
+    theTypeMesh.append( tr( "MT_QUADRILATERAL" ) );
+  }
+  if ( myIsOnGeometry && ( myMaxShapeDim == 3 || myMaxShapeDim == -1 ) )
+  {
+    theTypeMesh.append( tr( "MT_TETRAHEDRAL" ) );
+    theTypeMesh.append( tr( "MT_HEXAHEDRAL" ) );
+  }
+
+}
+//================================================================================
+/*!
+ * \brief Set available types of mesh
+  * \param theTypeMesh - List of available types of mesh
+ */
+//================================================================================
+void SMESHGUI_MeshOp::setAvailableMeshType( const QStringList& theTypeMesh )
+{
+  myDlg->setAvailableMeshType( theTypeMesh );
+}
+
+//================================================================================
+/*!
+ * \brief SLOT. Is called when the user select type of mesh
+  * \param theTabIndex - Index of current active tab
+  * \param theIndex - Index of current type of mesh
+ */
+//================================================================================
+void SMESHGUI_MeshOp::onAlgoSetByMeshType( const int theTabIndex, const int theIndex)
+{
+  int aDim;
+  if ( !myIsOnGeometry ) return;
+  THypDataList anAvailableAlgsData;
+  QStringList anAvailableAlgs;
+  QString anCompareType = "ANY";
+  bool isAvailableChoiceAlgo = false;
+  int anCurrentAvailableAlgo = 0;
+  bool isNone = true;
+  switch ( theIndex ) {
+    case MT_ANY:
+    {
+      for ( int dim = SMESH::DIM_2D; dim <= SMESH::DIM_3D; dim++ )
+      {
+        isNone = currentHyp( dim, Algo ) < 0;
+        isAvailableChoiceAlgo = false;
+        // retrieves a list of available algorithms from resources
+        availableHyps( dim, Algo, anAvailableAlgs, anAvailableAlgsData );
+        //return current algo in current tab
+        if ( !isNone && !myAvailableHypData[dim][Algo].empty() ){
+          for (int i = 0 ; i < anAvailableAlgsData.count(); i++)
+          {
+            HypothesisData* algoAny = anAvailableAlgsData.at(i);
+            HypothesisData* algoCur = myAvailableHypData[dim][Algo].at( currentHyp( dim, Algo ) );
+            QString tem = algoAny->Label;
+            if ( algoAny->Label == algoCur->Label ){
+              isAvailableChoiceAlgo = true;
+              anCurrentAvailableAlgo = i;
+              break;
+            }
+          }
+        }
+        else if ( !isNone ){
+          isAvailableChoiceAlgo = true;
+          anCurrentAvailableAlgo = currentHyp( dim, Algo );
+        }
+        //set new algorithm list and select the current algorithm
+        myAvailableHypData[dim][Algo] = anAvailableAlgsData;
+        myDlg->tab( dim )->setAvailableHyps( Algo, anAvailableAlgs );
+        if ( isAvailableChoiceAlgo )
+          setCurrentHyp( dim, Algo, anCurrentAvailableAlgo );
+      }
+      int aMaxShapeDim = ( myMaxShapeDim == -1 ) ? SMESH::DIM_3D : myMaxShapeDim;
+      for ( int i = SMESH::DIM_0D; i <= aMaxShapeDim; i++ ) {
+        myDlg->enableTab( i );
+      }
+      myDlg->setCurrentTab( theTabIndex );
+      myDlg->setHypoSets( SMESH::GetHypothesesSets( aMaxShapeDim, "" ) );
+    }
+    break;
+    case MT_TRIANGULAR:{
+      aDim = SMESH::DIM_2D;
+      anCompareType = "TRIA";
+    }
+    break;
+    case MT_QUADRILATERAL:{
+      aDim = SMESH::DIM_2D;
+      anCompareType = "QUAD";
+    }
+    break;
+    case MT_TETRAHEDRAL:{
+      aDim = SMESH::DIM_3D;
+      anCompareType = "TETRA";
+    }
+    break;
+    case MT_HEXAHEDRAL:{
+      aDim = SMESH::DIM_3D;
+      anCompareType = "HEXA";
+    }
+    break;
+    default:;
+  }
+  if ( anCompareType != "ANY" )
+  {
+    QString anCurrentAlgo;
+    bool isReqDisBound = true;
+    isNone = currentHyp( aDim, Algo ) < 0;
+    // retrieves a list of available algorithms from resources
+    availableHyps( aDim, Algo, anAvailableAlgs, anAvailableAlgsData );
+    // finding algorithm which is selected
+    if ( !isNone && !myAvailableHypData[aDim][Algo].empty() &&
+         myAvailableHypData[aDim][Algo].count() != anAvailableAlgsData.count() ){
+      anCurrentAlgo = myAvailableHypData[aDim][Algo].at( currentHyp( aDim, Algo ) )->Label;
+      isReqDisBound = myAvailableHypData[aDim][Algo].at( currentHyp( aDim, Algo ) )->InputTypes.isEmpty();
+    }
+    else if ( !isNone ){
+      anCurrentAlgo = anAvailableAlgsData.at( currentHyp( aDim, Algo ) )->Label;
+      isReqDisBound = anAvailableAlgsData.at( currentHyp( aDim, Algo ) )->InputTypes.isEmpty();
+    }
+    anAvailableAlgs.clear();
+    myAvailableHypData[aDim][Algo].clear();
+    // finding and adding algorithm depending on the type mesh
+    for ( int i = 0 ; i < anAvailableAlgsData.count(); i++ )
+    {
+      HypothesisData* algoIn = anAvailableAlgsData.at( i );
+      bool isAvailableAlgo = ( algoIn->OutputTypes.count() == 0 );
+      QStringList::const_iterator inElemType = algoIn->OutputTypes.begin();
+      for ( ; inElemType != algoIn->OutputTypes.end(); inElemType++ )
+      {
+        if ( *inElemType == anCompareType ){
+          isAvailableAlgo = true;
+          break;
+        }
+      }
+      if ( isAvailableAlgo || algoIn->OutputTypes.count()==0 ){
+        anAvailableAlgs.append( algoIn->Label );
+        myAvailableHypData[aDim][Algo].append( algoIn );
+      }
+      //algorithm will be active, if the chosen algorithm available in the current mesh type
+      if ( !isNone &&  isAvailableAlgo && algoIn->Label == anCurrentAlgo ){
+        isAvailableChoiceAlgo = true;
+        anCurrentAvailableAlgo = anAvailableAlgs.count() - 1 ;
+      }
+    }
+    //set new algorithm list and select the current algorithm
+    myDlg->tab( aDim )->setAvailableHyps( Algo, anAvailableAlgs );
+    if ( isAvailableChoiceAlgo )
+      setCurrentHyp( aDim, Algo, anCurrentAvailableAlgo );
+    int aMaxShapeDim = ( myMaxShapeDim == -1 ) ? SMESH::DIM_3D : myMaxShapeDim;
+    if ( isNone || isReqDisBound || !isAvailableChoiceAlgo ) {
+      for ( int i = SMESH::DIM_0D; i <= aMaxShapeDim; i++ ) {
+        if ( aDim != i ) {
+          myDlg->disableTab( i );
+          setCurrentHyp(i, Algo, -1);
+        }
+      }
+    }
+    else if ( !isNone ){
+      if ( aDim == SMESH::DIM_2D){
+        myDlg->disableTab( SMESH::DIM_3D );
+        setCurrentHyp( SMESH::DIM_3D, Algo, -1);
+      }
+      for ( int i = aMaxShapeDim; i > SMESH::DIM_0D; i-- )
+      {
+        isReqDisBound = ( currentHyp( i, Algo ) < 0 ) ? true : myAvailableHypData[i][Algo].at( currentHyp( i, Algo ) )->InputTypes.isEmpty();
+        if ( aMaxShapeDim != i && isReqDisBound) {
+          for (int j = i - 1; j >= SMESH::DIM_0D; j--){
+            myDlg->disableTab( j );
+            setCurrentHyp( j , Algo, -1 );
+          }
+          break;
+        }
+      }
+    }
+    myDlg->setHypoSets( SMESH::GetHypothesesSets( aDim, anCompareType ) );
+    myDlg->enableTab( aDim );
+    myDlg->setCurrentTab( aDim );
+  }
+}
index 91c39ed..5f56882 100644 (file)
@@ -48,6 +48,7 @@ class SMESHGUI_EXPORT SMESHGUI_MeshOp : public SMESHGUI_SelectionOp
       
 public:
   enum HypType{ Algo = 0, MainHyp, AddHyp, NbHypTypes };
+  enum MeshType{ MT_ANY = 0, MT_TRIANGULAR, MT_QUADRILATERAL, MT_TETRAHEDRAL, MT_HEXAHEDRAL };
 
   typedef std::pair<SMESH::SMESH_Hypothesis_var, QString> THypItem;
   typedef QList< THypItem > THypList;
@@ -83,6 +84,7 @@ protected slots:
   void                           processSet();
   void                           onHypoCreated( int );
   void                           onHypoEdited( int );
+  void                           onAlgoSetByMeshType( const int, const int );
 
 private:
   typedef QList<HypothesisData*> THypDataList; // typedef: list of hypothesis data
@@ -125,7 +127,8 @@ private:
   char*                          isSubmeshIgnored() const;
   _PTR(SObject)                  getSubmeshByGeom() const;
   void                           selectObject( _PTR(SObject) ) const;
-
+  void                           createMeshTypeList( QStringList& );
+  void                           setAvailableMeshType( const QStringList& );
 private:
   SMESHGUI_MeshDlg*              myDlg;
   SMESHGUI_ShapeByMeshOp*        myShapeByMeshOp;
@@ -136,13 +139,12 @@ private:
   TDim2Type2HypList              myExistingHyps; //!< all hypothesis of SMESH module
   TDim2Type2HypList              myObjHyps;      //!< hypothesis assigned to the current 
                                                  //   edited mesh/sub-mesh
-
   // hypdata corresponding to hypotheses present in myDlg
   THypDataList                   myAvailableHypData[4][NbHypTypes];
 
   bool                           myIgnoreAlgoSelection;
   HypothesesSet* myHypoSet;
-  int myDim, myType;
+  int myDim, myType, myMaxShapeDim;
 
   QString                        myObjectToSelect;
 };
index 174bcf1..df41ace 100644 (file)
@@ -5978,6 +5978,10 @@ Please specify them and try again</translation>
         <translation>Mesh</translation>
     </message>
     <message>
+        <source>MESH_TYPE</source>
+        <translation>Mesh type</translation>
+    </message>
+    <message>
         <source>NAME</source>
         <translation>Name</translation>
     </message>
@@ -6033,6 +6037,26 @@ Please specify it and try again</translation>
         <translation>Mesh is null</translation>
     </message>
     <message>
+        <source>MT_ANY</source>
+        <translation>Any</translation>
+    </message>
+    <message>
+        <source>MT_HEXAHEDRAL</source>
+        <translation>Hexahedral</translation>
+    </message>
+    <message>
+        <source>MT_TETRAHEDRAL</source>
+        <translation>Tetrahedral</translation>
+    </message>
+    <message>
+        <source>MT_TRIANGULAR</source>
+        <translation>Triangular</translation>
+    </message>
+    <message>
+        <source>MT_QUADRILATERAL</source>
+        <translation>Quadrilaterial</translation>
+    </message>
+    <message>
         <source>NAME_OF_MESH_IS_EMPTY</source>
         <translation>Name of mesh is empty
 Please enter valid name and try again</translation>
index 749b17e..b8ddd3e 100755 (executable)
@@ -5972,6 +5972,10 @@ Indiquez-les et essayez de nouveau</translation>
         <translation>Maillage</translation>
     </message>
     <message>
+        <source>MESH_TYPE</source>
+        <translation type="unfinished">Mesh type</translation>
+    </message>
+    <message>
         <source>NAME</source>
         <translation>Nom</translation>
     </message>
@@ -6027,6 +6031,26 @@ Spécifiez-le et essayez de nouveau</translation>
         <translation>Le maillage est nul</translation>
     </message>
     <message>
+        <source>MT_ANY</source>
+        <translation type="unfinished">Any</translation>
+    </message>
+    <message>
+        <source>MT_HEXAHEDRAL</source>
+        <translation type="unfinished">Hexahedral</translation>
+    </message>
+    <message>
+        <source>MT_TETRAHEDRAL</source>
+        <translation type="unfinished">Tetrahedral</translation>
+    </message>
+    <message>
+        <source>MT_TRIANGULAR</source>
+        <translation type="unfinished">Triangular</translation>
+    </message>
+    <message>
+        <source>MT_QUADRILATERAL</source>
+        <translation>Quadrilaterial</translation>
+    </message>
+    <message>
         <source>NAME_OF_MESH_IS_EMPTY</source>
         <translation>Le nom du maillage est vide
 Indiquez un nom valide et essayez de nouveau</translation>
index 7d2ed2a..fd40175 100644 (file)
 //
 #include "SMESH_Block.hxx"
 
+#include "SMDS_MeshNode.hxx"
+#include "SMDS_MeshVolume.hxx"
+#include "SMDS_VolumeTool.hxx"
+#include "SMESH_MeshAlgos.hxx"
+
 #include <BRepAdaptor_Curve.hxx>
 #include <BRepAdaptor_Curve2d.hxx>
 #include <BRepAdaptor_Surface.hxx>
 #include <math_Matrix.hxx>
 #include <math_Vector.hxx>
 
-#include "SMDS_MeshNode.hxx"
-#include "SMDS_MeshVolume.hxx"
-#include "SMDS_VolumeTool.hxx"
-#include "utilities.h"
+#include <utilities.h>
 
 #include <list>
 #include <limits>
@@ -309,24 +311,15 @@ gp_XYZ SMESH_Block::TFace::Point( const gp_XYZ& theParams ) const
 
 namespace
 {
+  inline
   bool isPntInTria( const gp_XY& p, const gp_XY& t0, const gp_XY& t1, const gp_XY& t2  )
   {
-    const double // matrix 2x2
-      T11 = t0.X()-t2.X(), T12 = t1.X()-t2.X(),
-      T21 = t0.Y()-t2.Y(), T22 = t1.Y()-t2.Y();
-    const double Tdet = T11*T22 - T12*T21; // matrix determinant
-    if ( Abs( Tdet ) < std::numeric_limits<double>::min() )
-      return false;
-    // matrix inverse
-    const double t11 = T22, t12 = -T12, t21 = -T21, t22 = T11;
-    // vector
-    const double r11 = p.X()-t2.X(), r12 = p.Y()-t2.Y();
-    // barycentric coordinates: mutiply matrix by vector
-    const double bc0 = (t11 * r11 + t12 * r12)/Tdet;
-    const double bc1 = (t21 * r11 + t22 * r12)/Tdet;
+    double bc0, bc1;
+    SMESH_MeshAlgos::GetBarycentricCoords( p, t0, t1, t2, bc0, bc1 );
     return ( bc0 >= 0. && bc1 >= 0. && bc0 + bc1 <= 1. );
   }
 
+  inline
   bool isPntInQuad( const gp_XY& p,
                     const gp_XY& q0, const gp_XY& q1, const gp_XY& q2, const gp_XY& q3 )
   {
index a2fcd0a..4fc6901 100644 (file)
@@ -68,20 +68,19 @@ class SMESHUtils_EXPORT SMESH_Block: public math_FunctionSetWithDerivatives
     // ----------------------------
     ID_NONE = 0,
 
-    ID_V000 = 1, ID_V100, ID_V010, ID_V110, ID_V001, ID_V101, ID_V011, ID_V111,
+    ID_V000 = 1, ID_V100, ID_V010, ID_V110, ID_V001, ID_V101, ID_V011, ID_V111, // 1-8
 
-    ID_Ex00, ID_Ex10, ID_Ex01, ID_Ex11,
-    ID_E0y0, ID_E1y0, ID_E0y1, ID_E1y1,
-    ID_E00z, ID_E10z, ID_E01z, ID_E11z,
+    ID_Ex00, ID_Ex10, ID_Ex01, ID_Ex11, // 9-12
+    ID_E0y0, ID_E1y0, ID_E0y1, ID_E1y1, // 13-16
+    ID_E00z, ID_E10z, ID_E01z, ID_E11z, // 17-20
 
-    ID_Fxy0, ID_Fxy1, ID_Fx0z, ID_Fx1z, ID_F0yz, ID_F1yz,
+    ID_Fxy0, ID_Fxy1, ID_Fx0z, ID_Fx1z, ID_F0yz, ID_F1yz, // 21-26
 
-    ID_Shell
-    };
-  enum { // to use TShapeID for indexing certain type subshapes
+    ID_Shell, // 27
 
-    ID_FirstV = ID_V000, ID_FirstE = ID_Ex00, ID_FirstF = ID_Fxy0
+    // to use TShapeID for indexing certain type subshapes
 
+    ID_FirstV = ID_V000, ID_FirstE = ID_Ex00, ID_FirstF = ID_Fxy0
   };
 
 
@@ -285,7 +284,7 @@ public:
                               std::list< int >  &       theNbEdgesInWires,
                               TopoDS_Vertex             theFirstVertex=TopoDS_Vertex(),
                               const bool                theShapeAnalysisAlgo=false);
-  // Return nb wires and a list of oredered edges.
+  // Return nb wires and a list of ordered edges.
   // It is used to assign indices to subshapes.
   // theFirstVertex may be NULL.
   // Always try to set a seam edge first
index 45acf33..2ecd29a 100644 (file)
@@ -1386,6 +1386,39 @@ double SMESH_MeshAlgos::GetDistance( const SMDS_MeshFace* face,
   return badDistance;
 }
 
+//================================================================================
+/*!
+ * \brief Returns barycentric coordinates of a point within a triangle.
+ *        A not returned bc2 = 1. - bc0 - bc1.
+ *        The point lies within the triangle if ( bc0 >= 0 && bc1 >= 0 && bc0+bc1 <= 1 )
+ */
+//================================================================================
+
+void SMESH_MeshAlgos::GetBarycentricCoords( const gp_XY& p,
+                                            const gp_XY& t0,
+                                            const gp_XY& t1,
+                                            const gp_XY& t2,
+                                            double &     bc0,
+                                            double &     bc1)
+{
+  const double // matrix 2x2
+    T11 = t0.X()-t2.X(), T12 = t1.X()-t2.X(),
+    T21 = t0.Y()-t2.Y(), T22 = t1.Y()-t2.Y();
+  const double Tdet = T11*T22 - T12*T21; // matrix determinant
+  if ( Abs( Tdet ) < std::numeric_limits<double>::min() )
+  {
+    bc0 = bc1 = 2.;
+    return;
+  }
+  // matrix inverse
+  const double t11 = T22, t12 = -T12, t21 = -T21, t22 = T11;
+  // vector
+  const double r11 = p.X()-t2.X(), r12 = p.Y()-t2.Y();
+  // barycentric coordinates: mutiply matrix by vector
+  bc0 = (t11 * r11 + t12 * r12)/Tdet;
+  bc1 = (t21 * r11 + t22 * r12)/Tdet;
+}
+
 //=======================================================================
 //function : FindFaceInSet
 //purpose  : Return a face having linked nodes n1 and n2 and which is
index a02ba10..dc4d8fc 100644 (file)
@@ -97,9 +97,16 @@ namespace SMESH_MeshAlgos
   /*!
    * \brief Return true if the point is IN or ON of the element
    */
-  SMESHUtils_EXPORT bool IsOut( const SMDS_MeshElement* element, const gp_Pnt& point, double tol );
+  SMESHUtils_EXPORT
+  bool IsOut( const SMDS_MeshElement* element, const gp_Pnt& point, double tol );
 
-  SMESHUtils_EXPORT double GetDistance( const SMDS_MeshFace* face, const gp_Pnt& point );
+  SMESHUtils_EXPORT
+  double GetDistance( const SMDS_MeshFace* face, const gp_Pnt& point );
+
+  SMESHUtils_EXPORT
+  void GetBarycentricCoords( const gp_XY& point,
+                             const gp_XY& t0, const gp_XY& t1, const gp_XY& t2,
+                             double &    bc0, double &    bc1);
 
   /*!
    * Return a face having linked nodes n1 and n2 and which is
@@ -107,35 +114,40 @@ namespace SMESH_MeshAlgos
    * - in elemSet provided that !elemSet.empty()
    * i1 and i2 optionally returns indices of n1 and n2
    */
-  SMESHUtils_EXPORT const SMDS_MeshElement* 
-    FindFaceInSet(const SMDS_MeshNode*    n1,
-                  const SMDS_MeshNode*    n2,
-                  const TIDSortedElemSet& elemSet,
-                  const TIDSortedElemSet& avoidSet,
-                  int*                    i1=0,
-                  int*                    i2=0);
+  SMESHUtils_EXPORT
+  const SMDS_MeshElement* FindFaceInSet(const SMDS_MeshNode*    n1,
+                                        const SMDS_MeshNode*    n2,
+                                        const TIDSortedElemSet& elemSet,
+                                        const TIDSortedElemSet& avoidSet,
+                                        int*                    i1=0,
+                                        int*                    i2=0);
   /*!
    * \brief Calculate normal of a mesh face
    */
-  SMESHUtils_EXPORT bool FaceNormal(const SMDS_MeshElement* F, gp_XYZ& normal, bool normalized=true);
+  SMESHUtils_EXPORT
+  bool FaceNormal(const SMDS_MeshElement* F, gp_XYZ& normal, bool normalized=true);
 
   /*!
    * \brief Return nodes common to two elements
    */
-  SMESHUtils_EXPORT std::vector< const SMDS_MeshNode*> GetCommonNodes(const SMDS_MeshElement* e1,
+  SMESHUtils_EXPORT
+  std::vector< const SMDS_MeshNode*> GetCommonNodes(const SMDS_MeshElement* e1,
                                                     const SMDS_MeshElement* e2);
 
   /*!
    * \brief Return SMESH_NodeSearcher. The caller is responsible for deleteing it
    */
-  SMESHUtils_EXPORT SMESH_NodeSearcher* GetNodeSearcher( SMDS_Mesh& mesh );
+  SMESHUtils_EXPORT
+  SMESH_NodeSearcher* GetNodeSearcher( SMDS_Mesh& mesh );
 
   /*!
    * \brief Return SMESH_ElementSearcher. The caller is responsible for deleting it
    */
-  SMESHUtils_EXPORT SMESH_ElementSearcher* GetElementSearcher( SMDS_Mesh& mesh );
-  SMESHUtils_EXPORT SMESH_ElementSearcher* GetElementSearcher( SMDS_Mesh& mesh,
-                                                               SMDS_ElemIteratorPtr elemIt );
+  SMESHUtils_EXPORT
+  SMESH_ElementSearcher* GetElementSearcher( SMDS_Mesh& mesh );
+  SMESHUtils_EXPORT
+  SMESH_ElementSearcher* GetElementSearcher( SMDS_Mesh& mesh,
+                                             SMDS_ElemIteratorPtr elemIt );
 }
 
 #endif
index 2555cc3..6e5bda2 100644 (file)
 
 #include "SMESH_Utils.hxx"
 
-#include <SMDS_MeshNode.hxx>
+#include "SMDS_MeshNode.hxx"
 
 #include <gp_XYZ.hxx>
+#include <gp_XY.hxx>
 
 #include <map>
 #include <list>
@@ -42,7 +43,7 @@ typedef std::map<const SMDS_MeshElement*,
                  std::list<const SMDS_MeshElement*>, TIDCompare > TElemOfElemListMap;
 typedef std::map<const SMDS_MeshElement*,
                  std::list<const SMDS_MeshNode*>,    TIDCompare > TElemOfNodeListMap;
-typedef std::map<const SMDS_MeshNode*, const SMDS_MeshNode*> TNodeNodeMap;
+typedef std::map<const SMDS_MeshNode*, const SMDS_MeshNode*>      TNodeNodeMap;
 
 //!< Set of elements sorted by ID, to be used to assure predictability of edition
 typedef std::set< const SMDS_MeshElement*, TIDCompare >      TIDSortedElemSet;
@@ -50,8 +51,8 @@ typedef std::set< const SMDS_MeshNode*,    TIDCompare >      TIDSortedNodeSet;
 
 typedef std::pair< const SMDS_MeshNode*, const SMDS_MeshNode* >   NLink;
 
-struct faceQuadStruct; // defined in StdMeshers_Quadrangle_2D.hxx
-typedef boost::shared_ptr<faceQuadStruct> TFaceQuadStructPtr;
+struct FaceQuadStruct; // defined in StdMeshers_Quadrangle_2D.hxx
+typedef boost::shared_ptr<FaceQuadStruct> TFaceQuadStructPtr;
 
 
 namespace SMESHUtils
@@ -137,6 +138,10 @@ typedef struct uvPtStruct
   double x, y; // 2d parameter, normalized [0,1]
   const SMDS_MeshNode * node;
 
+  uvPtStruct(): node(NULL) {}
+
+  inline gp_XY UV() const { return gp_XY( u, v ); }
+
   struct NodeAccessor // accessor to iterate on nodes in UVPtStructVec
   {
     static const SMDS_MeshNode* value(std::vector< uvPtStruct >::const_iterator it)
index 929da0d..7a59ad8 100644 (file)
@@ -2075,6 +2075,7 @@ bool _pyMesh::NeedMeshAccess( const Handle(_pyCommand)& theCommand )
         "GetSubMeshElementsId","GetSubMeshNodesId","GetSubMeshElementType","Dump","GetNodeXYZ",
         "GetNodeInverseElements","GetShapeID","GetShapeIDForElem","GetElemNbNodes",
         "GetElemNode","IsMediumNode","IsMediumNodeOfAnyElem","ElemNbEdges","ElemNbFaces",
+        "GetElemFaceNodes", "GetFaceNormal", "FindElementByNodes",
         "IsPoly","IsQuadratic","BaryCenter","GetHypothesisList", "SetAutoColor", "GetAutoColor",
         "Clear", "ConvertToStandalone", "GetMeshOrder", "SetMeshOrder"
         ,"" }; // <- mark of end
index 43ffd02..78d9bf9 100644 (file)
@@ -240,13 +240,20 @@ namespace SMESH
   template<class TArray>
   void DumpArray(const TArray& theArray, TPythonDump & theStream)
   {
-    theStream << "[ ";
-    for (int i = 1; i <= theArray.length(); i++) {
-      theStream << theArray[i-1];
-      if ( i < theArray.length() )
-        theStream << ", ";
+    if ( theArray.length() == 0 )
+    {
+      theStream << "[]";
+    }
+    else
+    {
+      theStream << "[ ";
+      for (int i = 1; i <= theArray.length(); i++) {
+        theStream << theArray[i-1];
+        if ( i < theArray.length() )
+          theStream << ", ";
+      }
+      theStream << " ]";
     }
-    theStream << " ]";
   }
 
   TPythonDump&
@@ -264,6 +271,13 @@ namespace SMESH
   }
 
   TPythonDump&
+  TPythonDump::operator<<(const SMESH::nodes_array& theArg)
+  {
+    DumpArray( theArg, *this );
+    return *this;
+  }
+
+  TPythonDump&
   TPythonDump::operator<<(const SMESH::string_array& theArray)
   {
     myStream << "[ ";
@@ -517,6 +531,11 @@ namespace SMESH
     DumpArray( *theList, *this );
     return *this;
   }
+  TPythonDump& TPythonDump::operator<<(const GEOM::ListOfGO& theList)
+  {
+    DumpArray( theList, *this );
+    return *this;
+  }
   TPythonDump& TPythonDump::operator<<(const SMESH::ListOfIDSources& theList)
   {
     DumpArray( theList, *this );
index 812d3c0..f6beb92 100644 (file)
@@ -44,6 +44,7 @@
 #include "SMESH_Gen_i.hxx"
 #include "SMESH_Group.hxx"
 #include "SMESH_Group_i.hxx"
+#include "SMESH_MeshAlgos.hxx"
 #include "SMESH_MeshEditor.hxx"
 #include "SMESH_MeshEditor_i.hxx"
 #include "SMESH_MeshPartDS.hxx"
@@ -4024,6 +4025,32 @@ SMESH::long_array* SMESH_Mesh_i::GetElemFaceNodes(CORBA::Long  elemId,
 }
 
 //=======================================================================
+//function : GetElemFaceNodes
+//purpose  : Returns three components of normal of given mesh face.
+//=======================================================================
+
+SMESH::double_array* SMESH_Mesh_i::GetFaceNormal(CORBA::Long  elemId)
+{
+  if ( _preMeshInfo )
+    _preMeshInfo->FullLoadFromFile();
+
+  SMESH::double_array_var aResult = new SMESH::double_array();
+
+  if ( SMESHDS_Mesh* mesh = _impl->GetMeshDS() )
+  {
+    gp_XYZ normal;
+    if ( SMESH_MeshAlgos::FaceNormal( mesh->FindElement(elemId), normal, /*normalized=*/true ))
+    {
+      aResult->length( 3 );
+      aResult[ 0 ] = normal.X();
+      aResult[ 1 ] = normal.Y();
+      aResult[ 2 ] = normal.Z();
+    }
+  }
+  return aResult._retn();
+}
+
+//=======================================================================
 //function : FindElementByNodes
 //purpose  : Returns an element based on all given nodes.
 //=======================================================================
index c52699b..2615bef 100644 (file)
@@ -533,6 +533,11 @@ public:
   SMESH::long_array* GetElemFaceNodes(CORBA::Long elemId, CORBA::Short faceIndex);
 
   /*!
+   * Returns three components of normal of given mesh face (or an empty array in KO case)
+   */
+  SMESH::double_array* GetFaceNormal(CORBA::Long faceId);
+
+  /*!
    * Returns an element based on all given nodes.
    */
   CORBA::Long FindElementByNodes(const SMESH::long_array& nodes);
index 595e732..1f9d5a6 100644 (file)
@@ -27,6 +27,7 @@
 
 #include <SALOMEconfig.h>
 #include CORBA_SERVER_HEADER(SMESH_Mesh)
+#include CORBA_SERVER_HEADER(GEOM_Gen)
 #include CORBA_SERVER_HEADER(SALOMEDS)
 
 #include <TCollection_AsciiString.hxx>
@@ -163,6 +164,9 @@ namespace SMESH
     operator<<(const SMESH::string_array& theArg);
 
     TPythonDump&
+    operator<<(const SMESH::nodes_array& theArg);
+
+    TPythonDump&
     operator<<(SMESH::SMESH_Hypothesis_ptr theArg);
 
     TPythonDump&
@@ -217,6 +221,9 @@ namespace SMESH
     operator<<(const SMESH::ListOfGroups * theList);
 
     TPythonDump&
+    operator<<(const GEOM::ListOfGO& theList);
+
+    TPythonDump&
     operator<<(const SMESH::ListOfIDSources& theList);
 
     static const char* SMESHGenName() { return "smeshgen"; }
index e28659b..4e5b2d2 100644 (file)
@@ -199,7 +199,8 @@ class StdMeshersBuilder_Segment(Mesh_Algorithm):
         hyp.SetDeflection(deflection)
         return hyp
 
-    ## Defines "Arithmetic1D" hypothesis to cut an edge in several segments with increasing arithmetic length
+    ## Defines "Arithmetic1D" hypothesis to cut an edge in several segments with a length
+    #  that changes in arithmetic progression
     #  @param start defines the length of the first segment
     #  @param end   defines the length of the last  segment
     #  @param reversedEdges is a list of edges to mesh using reversed orientation.
@@ -226,6 +227,32 @@ class StdMeshersBuilder_Segment(Mesh_Algorithm):
         hyp.SetObjectEntry( entry )
         return hyp
 
+    ## Defines "GeometricProgression" hypothesis to cut an edge in several
+    #  segments with a length that changes in Geometric progression
+    #  @param start defines the length of the first segment
+    #  @param ratio defines the common ratio of the geometric progression
+    #  @param reversedEdges is a list of edges to mesh using reversed orientation.
+    #                       A list item can also be a tuple (edge, 1st_vertex_of_edge)
+    #  @param UseExisting if ==true - searches for an existing hypothesis created with
+    #                     the same parameters, else (default) - creates a new one
+    #  @return an instance of StdMeshers_Geometric1D hypothesis
+    #  @ingroup l3_hypos_1dhyps
+    def GeometricProgression(self, start, ratio, reversedEdges=[], UseExisting=0):
+        reversedEdgeInd = self.ReversedEdgeIndices(reversedEdges)
+        entry = self.MainShapeEntry()
+        from salome.smesh.smeshBuilder import IsEqual
+        compFun = lambda hyp, args: ( IsEqual(hyp.GetLength(1), args[0]) and \
+                                      IsEqual(hyp.GetLength(0), args[1]) and \
+                                      hyp.GetReversedEdges() == args[2]  and \
+                                      (not args[2] or hyp.GetObjectEntry() == args[3]))
+        hyp = self.Hypothesis("GeometricProgression", [start, ratio, reversedEdgeInd, entry],
+                              UseExisting=UseExisting, CompareMethod=compFun)
+        hyp.SetStartLength( start )
+        hyp.SetCommonRatio( ratio )
+        hyp.SetReversedEdges( reversedEdgeInd )
+        hyp.SetObjectEntry( entry )
+        return hyp
+
     ## Defines "FixedPoints1D" hypothesis to cut an edge using parameter
     # on curve from 0 to 1 (additionally it is neecessary to check
     # orientation of edges and create list of reversed edges if it is
@@ -237,7 +264,7 @@ class StdMeshersBuilder_Segment(Mesh_Algorithm):
     #                       A list item can also be a tuple (edge, 1st_vertex_of_edge)
     #  @param UseExisting if ==true - searches for an existing hypothesis created with
     #                     the same parameters, else (default) - creates a new one
-    #  @return an instance of StdMeshers_Arithmetic1D hypothesis
+    #  @return an instance of StdMeshers_FixedPoints1D hypothesis
     #  @ingroup l3_hypos_1dhyps
     def FixedPoints1D(self, points, nbSegs=[1], reversedEdges=[], UseExisting=0):
         if not isinstance(reversedEdges,list): #old version script, before adding reversedEdges
@@ -295,12 +322,23 @@ class StdMeshersBuilder_Segment(Mesh_Algorithm):
         hyp.SetDeflection(d)
         return hyp
 
-    ## Defines "Propagation" hypothesis that propagates all other hypotheses on all other edges that are at
-    #  the opposite side in case of quadrangular faces
+    ## Defines "Propagation" hypothesis that propagates 1D hypotheses
+    #  from an edge where this hypothesis is assigned to
+    #  on all other edges that are at the opposite side in case of quadrangular faces
+    #  This hypothesis should be assigned to an edge to propagate a hypothesis from.
     #  @ingroup l3_hypos_additi
     def Propagation(self):
         return self.Hypothesis("Propagation", UseExisting=1, CompareMethod=self.CompareEqualHyp)
 
+    ## Defines "Propagation of Node Distribution" hypothesis that propagates
+    #  distribution of nodes from an edge where this hypothesis is assigned to,
+    #  to opposite edges of quadrangular faces, so that number of segments on all these
+    #  edges will be the same, as well as relations between segment lengths. 
+    #  @ingroup l3_hypos_additi
+    def PropagationOfDistribution(self):
+        return self.Hypothesis("PropagOfDistribution", UseExisting=1,
+                               CompareMethod=self.CompareEqualHyp)
+
     ## Defines "AutomaticLength" hypothesis
     #  @param fineness for the fineness [0-1]
     #  @param UseExisting if ==true - searches for an existing hypothesis created with the
@@ -552,28 +590,56 @@ class StdMeshersBuilder_Quadrangle(Mesh_Algorithm):
     #                    same number of segments, the other pair must have an even difference
     #                    between the numbers of segments on the sides.
     #  @param triangleVertex: vertex of a trilateral geometrical face, around which triangles
-    #                  will be created while other elements will be quadrangles.
-    #                  Vertex can be either a GEOM_Object or a vertex ID within the
-    #                  shape to mesh
-    #  @param UseExisting: if ==true - searches for the existing hypothesis created with
-    #                  the same parameters, else (default) - creates a new one
+    #                    will be created while other elements will be quadrangles.
+    #                    Vertex can be either a GEOM_Object or a vertex ID within the
+    #                    shape to mesh
+    #  @param enfVertices: list of shapes defining positions where nodes (enforced nodes)
+    #                    must be created by the mesher. Shapes can be of any type,
+    #                    vertices of given shapes define positions of enforced nodes.
+    #                    Only vertices successfully projected to the face are used.
+    #  @param enfPoint: list of points giving positions of enforced nodes.
+    #                    Point can be defined either as SMESH.PointStruct's
+    #                    ([SMESH.PointStruct(x1,y1,z1), SMESH.PointStruct(x2,y2,z2),...])
+    #                    or triples of values ([[x1,y1,z1], [x2,y2,z2], ...]).
+    #                    In the case if the defined QuadrangleParameters() refer to a sole face,
+    #                    all given points must lie on this face, else the mesher fails.
+    #  @param UseExisting: if \c True - searches for the existing hypothesis created with
+    #                    the same parameters, else (default) - creates a new one
     #  @ingroup l3_hypos_quad
-    def QuadrangleParameters(self, quadType=StdMeshers.QUAD_STANDARD, triangleVertex=0, UseExisting=0):
-        import GEOM
+    def QuadrangleParameters(self, quadType=StdMeshers.QUAD_STANDARD, triangleVertex=0,
+                             enfVertices=[],enfPoints=[],UseExisting=0):
+        import GEOM, SMESH
         vertexID = triangleVertex
         if isinstance( triangleVertex, GEOM._objref_GEOM_Object ):
             vertexID = self.mesh.geompyD.GetSubShapeID( self.mesh.geom, triangleVertex )
+        if isinstance( enfVertices, int ) and not enfPoints and not UseExisting:
+            # a call of old syntax, before inserting enfVertices and enfPoints before UseExisting
+            UseExisting, enfVertices = enfVertices, []
+        pStructs, xyz = [], []
+        for p in enfPoints:
+            if isinstance( p, SMESH.PointStruct ):
+                xyz.append(( p.x, p.y, p.z ))
+                pStructs.append( p )
+            else:
+                xyz.append(( p[0], p[1], p[2] ))
+                pStructs.append( SMESH.PointStruct( p[0], p[1], p[2] ))
         if not self.params:
             compFun = lambda hyp,args: \
                       hyp.GetQuadType() == args[0] and \
-                      ( hyp.GetTriaVertex()==args[1] or ( hyp.GetTriaVertex()<1 and args[1]<1))
-            self.params = self.Hypothesis("QuadrangleParams", [quadType,vertexID],
+                      (hyp.GetTriaVertex()==args[1] or ( hyp.GetTriaVertex()<1 and args[1]<1)) and \
+                      ((hyp.GetEnforcedNodes()) == (args[2],args[3])) # True w/o enfVertices only
+            entries = [ shape.GetStudyEntry() for shape in enfVertices ]
+            self.params = self.Hypothesis("QuadrangleParams", [quadType,vertexID,entries,xyz],
                                           UseExisting = UseExisting, CompareMethod=compFun)
             pass
         if self.params.GetQuadType() != quadType:
             self.params.SetQuadType(quadType)
         if vertexID > 0:
             self.params.SetTriaVertex( vertexID )
+        from salome.smesh.smeshBuilder import AssureGeomPublished
+        for v in enfVertices:
+            AssureGeomPublished( self.mesh, v )
+        self.params.SetEnforcedNodes( enfVertices, pStructs )
         return self.params
 
     ## Defines "QuadrangleParams" hypothesis with a type of quadrangulation that only
@@ -980,7 +1046,8 @@ class StdMeshersBuilder_Prism3D(Mesh_Algorithm):
         return hyp
 
     ## Defines "Arithmetic1D" hypothesis, specifying the distribution of segments
-    #  to build between the inner and the outer shells with a length that changes in arithmetic progression
+    #  to build between the inner and the outer shells with a length that changes
+    #  in arithmetic progression
     #  @param start  the length of the first segment
     #  @param end    the length of the last  segment
     def Arithmetic1D(self, start, end ):
@@ -992,6 +1059,20 @@ class StdMeshersBuilder_Prism3D(Mesh_Algorithm):
         hyp.SetLength(end  , 0)
         return hyp
 
+    ## Defines "GeometricProgression" hypothesis, specifying the distribution of segments
+    #  to build between the inner and the outer shells with a length that changes
+    #  in Geometric progression
+    #  @param start  the length of the first segment
+    #  @param ratio  the common ratio of the geometric progression
+    def GeometricProgression(self, start, ratio ):
+        if self.algoType != "RadialPrism_3D":
+            print "Prism_3D algorith doesn't support any hyposesis"
+            return None
+        hyp = self.OwnHypothesis("GeometricProgression", [start, ratio])
+        hyp.SetStartLength( start )
+        hyp.SetCommonRatio( ratio )
+        return hyp
+
     ## Defines "StartEndLength" hypothesis, specifying distribution of segments
     #  to build between the inner and the outer shells as geometric length increasing
     #  @param start for the length of the first segment
@@ -1148,6 +1229,16 @@ class StdMeshersBuilder_RadialQuadrangle1D2D(Mesh_Algorithm):
         hyp.SetLength(end  , 0)
         return hyp
 
+    ## Defines "GeometricProgression" hypothesis, specifying the distribution of segments
+    #  with a length that changes in Geometric progression
+    #  @param start  the length of the first segment
+    #  @param ratio  the common ratio of the geometric progression
+    def GeometricProgression(self, start, ratio ):
+        hyp = self.OwnHypothesis("GeometricProgression", [start, ratio])
+        hyp.SetStartLength( start )
+        hyp.SetCommonRatio( ratio )
+        return hyp
+
     ## Defines "StartEndLength" hypothesis, specifying distribution of segments
     #  as geometric length increasing
     #  @param start for the length of the first segment
index 5e68dc2..a5aae93 100644 (file)
@@ -2418,6 +2418,12 @@ class Mesh:
     def GetElemFaceNodes(self,elemId, faceIndex):
         return self.mesh.GetElemFaceNodes(elemId, faceIndex)
 
+    ## Returns three components of normal of given mesh face
+    #  (or an empty array in KO case)
+    #  @ingroup l1_meshinfo
+    def GetFaceNormal(self, faceId):
+        return self.mesh.GetFaceNormal(faceId)
+
     ## Returns an element based on all given nodes.
     #  @ingroup l1_meshinfo
     def FindElementByNodes(self,nodes):
@@ -3181,7 +3187,8 @@ class Mesh:
     #  Note that nodes built on edges and boundary nodes are always fixed.
     #  @param MaxNbOfIterations the maximum number of iterations
     #  @param MaxAspectRatio varies in range [1.0, inf]
-    #  @param Method is Laplacian(LAPLACIAN_SMOOTH) or Centroidal(CENTROIDAL_SMOOTH)
+    #  @param Method is either Laplacian (SMESH.SMESH_MeshEditor.LAPLACIAN_SMOOTH)
+    #         or Centroidal (SMESH.SMESH_MeshEditor.CENTROIDAL_SMOOTH)
     #  @return TRUE in case of success, FALSE otherwise.
     #  @ingroup l2_modif_smooth
     def Smooth(self, IDsOfElements, IDsOfFixedNodes,
@@ -3199,7 +3206,8 @@ class Mesh:
     #  Note that nodes built on edges and boundary nodes are always fixed.
     #  @param MaxNbOfIterations the maximum number of iterations
     #  @param MaxAspectRatio varies in range [1.0, inf]
-    #  @param Method is Laplacian(LAPLACIAN_SMOOTH) or Centroidal(CENTROIDAL_SMOOTH)
+    #  @param Method is either Laplacian (SMESH.SMESH_MeshEditor.LAPLACIAN_SMOOTH)
+    #         or Centroidal (SMESH.SMESH_MeshEditor.CENTROIDAL_SMOOTH)
     #  @return TRUE in case of success, FALSE otherwise.
     #  @ingroup l2_modif_smooth
     def SmoothObject(self, theObject, IDsOfFixedNodes,
@@ -3215,7 +3223,8 @@ class Mesh:
     #  Note that nodes built on edges and boundary nodes are always fixed.
     #  @param MaxNbOfIterations the maximum number of iterations
     #  @param MaxAspectRatio varies in range [1.0, inf]
-    #  @param Method is Laplacian(LAPLACIAN_SMOOTH) or Centroidal(CENTROIDAL_SMOOTH)
+    #  @param Method is either Laplacian (SMESH.SMESH_MeshEditor.LAPLACIAN_SMOOTH)
+    #         or Centroidal (SMESH.SMESH_MeshEditor.CENTROIDAL_SMOOTH)
     #  @return TRUE in case of success, FALSE otherwise.
     #  @ingroup l2_modif_smooth
     def SmoothParametric(self, IDsOfElements, IDsOfFixedNodes,
@@ -3233,7 +3242,8 @@ class Mesh:
     #  Note that nodes built on edges and boundary nodes are always fixed.
     #  @param MaxNbOfIterations the maximum number of iterations
     #  @param MaxAspectRatio varies in range [1.0, inf]
-    #  @param Method Laplacian(LAPLACIAN_SMOOTH) or Centroidal(CENTROIDAL_SMOOTH)
+    #  @param Method is either Laplacian (SMESH.SMESH_MeshEditor.LAPLACIAN_SMOOTH)
+    #         or Centroidal (SMESH.SMESH_MeshEditor.CENTROIDAL_SMOOTH)
     #  @return TRUE in case of success, FALSE otherwise.
     #  @ingroup l2_modif_smooth
     def SmoothParametricObject(self, theObject, IDsOfFixedNodes,
index a81f825..0d018d5 100644 (file)
@@ -267,21 +267,26 @@ class Mesh_Algorithm:
     #  @param thickness total thickness of layers of prisms
     #  @param numberOfLayers number of layers of prisms
     #  @param stretchFactor factor (>1.0) of growth of layer thickness towards inside of mesh
-    #  @param ignoreFaces list of geometrical faces (or their ids) not to generate layers on
+    #  @param faces list of geometrical faces (or their ids).
+    #         Viscous layers are either generated on these faces or not, depending on
+    #         the value of \a isFacesToIgnore parameter.
+    #  @param isFacesToIgnore if \c True, the Viscous layers are not generated on the
+    #         faces specified by the previous parameter (\a faces).
     #  @ingroup l3_hypos_additi
-    def ViscousLayers(self, thickness, numberOfLayers, stretchFactor, ignoreFaces=[]):
+    def ViscousLayers(self, thickness, numberOfLayers, stretchFactor,
+                      faces=[], isFacesToIgnore=True ):
         if not isinstance(self.algo, SMESH._objref_SMESH_3D_Algo):
             raise TypeError, "ViscousLayers are supported by 3D algorithms only"
         if not "ViscousLayers" in self.GetCompatibleHypothesis():
             raise TypeError, "ViscousLayers are not supported by %s"%self.algo.GetName()
-        if ignoreFaces and isinstance( ignoreFaces[0], geomBuilder.GEOM._objref_GEOM_Object ):
-            ignoreFaces = [ self.mesh.geompyD.GetSubShapeID(self.mesh.geom, f) for f in ignoreFaces ]
+        if faces and isinstance( faces[0], geomBuilder.GEOM._objref_GEOM_Object ):
+            faces = [ self.mesh.geompyD.GetSubShapeID(self.mesh.geom, f) for f in faces ]
         hyp = self.Hypothesis("ViscousLayers",
-                              [thickness, numberOfLayers, stretchFactor, ignoreFaces])
+                              [thickness, numberOfLayers, stretchFactor, faces])
         hyp.SetTotalThickness(thickness)
         hyp.SetNumberLayers(numberOfLayers)
         hyp.SetStretchFactor(stretchFactor)
-        hyp.SetIgnoreFaces(ignoreFaces)
+        hyp.SetFaces(faces, isFacesToIgnore)
         return hyp
 
     ## Defines "ViscousLayers2D" hypothesis to give parameters of layers of quadrilateral
@@ -290,9 +295,9 @@ class Mesh_Algorithm:
     #  @param thickness total thickness of layers of quadrilaterals
     #  @param numberOfLayers number of layers
     #  @param stretchFactor factor (>1.0) of growth of layer thickness towards inside of mesh
-    #  @param edges list of geometrical edge (or their ids).
+    #  @param edges list of geometrical edges (or their ids).
     #         Viscous layers are either generated on these edges or not, depending on
-    #         the values of \a isEdgesToIgnore parameter.
+    #         the value of \a isEdgesToIgnore parameter.
     #  @param isEdgesToIgnore if \c True, the Viscous layers are not generated on the
     #         edges specified by the previous parameter (\a edges).
     #  @ingroup l3_hypos_additi
@@ -313,7 +318,7 @@ class Mesh_Algorithm:
         hyp.SetEdges(edges, isEdgesToIgnore)
         return hyp
 
-    ## Transform a list of ether edges or tuples (edge, 1st_vertex_of_edge)
+    ## Transform a list of either edges or tuples (edge, 1st_vertex_of_edge)
     #  into a list acceptable to SetReversedEdges() of some 1D hypotheses
     #  @ingroup l3_hypos_1dhyps
     def ReversedEdgeIndices(self, reverseList):
index de84b50..2dde1b2 100644 (file)
@@ -75,8 +75,10 @@ ENDIF(SALOME_SMESH_ENABLE_MEFISTO)
 # header files / no moc processing
 SET(StdMeshers_HEADERS
   StdMeshers_LocalLength.hxx
+  StdMeshers_Reversible1D.hxx
   StdMeshers_StartEndLength.hxx
   StdMeshers_Arithmetic1D.hxx
+  StdMeshers_Geometric1D.hxx
   StdMeshers_FixedPoints1D.hxx
   StdMeshers_NumberOfSegments.hxx
   StdMeshers_Deflection1D.hxx
@@ -136,8 +138,10 @@ ENDIF(SALOME_SMESH_ENABLE_MEFISTO)
 # sources / static
 SET(StdMeshers_SOURCES
   StdMeshers_LocalLength.cxx
+  StdMeshers_Reversible1D.cxx
   StdMeshers_StartEndLength.cxx
   StdMeshers_Arithmetic1D.cxx
+  StdMeshers_Geometric1D.cxx
   StdMeshers_FixedPoints1D.cxx
   StdMeshers_NumberOfSegments.cxx
   StdMeshers_Deflection1D.cxx
index fa4aa90..369d4e4 100644 (file)
 
 #include "utilities.h"
 
-#include <Precision.hxx>
-#include <Bnd_Box.hxx>
-
 #include <limits>
 
+#include <Bnd_Box.hxx>
+#include <Precision.hxx>
+#include <gp_Vec.hxx>
+
 using namespace std;
 
 //=======================================================================
@@ -48,10 +49,23 @@ StdMeshers_CartesianParameters3D::StdMeshers_CartesianParameters3D(int         h
                                                                    int         studyId,
                                                                    SMESH_Gen * gen)
   : SMESH_Hypothesis(hypId, studyId, gen),
-    _sizeThreshold( 4.0 ) // default according to the customer specification
+    _sizeThreshold( 4.0 ), // default according to the customer specification
+    _toAddEdges( false )
 {
   _name = "CartesianParameters3D"; // used by "Cartesian_3D"
   _param_algo_dim = 3; // 3D
+
+  _axisDirs[0] = 1.;
+  _axisDirs[1] = 0.;
+  _axisDirs[2] = 0.;
+
+  _axisDirs[3] = 0.;
+  _axisDirs[4] = 1.;
+  _axisDirs[5] = 0.;
+
+  _axisDirs[6] = 0.;
+  _axisDirs[7] = 0.;
+  _axisDirs[8] = 1.;
 }
 
 
@@ -309,6 +323,44 @@ void StdMeshers_CartesianParameters3D::GetCoordinates(std::vector<double>& xNode
 }
 
 //=======================================================================
+//function : SetAxisDirs
+//purpose  : Sets directions of axes
+//=======================================================================
+
+void StdMeshers_CartesianParameters3D::SetAxisDirs(const double* the9DirComps)
+  throw ( SALOME_Exception )
+{
+  gp_Vec x( the9DirComps[0],
+            the9DirComps[1],
+            the9DirComps[2] );
+  gp_Vec y( the9DirComps[3],
+            the9DirComps[4],
+            the9DirComps[5] );
+  gp_Vec z( the9DirComps[6],
+            the9DirComps[7],
+            the9DirComps[8] );
+  if ( x.Magnitude() < RealSmall() ||
+       y.Magnitude() < RealSmall() ||
+       z.Magnitude() < RealSmall() )
+    throw SALOME_Exception("Zero magnitude of axis direction");
+
+  if ( x.IsParallel( y, M_PI / 180. ) ||
+       x.IsParallel( z, M_PI / 180. ) ||
+       y.IsParallel( z, M_PI / 180. ))
+    throw SALOME_Exception("Parallel axis directions");
+
+  bool isChanged = false;
+  for ( int i = 0; i < 9; ++i )
+  {
+    if ( Abs( _axisDirs[i] - the9DirComps[i] ) > 1e-7 )
+      isChanged = true;
+    _axisDirs[i] = the9DirComps[i];
+  }
+  if ( isChanged )
+    NotifySubMeshesHypothesisModification();
+}
+
+//=======================================================================
 //function : GetGrid
 //purpose  : Return coordinates of node positions along the three axes
 //=======================================================================
@@ -333,6 +385,33 @@ double StdMeshers_CartesianParameters3D::GetSizeThreshold() const
 }
 
 //=======================================================================
+//function : SetToAddEdges
+//purpose  : Enables implementation of geometrical edges into the mesh. If this feature
+//           is disabled, sharp edges of the shape are lost ("smoothed") in the mesh if
+//           they don't coincide with the grid lines
+//=======================================================================
+
+void StdMeshers_CartesianParameters3D::SetToAddEdges(bool toAdd)
+{
+  if ( _toAddEdges != toAdd )
+  {
+    _toAddEdges = toAdd;
+    NotifySubMeshesHypothesisModification();
+  }
+}
+
+//=======================================================================
+//function : GetToAddEdges
+//purpose  : Returns true if implementation of geometrical edges into the
+//           mesh is enabled
+//=======================================================================
+
+bool StdMeshers_CartesianParameters3D::GetToAddEdges() const
+{
+  return _toAddEdges;
+}
+
+//=======================================================================
 //function : IsDefined
 //purpose  : Return true if parameters are well defined
 //=======================================================================
@@ -369,6 +448,7 @@ std::ostream & StdMeshers_CartesianParameters3D::SaveTo(std::ostream & save)
     for ( size_t j = 0; j < _spaceFunctions[i].size(); ++j )
       save << _spaceFunctions[i][j] << " ";
   }
+  save << _toAddEdges << " ";
 
   return save;
 }
@@ -382,7 +462,7 @@ std::istream & StdMeshers_CartesianParameters3D::LoadFrom(std::istream & load)
 {
   bool ok;
 
-  ok = (load >> _sizeThreshold  );
+  ok = ( load >> _sizeThreshold );
   for ( int ax = 0; ax < 3; ++ax )
   {
     if (ok)
@@ -419,6 +499,9 @@ std::istream & StdMeshers_CartesianParameters3D::LoadFrom(std::istream & load)
       }
     }
   }
+
+  load >> _toAddEdges;
+
   return load;
 }
 
index adb22f0..9897095 100644 (file)
@@ -100,6 +100,10 @@ public:
                       std::vector<double>& yNodes,
                       std::vector<double>& zNodes,
                       const Bnd_Box&       bndBox) const throw ( SALOME_Exception );
+
+  void SetAxisDirs(const double* the9DirComps) throw ( SALOME_Exception );
+  const double* GetAxisDirs() const { return _axisDirs; }
+
   /*!
    * Set size threshold. A polyhedral cell got by cutting an initial
    * hexahedron by geometry boundary is considered small and is removed if
@@ -112,6 +116,14 @@ public:
   double GetSizeThreshold() const;
 
   /*!
+   * \brief Enables implementation of geometrical edges into the mesh. If this feature
+   *        is disabled, sharp edges of the shape are lost ("smoothed") in the mesh if
+   *        they don't coincide with the grid lines
+   */
+  void SetToAddEdges(bool toAdd);
+  bool GetToAddEdges() const;
+
+  /*!
    * \brief Return true if parameters are well defined
    */
   bool IsDefined() const;
@@ -138,7 +150,10 @@ public:
   std::vector<std::string> _spaceFunctions[3];
   std::vector<double>      _internalPoints[3];
 
+  double _axisDirs[9];
+
   double _sizeThreshold;
+  bool   _toAddEdges;
 };
 
 #endif
index 83447b1..36b2e2d 100644 (file)
 #include "SMESH_subMeshEventListener.hxx"
 #include "StdMeshers_CartesianParameters3D.hxx"
 
-#include "utilities.h"
-#include "Utils_ExceptHandlers.hxx"
+#include <utilities.h>
+#include <Utils_ExceptHandlers.hxx>
 #include <Basics_OCCTVersion.hxx>
 
+#include <BRepAdaptor_Curve.hxx>
 #include <BRepAdaptor_Surface.hxx>
 #include <BRepBndLib.hxx>
 #include <BRepBuilderAPI_Copy.hxx>
@@ -44,6 +45,7 @@
 #include <BRep_Tool.hxx>
 #include <Bnd_Box.hxx>
 #include <ElSLib.hxx>
+#include <GCPnts_UniformDeflection.hxx>
 #include <Geom2d_BSplineCurve.hxx>
 #include <Geom2d_BezierCurve.hxx>
 #include <Geom2d_TrimmedCurve.hxx>
@@ -63,7 +65,6 @@
 #include <TopExp.hxx>
 #include <TopExp_Explorer.hxx>
 #include <TopLoc_Location.hxx>
-#include <TopTools_MapIteratorOfMapOfShape.hxx>
 #include <TopTools_MapOfShape.hxx>
 #include <TopoDS.hxx>
 #include <TopoDS_Face.hxx>
@@ -76,7 +77,7 @@
 #include <gp_Sphere.hxx>
 #include <gp_Torus.hxx>
 
-//#undef WITH_TBB
+#undef WITH_TBB
 #ifdef WITH_TBB
 #include <tbb/parallel_for.h>
 //#include <tbb/enumerable_thread_specific.h>
@@ -150,6 +151,8 @@ bool StdMeshers_Cartesian_3D::CheckHypothesis (SMESH_Mesh&          aMesh,
 
 namespace
 {
+  typedef int TGeomID;
+
   //=============================================================================
   // Definitions of internal utils
   // --------------------------------------------------------------------------
@@ -161,17 +164,40 @@ namespace
   };
   // --------------------------------------------------------------------------
   /*!
-   * \brief Data of intersection between a GridLine and a TopoDS_Face
+   * \brief Common data of any intersection between a Grid and a shape
    */
-  struct IntersectionPoint
+  struct B_IntersectPoint
   {
-    double                       _paramOnLine;
-    mutable Transition           _transition;
     mutable const SMDS_MeshNode* _node;
-    mutable size_t               _indexOnLine;
+    mutable vector< TGeomID >    _faceIDs;
+
+    B_IntersectPoint(): _node(NULL) {}
+    void Add( const vector< TGeomID >& fIDs, const SMDS_MeshNode* n=0 ) const;
+    bool HasCommonFace( const B_IntersectPoint * other ) const;
+    bool IsOnFace( int faceID ) const;
+    virtual ~B_IntersectPoint() {}
+  };
+  // --------------------------------------------------------------------------
+  /*!
+   * \brief Data of intersection between a GridLine and a TopoDS_Face
+   */
+  struct F_IntersectPoint : public B_IntersectPoint
+  {
+    double             _paramOnLine;
+    mutable Transition _transition;
+    mutable size_t     _indexOnLine;
 
-    IntersectionPoint(): _node(0) {}
-    bool operator< ( const IntersectionPoint& o ) const { return _paramOnLine < o._paramOnLine; }
+    bool operator< ( const F_IntersectPoint& o ) const { return _paramOnLine < o._paramOnLine; }
+  };
+  // --------------------------------------------------------------------------
+  /*!
+   * \brief Data of intersection between GridPlanes and a TopoDS_EDGE
+   */
+  struct E_IntersectPoint : public B_IntersectPoint
+  {
+    gp_Pnt  _point;
+    double  _uvw[3];
+    TGeomID _shapeID;
   };
   // --------------------------------------------------------------------------
   /*!
@@ -181,10 +207,23 @@ namespace
   {
     gp_Lin _line;
     double _length; // line length
-    multiset< IntersectionPoint > _intPoints;
+    multiset< F_IntersectPoint > _intPoints;
 
     void RemoveExcessIntPoints( const double tol );
-    bool GetIsOutBefore( multiset< IntersectionPoint >::iterator ip, bool prevIsOut );
+    bool GetIsOutBefore( multiset< F_IntersectPoint >::iterator ip, bool prevIsOut );
+  };
+  // --------------------------------------------------------------------------
+  /*!
+   * \brief Planes of the grid used to find intersections of an EDGE with a hexahedron
+   */
+  struct GridPlanes
+  {
+    double _factor;
+    gp_XYZ _uNorm, _vNorm, _zNorm;
+    vector< gp_XYZ > _origins; // origin points of all planes in one direction
+    vector< double > _zProjs;  // projections of origins to _zNorm
+
+    gp_XY GetUV( const gp_Pnt& p, const gp_Pnt& origin );
   };
   // --------------------------------------------------------------------------
   /*!
@@ -234,11 +273,15 @@ namespace
   struct Grid
   {
     vector< double >   _coords[3]; // coordinates of grid nodes
+    gp_XYZ             _axes  [3]; // axis directions
     vector< GridLine > _lines [3]; //  in 3 directions
     double             _tol, _minCellSize;
 
-    vector< const SMDS_MeshNode* > _nodes; // mesh nodes at grid nodes
-    vector< bool >                 _isBndNode; // is mesh node at intersection with geometry
+    vector< const SMDS_MeshNode* >    _nodes; // mesh nodes at grid nodes
+    vector< const F_IntersectPoint* > _gridIntP; // grid node intersection with geometry
+
+    list< E_IntersectPoint >          _edgeIntP; // intersections with EDGEs
+    TopTools_IndexedMapOfShape        _shapes;
 
     size_t CellIndex( size_t i, size_t j, size_t k ) const
     {
@@ -257,6 +300,7 @@ namespace
     void SetCoordinates(const vector<double>& xCoords,
                         const vector<double>& yCoords,
                         const vector<double>& zCoords,
+                        const double*         axesDirs,
                         const TopoDS_Shape&   shape );
     void ComputeNodes(SMESH_MesherHelper& helper);
   };
@@ -302,10 +346,11 @@ namespace
   struct FaceGridIntersector
   {
     TopoDS_Face _face;
+    TGeomID     _faceID;
     Grid*       _grid;
     Bnd_Box     _bndBox;
     __IntCurvesFace_Intersector* _surfaceInt;
-    vector< std::pair< GridLine*, IntersectionPoint > > _intersections;
+    vector< std::pair< GridLine*, F_IntersectPoint > > _intersections;
 
     FaceGridIntersector(): _grid(0), _surfaceInt(0) {}
     void Intersect();
@@ -314,7 +359,12 @@ namespace
     void StoreIntersections()
     {
       for ( size_t i = 0; i < _intersections.size(); ++i )
-        _intersections[i].first->_intPoints.insert( _intersections[i].second );
+      {
+        multiset< F_IntersectPoint >::iterator ip = 
+          _intersections[i].first->_intPoints.insert( _intersections[i].second );
+        ip->_faceIDs.reserve( 1 );
+        ip->_faceIDs.push_back( _faceID );
+      }
     }
     const Bnd_Box& GetFaceBndBox()
     {
@@ -352,7 +402,7 @@ namespace
     gp_Torus    _torus;
     __IntCurvesFace_Intersector* _surfaceInt;
 
-    vector< IntersectionPoint > _intPoints;
+    vector< F_IntersectPoint > _intPoints;
 
     void IntersectWithPlane   (const GridLine& gridLine);
     void IntersectWithCylinder(const GridLine& gridLine);
@@ -381,20 +431,54 @@ namespace
     struct _Face;
     struct _Link;
     // --------------------------------------------------------------------------------
-    struct _Node //!< node either at a hexahedron corner or at GridLine intersection
+    struct _Node //!< node either at a hexahedron corner or at intersection
     {
-      const SMDS_MeshNode*     _node; // mesh node at hexahedron corner
-      const IntersectionPoint* _intPoint;
-
-      _Node(const SMDS_MeshNode* n=0, const IntersectionPoint* ip=0):_node(n), _intPoint(ip) {} 
-      const SMDS_MeshNode* Node() const { return _intPoint ? _intPoint->_node : _node; }
-      //bool IsCorner() const { return _node; }
+      const SMDS_MeshNode*       _node; // mesh node at hexahedron corner
+      const B_IntersectPoint* _intPoint;
+
+      _Node(const SMDS_MeshNode* n=0, const B_IntersectPoint* ip=0):_node(n), _intPoint(ip) {} 
+      const SMDS_MeshNode*    Node() const
+      { return ( _intPoint && _intPoint->_node ) ? _intPoint->_node : _node; }
+      const F_IntersectPoint* FaceIntPnt() const
+      { return static_cast< const F_IntersectPoint* >( _intPoint ); }
+      const E_IntersectPoint* EdgeIntPnt() const
+      { return static_cast< const E_IntersectPoint* >( _intPoint ); }
+      void Add( const E_IntersectPoint* ip )
+      {
+        if ( !_intPoint ) {
+          _intPoint = ip;
+        }
+        else if ( !_intPoint->_node ) {
+          ip->Add( _intPoint->_faceIDs );
+          _intPoint = ip;
+        }
+        else  {
+          _intPoint->Add( ip->_faceIDs );
+        }
+      }
+      bool IsLinked( const B_IntersectPoint* other ) const
+      {
+        return _intPoint && _intPoint->HasCommonFace( other );
+      }
+      bool IsOnFace( int faceID ) const // returns true if faceID is found
+      {
+        return _intPoint ? _intPoint->IsOnFace( faceID ) : false;
+      }
+      gp_Pnt Point() const
+      {
+        if ( const SMDS_MeshNode* n = Node() )
+          return SMESH_TNodeXYZ( n );
+        if ( const E_IntersectPoint* eip =
+             dynamic_cast< const E_IntersectPoint* >( _intPoint ))
+          return eip->_point;
+        return gp_Pnt( 1e100, 0, 0 );
+      }
     };
     // --------------------------------------------------------------------------------
     struct _Link // link connecting two _Node's
     {
       _Node* _nodes[2];
-      vector< _Node _intNodes; // _Node's at GridLine intersections
+      vector< _Node > _intNodes; // _Node's at GridLine intersections
       vector< _Link > _splits;
       vector< _Face*> _faces;
     };
@@ -412,21 +496,38 @@ namespace
       }
       _Node* FirstNode() const { return _link->_nodes[ _reverse ]; }
       _Node* LastNode() const { return _link->_nodes[ !_reverse ]; }
+      operator bool() const { return _link; }
+      vector< TGeomID > GetNotUsedFace(const set<TGeomID>& usedIDs ) const // returns a supporting FACEs
+      {
+        vector< TGeomID > faces;
+        const B_IntersectPoint *ip0, *ip1;
+        if (( ip0 = _link->_nodes[0]->_intPoint ) &&
+            ( ip1 = _link->_nodes[1]->_intPoint ))
+        {
+          for ( size_t i = 0; i < ip0->_faceIDs.size(); ++i )
+            if ( ip1->IsOnFace ( ip0->_faceIDs[i] ) &&
+                 !usedIDs.count( ip0->_faceIDs[i] ) )
+              faces.push_back( ip0->_faceIDs[i] );
+        }
+        return faces;
+      }
     };
     // --------------------------------------------------------------------------------
     struct _Face
     {
-      vector< _OrientedLink > _links;
-      vector< _Link >         _polyLinks; // links added to close a polygonal face
+      vector< _OrientedLink > _links;       // links on GridLine's
+      vector< _Link >         _polyLinks;   // links added to close a polygonal face
+      vector< _Node >         _edgeNodes;   // nodes at intersection with EDGEs
     };
     // --------------------------------------------------------------------------------
     struct _volumeDef // holder of nodes of a volume mesh element
     {
-      vector< const SMDS_MeshNode* > _nodes;
-      vector< int >                  _quantities;
+      //vector< const SMDS_MeshNode* > _nodes;
+      vector< _Node* > _nodes;
+      vector< int >    _quantities;
       typedef boost::shared_ptr<_volumeDef> Ptr;
-      void set( const vector< const SMDS_MeshNode* >& nodes,
-                const vector< int > quant = vector< int >() )
+      void set( const vector< _Node* >& nodes,
+                const vector< int >&    quant = vector< int >() )
       { _nodes = nodes; _quantities = quant; }
       // static Ptr New( const vector< const SMDS_MeshNode* >& nodes,
       //                 const vector< int > quant = vector< int >() )
@@ -440,13 +541,19 @@ namespace
 
     // topology of a hexahedron
     int   _nodeShift[8];
-    _Node _hexNodes[8];
-    _Link _hexLinks[12];
-    _Face _hexQuads[6];
+    _Node _hexNodes [8];
+    _Link _hexLinks [12];
+    _Face _hexQuads [6];
 
     // faces resulted from hexahedron intersection
     vector< _Face > _polygons;
 
+    // intresections with EDGEs
+    vector< const E_IntersectPoint* > _edgeIntPnts;
+
+    // nodes inside the hexahedron (at VERTEXes)
+    vector< _Node > _vertexNodes;
+
     // computed volume elements
     //vector< _volumeDef::Ptr > _volumeDefs;
     _volumeDef _volumeDefs;
@@ -459,7 +566,8 @@ namespace
 
   public:
     Hexahedron(const double sizeThreshold, Grid* grid);
-    int MakeElements(SMESH_MesherHelper& helper);
+    int MakeElements(SMESH_MesherHelper&                      helper,
+                     const map< TGeomID, vector< TGeomID > >& edge2faceIDsMap);
     void ComputeElements();
     void Init() { init( _i, _j, _k ); }
 
@@ -467,6 +575,17 @@ namespace
     Hexahedron(const Hexahedron& other );
     void init( size_t i, size_t j, size_t k );
     void init( size_t i );
+    void addEdges(SMESH_MesherHelper&                      helper,
+                  vector< Hexahedron* >&                   intersectedHex,
+                  const map< TGeomID, vector< TGeomID > >& edge2faceIDsMap);
+    gp_Pnt findIntPoint( double u1, double proj1, double u2, double proj2,
+                         double proj, BRepAdaptor_Curve& curve,
+                         const gp_XYZ& axis, const gp_XYZ& origin );
+    int  getEntity( const E_IntersectPoint* ip, int* facets, int& sub );
+    bool addIntersection( const E_IntersectPoint& ip,
+                          vector< Hexahedron* >&  hexes,
+                          int ijk[], int dIJK[] );
+    bool findChain( _Node* n1, _Node* n2, _Face& quad, vector<_Node*>& chainNodes );
     int  addElements(SMESH_MesherHelper& helper);
     bool isInHole() const;
     bool checkPolyhedronSize() const;
@@ -474,8 +593,17 @@ namespace
     bool addTetra();
     bool addPenta();
     bool addPyra ();
+    _Node* FindEqualNode( vector< _Node >&        nodes,
+                          const E_IntersectPoint* ip,
+                          const double            tol2 )
+    {
+      for ( size_t i = 0; i < nodes.size(); ++i )
+        if ( nodes[i].Point().SquareDistance( ip->_point ) <= tol2 )
+          return & nodes[i];
+      return 0;
+    }
   };
+
 #ifdef WITH_TBB
   // --------------------------------------------------------------------------
   /*!
@@ -507,11 +635,34 @@ namespace
         _faceVec[i].Intersect();
     }
   };
-
 #endif
+
   //=============================================================================
   // Implementation of internal utils
   //=============================================================================
+  /*!
+   * \brief adjust \a i to have \a val between values[i] and values[i+1]
+   */
+  inline void locateValue( int & i, double val, const vector<double>& values,
+                           int& di, double tol )
+  {
+    val += values[0]; // input \a val is measured from 0.
+    if ( i > values.size()-2 )
+      i = values.size()-2;
+    else
+      while ( i+2 < values.size() && val > values[ i+1 ])
+        ++i;
+    while ( i > 0 && val < values[ i ])
+      --i;
+
+    if ( i > 0 && val - values[ i ] < tol )
+      di = -1;
+    else if ( i+2 < values.size() && values[ i+1 ] - val < tol )
+      di = 1;
+    else
+      di = 0;
+  }
+  //=============================================================================
   /*
    * Remove coincident intersection points
    */
@@ -520,15 +671,16 @@ namespace
     if ( _intPoints.size() < 2 ) return;
 
     set< Transition > tranSet;
-    multiset< IntersectionPoint >::iterator ip1, ip2 = _intPoints.begin();
+    multiset< F_IntersectPoint >::iterator ip1, ip2 = _intPoints.begin();
     while ( ip2 != _intPoints.end() )
     {
       tranSet.clear();
       ip1 = ip2++;
-      while ( ip2->_paramOnLine - ip1->_paramOnLine <= tol  && ip2 != _intPoints.end())
+      while ( ip2 != _intPoints.end() && ip2->_paramOnLine - ip1->_paramOnLine <= tol )
       {
         tranSet.insert( ip1->_transition );
         tranSet.insert( ip2->_transition );
+        ip2->Add( ip1->_faceIDs );
         _intPoints.erase( ip1 );
         ip1 = ip2++;
       }
@@ -547,7 +699,7 @@ namespace
   /*
    * Return "is OUT" state for nodes before the given intersection point
    */
-  bool GridLine::GetIsOutBefore( multiset< IntersectionPoint >::iterator ip, bool prevIsOut )
+  bool GridLine::GetIsOutBefore( multiset< F_IntersectPoint >::iterator ip, bool prevIsOut )
   {
     if ( ip->_transition == Trans_IN )
       return true;
@@ -558,7 +710,7 @@ namespace
       // singularity point (apex of a cone)
       if ( _intPoints.size() == 1 || ip == _intPoints.begin() )
         return true;
-      multiset< IntersectionPoint >::iterator ipBef = ip, ipAft = ++ip;
+      multiset< F_IntersectPoint >::iterator ipBef = ip, ipAft = ++ip;
       if ( ipAft == _intPoints.end() )
         return false;
       --ipBef;
@@ -570,12 +722,64 @@ namespace
   }
   //================================================================================
   /*
+   * Returns parameters of a point in i-th plane
+   */
+  gp_XY GridPlanes::GetUV( const gp_Pnt& p, const gp_Pnt& origin )
+  {
+    gp_Vec v( origin, p );
+    return gp_XY( v.Dot( _uNorm ) * _factor,
+                  v.Dot( _vNorm ) * _factor );
+  }
+  //================================================================================
+  /*
+   * Adds face IDs
+   */
+  void B_IntersectPoint::Add( const vector< TGeomID >& fIDs,
+                              const SMDS_MeshNode*     n) const
+  {
+    if ( _faceIDs.empty() )
+      _faceIDs = fIDs;
+    else
+      for ( size_t i = 0; i < fIDs.size(); ++i )
+      {
+        vector< TGeomID >::iterator it =
+          std::find( _faceIDs.begin(), _faceIDs.end(), fIDs[i] );
+        if ( it == _faceIDs.end() )
+          _faceIDs.push_back( fIDs[i] );
+      }
+    if ( !_node )
+      _node = n;
+  }
+  //================================================================================
+  /*
+   * Returns \c true if \a other B_IntersectPoint holds the same face ID
+   */
+  bool B_IntersectPoint::HasCommonFace( const B_IntersectPoint * other ) const
+  {
+    if ( other )
+      for ( size_t i = 0; i < other->_faceIDs.size(); ++i )
+        if ( IsOnFace( other->_faceIDs[i] ) )
+          return true;
+    return false;
+  }
+  //================================================================================
+  /*
+   * Returns \c true if \a faceID in in this->_faceIDs
+   */
+  bool B_IntersectPoint::IsOnFace( int faceID ) const // returns true if faceID is found
+  {
+    vector< TGeomID >::const_iterator it =
+      std::find( _faceIDs.begin(), _faceIDs.end(), faceID );
+    return ( it != _faceIDs.end() );
+  }
+  //================================================================================
+  /*
    * Return an iterator on GridLine's in a given direction
    */
   LineIndexer Grid::GetLineIndexer(size_t iDir) const
   {
     const size_t indices[] = { 1,2,0, 0,2,1, 0,1,2 };
-    const string s[] = { "X", "Y", "Z" };
+    const string s      [] = { "X", "Y", "Z" };
     LineIndexer li( _coords[0].size(),  _coords[1].size(),    _coords[2].size(),
                     indices[iDir*3],    indices[iDir*3+1],    indices[iDir*3+2],
                     s[indices[iDir*3]], s[indices[iDir*3+1]], s[indices[iDir*3+2]]);
@@ -588,11 +792,21 @@ namespace
   void Grid::SetCoordinates(const vector<double>& xCoords,
                             const vector<double>& yCoords,
                             const vector<double>& zCoords,
+                            const double*         axesDirs,
                             const TopoDS_Shape&   shape)
   {
     _coords[0] = xCoords;
     _coords[1] = yCoords;
     _coords[2] = zCoords;
+    _axes[0].SetCoord( axesDirs[0],
+                       axesDirs[1],
+                       axesDirs[2]);
+    _axes[1].SetCoord( axesDirs[3],
+                       axesDirs[4],
+                       axesDirs[5]);
+    _axes[2].SetCoord( axesDirs[6],
+                       axesDirs[7],
+                       axesDirs[8]);
 
     // compute tolerance
     _minCellSize = Precision::Infinite();
@@ -649,7 +863,7 @@ namespace
     const size_t nbGridNodes = _coords[0].size() * _coords[1].size() * _coords[2].size();
     vector< bool > isNodeOut( nbGridNodes, false );
     _nodes.resize( nbGridNodes, 0 );
-    _isBndNode.resize( nbGridNodes, false );
+    _gridIntP.resize( nbGridNodes, NULL );
 
     for ( int iDir = 0; iDir < 3; ++iDir ) // loop on 3 line directions
     {
@@ -669,8 +883,8 @@ namespace
 
         GridLine& line = _lines[ iDir ][ li.LineIndex() ];
         line.RemoveExcessIntPoints( _tol );
-        multiset< IntersectionPoint >& intPnts = _lines[ iDir ][ li.LineIndex() ]._intPoints;
-        multiset< IntersectionPoint >::iterator ip = intPnts.begin();
+        multiset< F_IntersectPoint >& intPnts = _lines[ iDir ][ li.LineIndex() ]._intPoints;
+        multiset< F_IntersectPoint >::iterator ip = intPnts.begin();
 
         bool isOut = true;
         const double* nodeCoord = & coords[0], *coord0 = nodeCoord, *coordEnd = coord0 + coords.size();
@@ -706,14 +920,18 @@ namespace
           else
           {
             int nodeIndex = nIndex0 + nShift * ( nodeCoord-coord0 );
-            if ( ! _nodes[ nodeIndex ] )
+            if ( !_nodes[ nodeIndex ] )
             {
               li.SetIndexOnLine( nodeCoord-coord0 );
               double xyz[3] = { _coords[0][ li.I() ], _coords[1][ li.J() ], _coords[2][ li.K() ]};
-              _nodes[ nodeIndex ] = helper.AddNode( xyz[0], xyz[1], xyz[2] );
-              _isBndNode[ nodeIndex ] = true;
+              _nodes   [ nodeIndex ] = helper.AddNode( xyz[0], xyz[1], xyz[2] );
+              _gridIntP[ nodeIndex ] = & * ip;
             }
-            //ip->_node = _nodes[ nodeIndex ];
+            if ( _gridIntP[ nodeIndex ] )
+              _gridIntP[ nodeIndex ]->Add( ip->_faceIDs );
+            else
+              _gridIntP[ nodeIndex ] = & * ip;
+            // ip->_node        = _nodes[ nodeIndex ]; -- to differ from ip on links
             ip->_indexOnLine = nodeCoord-coord0;
             if ( ++nodeCoord < coordEnd )
               nodeParam = *nodeCoord - *coord0;
@@ -744,7 +962,7 @@ namespace
       LineIndexer li = GetLineIndexer( iDir );
       for ( ; li.More(); ++li )
       {
-        multiset< IntersectionPoint >& intPnts = _lines[ iDir ][ li.LineIndex() ]._intPoints;
+        multiset< F_IntersectPoint >& intPnts = _lines[ iDir ][ li.LineIndex() ]._intPoints;
         if ( intPnts.empty() ) continue;
         if ( intPnts.size() == 1 )
         {
@@ -944,7 +1162,7 @@ namespace
   {
     if ( !toClassify || UVIsOnFace() )
     {
-      IntersectionPoint p;
+      F_IntersectPoint p;
       p._paramOnLine = _w;
       p._transition  = _transition;
       _intPoints.push_back( p );
@@ -1291,16 +1509,18 @@ namespace
     _origNodeInd   = _grid->NodeIndex( i,j,k );
     for ( int iN = 0; iN < 8; ++iN )
     {
-      _hexNodes[iN]._node = _grid->_nodes[ _origNodeInd + _nodeShift[iN] ];
+      _hexNodes[iN]._node     = _grid->_nodes   [ _origNodeInd + _nodeShift[iN] ];
+      _hexNodes[iN]._intPoint = _grid->_gridIntP[ _origNodeInd + _nodeShift[iN] ];
       _nbCornerNodes += bool( _hexNodes[iN]._node );
-      _nbBndNodes    += _grid->_isBndNode[ _origNodeInd + _nodeShift[iN] ];
+      _nbBndNodes    += bool( _hexNodes[iN]._intPoint );
     }
 
     _sideLength[0] = _grid->_coords[0][i+1] - _grid->_coords[0][i];
     _sideLength[1] = _grid->_coords[1][j+1] - _grid->_coords[1][j];
     _sideLength[2] = _grid->_coords[2][k+1] - _grid->_coords[2][k];
 
-    if ( _nbCornerNodes < 8 && _nbIntNodes + _nbCornerNodes > 3)
+    if ( _nbIntNodes + _edgeIntPnts.size() > 0 &&
+         _nbIntNodes + _nbCornerNodes + _edgeIntPnts.size() > 3)
     {
       _Link split;
       // create sub-links (_splits) by splitting links with _intNodes
@@ -1309,21 +1529,120 @@ namespace
         _Link& link = _hexLinks[ iLink ];
         link._splits.clear();
         split._nodes[ 0 ] = link._nodes[0];
-        for ( size_t i = 0; i < link._intNodes.size(); ++ i )
+        bool isOut = ( ! link._nodes[0]->Node() );
+        //int iEnd = link._intNodes.size() - bool( link._nodes[1]->_intPoint );
+        for ( size_t i = 0; i < link._intNodes.size(); ++i )
         {
-          if ( split._nodes[ 0 ]->Node() )
+          if ( link._intNodes[i].Node() )
           {
-            split._nodes[ 1 ] = &link._intNodes[i];
-            link._splits.push_back( split );
+            if ( split._nodes[ 0 ]->Node() && !isOut )
+            {
+              split._nodes[ 1 ] = &link._intNodes[i];
+              link._splits.push_back( split );
+            }
+            split._nodes[ 0 ] = &link._intNodes[i];
+          }
+          switch ( link._intNodes[i].FaceIntPnt()->_transition ) {
+          case Trans_OUT: isOut = true; break;
+          case Trans_IN : isOut = false; break;
+          default:; // isOut remains the same
           }
-          split._nodes[ 0 ] = &link._intNodes[i];
         }
-        if ( link._nodes[ 1 ]->Node() && split._nodes[ 0 ]->Node() )
+        if ( link._nodes[ 1 ]->Node() && split._nodes[ 0 ]->Node() && !isOut )
         {
           split._nodes[ 1 ] = link._nodes[1];
           link._splits.push_back( split );
         }
       }
+
+      // Create _Node's at intersections with EDGEs.
+
+      const double tol2 = _grid->_tol * _grid->_tol;
+      int facets[3], nbFacets, subEntity;
+
+      for ( size_t iP = 0; iP < _edgeIntPnts.size(); ++iP )
+      {
+        nbFacets = getEntity( _edgeIntPnts[iP], facets, subEntity );
+        _Node* equalNode = 0;
+        switch( nbFacets ) {
+        case 1: // in a _Face
+        {
+          _Face& quad = _hexQuads[ facets[0] - SMESH_Block::ID_FirstF ];
+          equalNode = FindEqualNode( quad._edgeNodes, _edgeIntPnts[ iP ], tol2 );
+          if ( equalNode ) {
+            equalNode->Add( _edgeIntPnts[ iP ] );
+          }
+          else {
+            quad._edgeNodes.push_back( _Node( 0, _edgeIntPnts[ iP ]));
+            ++_nbIntNodes;
+          }
+          break;
+        }
+        case 2: // on a _Link
+        {
+          _Link& link = _hexLinks[ subEntity - SMESH_Block::ID_FirstE ];
+          if ( link._splits.size() > 0 )
+          {
+            equalNode = FindEqualNode( link._intNodes, _edgeIntPnts[ iP ], tol2 );
+            if ( equalNode )
+              equalNode->Add( _edgeIntPnts[ iP ] );
+          }
+          else
+          {
+            for ( int iF = 0; iF < 2; ++iF )
+            {
+              _Face& quad = _hexQuads[ facets[iF] - SMESH_Block::ID_FirstF ];
+              equalNode = FindEqualNode( quad._edgeNodes, _edgeIntPnts[ iP ], tol2 );
+              if ( equalNode ) {
+                equalNode->Add( _edgeIntPnts[ iP ] );
+              }
+              else {
+                quad._edgeNodes.push_back( _Node( 0, _edgeIntPnts[ iP ]));
+                ++_nbIntNodes;
+              }
+            }
+          }
+          break;
+        }
+        case 3: // at a corner
+        {
+          _Node& node = _hexNodes[ subEntity - SMESH_Block::ID_FirstV ];
+          if ( node.Node() > 0 )
+          {
+            if ( node._intPoint )
+              node._intPoint->Add( _edgeIntPnts[ iP ]->_faceIDs, _edgeIntPnts[ iP ]->_node );
+          }
+          else
+          {
+            for ( int iF = 0; iF < 3; ++iF )
+            {
+              _Face& quad = _hexQuads[ facets[iF] - SMESH_Block::ID_FirstF ];
+              equalNode = FindEqualNode( quad._edgeNodes, _edgeIntPnts[ iP ], tol2 );
+              if ( equalNode ) {
+                equalNode->Add( _edgeIntPnts[ iP ] );
+              }
+              else {
+                quad._edgeNodes.push_back( _Node( 0, _edgeIntPnts[ iP ]));
+                ++_nbIntNodes;
+              }
+            }
+          }
+          break;
+        }
+        default: // inside a hex
+        {
+          equalNode = FindEqualNode( _vertexNodes, _edgeIntPnts[ iP ], tol2 );
+          if ( equalNode ) {
+            equalNode->Add( _edgeIntPnts[ iP ] );
+          }
+          else {
+            _vertexNodes.push_back( _Node( 0, _edgeIntPnts[iP] ));
+            ++_nbIntNodes;
+          }
+        }
+        } // switch( nbFacets )
+
+      } // loop on _edgeIntPnts
     }
   }
   //================================================================================
@@ -1351,87 +1670,173 @@ namespace
     if ( _nbCornerNodes + _nbIntNodes < 4 )
       return;
 
-    if ( _nbBndNodes == _nbCornerNodes && isInHole() )
+    if ( _nbBndNodes == _nbCornerNodes && _nbIntNodes == 0 && isInHole() )
       return;
 
     _polygons.clear();
-
-    vector<const SMDS_MeshNode* > polyhedraNodes;
-    vector<int>                   quantities;
+    _polygons.reserve( 10 );
 
     // create polygons from quadrangles and get their nodes
 
-    vector<_Node*> nodes;
-    nodes.reserve( _nbCornerNodes + _nbIntNodes );
-
     _Link polyLink;
-    polyLink._faces.reserve( 1 );
+    vector< _OrientedLink > splits;
+    vector<_Node*> chainNodes;
+
+    bool hasEdgeIntersections = !_edgeIntPnts.empty();
 
     for ( int iF = 0; iF < 6; ++iF ) // loop on 6 sides of a hexahedron
     {
-      const _Face& quad = _hexQuads[ iF ] ;
+      _Face& quad = _hexQuads[ iF ] ;
 
       _polygons.resize( _polygons.size() + 1 );
-      _Face& polygon = _polygons.back();
-      polygon._links.clear();
-      polygon._polyLinks.clear(); polygon._polyLinks.reserve( 10 );
+      _Face* polygon = &_polygons.back();
+      polygon->_polyLinks.reserve( 20 );
 
-      // add splits of a link to a polygon and collect info on nodes
-      //int nbIn = 0, nbOut = 0, nbCorners = 0;
-      nodes.clear();
+      splits.clear();
       for ( int iE = 0; iE < 4; ++iE ) // loop on 4 sides of a quadrangle
+        for ( int iS = 0; iS < quad._links[ iE ].NbResultLinks(); ++iS )
+          splits.push_back( quad._links[ iE ].ResultLink( iS ));
+
+      // add splits of links to a polygon and add _polyLinks to make
+      // polygon's boundary closed
+
+      int nbSplits = splits.size();
+      if ( nbSplits < 2 && quad._edgeNodes.empty() )
+        nbSplits = 0;
+
+      if ( nbSplits == 0 && !quad._edgeNodes.empty() )
+      {
+        // make _vertexNodes from _edgeNodes of an empty quad
+        const double tol2 = _grid->_tol * _grid->_tol;
+        for ( size_t iP = 0; iP < quad._edgeNodes.size(); ++iP )
+        {
+          _Node* equalNode =
+            FindEqualNode( _vertexNodes, quad._edgeNodes[ iP ].EdgeIntPnt(), tol2 );
+          if ( equalNode )
+            equalNode->Add( quad._edgeNodes[ iP ].EdgeIntPnt() );
+          else
+            _vertexNodes.push_back( quad._edgeNodes[ iP ]);
+        }
+      }
+
+      while ( nbSplits > 0 )
       {
-        int nbSpits = quad._links[ iE ].NbResultLinks();
-        for ( int iS = 0; iS < nbSpits; ++iS )
+        size_t iS = 0;
+        while ( !splits[ iS ] )
+          ++iS;
+
+        if ( !polygon->_links.empty() )
+        {
+          _polygons.resize( _polygons.size() + 1 );
+          polygon = &_polygons.back();
+          polygon->_polyLinks.reserve( 20 );
+        }
+        polygon->_links.push_back( splits[ iS ] );
+        splits[ iS++ ]._link = 0;
+        --nbSplits;
+
+        _Node* nFirst = polygon->_links.back().FirstNode();
+        _Node *n1,*n2 = polygon->_links.back().LastNode();
+        for ( ; nFirst != n2 && iS < splits.size(); ++iS )
         {
-          _OrientedLink split = quad._links[ iE ].ResultLink( iS );
-          _Node* n = split.FirstNode();
-          if ( !polygon._links.empty() )
+          _OrientedLink& split = splits[ iS ];
+          if ( !split ) continue;
+
+          n1 = split.FirstNode();
+          if ( n1 != n2 )
           {
-            _Node* nPrev = polygon._links.back().LastNode();
-            if ( nPrev != n )
+            // try to connect to intersections with EDGES
+            if ( quad._edgeNodes.size() > 0  &&
+                 findChain( n2, n1, quad, chainNodes ))
             {
-              polyLink._nodes[0] = nPrev;
-              polyLink._nodes[1] = n;
-              polygon._polyLinks.push_back( polyLink );
-              polygon._links.push_back( _OrientedLink( &polygon._polyLinks.back() ));
-              nodes.push_back( nPrev );
+              for ( size_t i = 1; i < chainNodes.size(); ++i )
+              {
+                polyLink._nodes[0] = chainNodes[i-1];
+                polyLink._nodes[1] = chainNodes[i];
+                polygon->_polyLinks.push_back( polyLink );
+                polygon->_links.push_back( _OrientedLink( &polygon->_polyLinks.back() ));
+              }
             }
+            // try to connect to a split ending on the same FACE
+            else
+            {
+              _OrientedLink foundSplit;
+              for ( int i = iS; i < splits.size() && !foundSplit; ++i )
+                if (( foundSplit = splits[ i ]) &&
+                    ( n2->IsLinked( foundSplit.FirstNode()->_intPoint )))
+                {
+                  polyLink._nodes[0] = n2;
+                  polyLink._nodes[1] = foundSplit.FirstNode();
+                  polygon->_polyLinks.push_back( polyLink );
+                  polygon->_links.push_back( _OrientedLink( &polygon->_polyLinks.back() ));
+                  iS = i - 1;
+                }
+                else
+                {
+                  foundSplit._link = 0;
+                }
+              if ( foundSplit )
+              {
+                n2 = foundSplit.FirstNode();
+                continue;
+              }
+              else
+              {
+                if ( n2->IsLinked( nFirst->_intPoint ))
+                  break;
+                polyLink._nodes[0] = n2;
+                polyLink._nodes[1] = n1;
+                polygon->_polyLinks.push_back( polyLink );
+                polygon->_links.push_back( _OrientedLink( &polygon->_polyLinks.back() ));
+              }
+            }
+          }
+          polygon->_links.push_back( split );
+          split._link = 0;
+          --nbSplits;
+          n2 = polygon->_links.back().LastNode();
+
+        } // loop on splits
+
+        if ( nFirst != n2 ) // close a polygon
+        {
+          findChain( n2, nFirst, quad, chainNodes );
+          for ( size_t i = 1; i < chainNodes.size(); ++i )
+          {
+            polyLink._nodes[0] = chainNodes[i-1];
+            polyLink._nodes[1] = chainNodes[i];
+            polygon->_polyLinks.push_back( polyLink );
+            polygon->_links.push_back( _OrientedLink( &polygon->_polyLinks.back() ));
           }
-          polygon._links.push_back( split );
-          nodes.push_back( n );
         }
-      }
-      if ( polygon._links.size() > 1 )
-      {
-        _Node* n1 = polygon._links.back().LastNode();
-        _Node* n2 = polygon._links.front().FirstNode();
-        if ( n1 != n2 )
+
+        if ( polygon->_links.size() < 3 && nbSplits > 0 )
         {
-          polyLink._nodes[0] = n1;
-          polyLink._nodes[1] = n2;
-          polygon._polyLinks.push_back( polyLink );
-          polygon._links.push_back( _OrientedLink( &polygon._polyLinks.back() ));
-          nodes.push_back( n1 );
+          polygon->_polyLinks.clear();
+          polygon->_links.clear();
         }
-        // add polygon to its links
-        for ( size_t iL = 0; iL < polygon._links.size(); ++iL )
-          polygon._links[ iL ]._link->_faces.push_back( &polygon );
-        // store polygon nodes
-        quantities.push_back( nodes.size() );
-        for ( size_t i = 0; i < nodes.size(); ++i )
-          polyhedraNodes.push_back( nodes[i]->Node() );
-      }
-      else
-      {
-        _polygons.resize( _polygons.size() - 1 );
-      }
-    }
+      } // while ( nbSplits > 0 )
+
+      if ( polygon->_links.size() < 3 )
+        _polygons.pop_back();
+
+    }  // loop on 6 sides of a hexahedron
 
     // create polygons closing holes in a polyhedron
 
+    // add polygons to their links
+    for ( size_t iP = 0; iP < _polygons.size(); ++iP )
+    {
+      _Face& polygon = _polygons[ iP ];
+      for ( size_t iL = 0; iL < polygon._links.size(); ++iL )
+      {
+        polygon._links[ iL ]._link->_faces.reserve( 2 );
+        polygon._links[ iL ]._link->_faces.push_back( &polygon );
+      }
+    }
     // find free links
     vector< _OrientedLink* > freeLinks;
+    freeLinks.reserve(20);
     for ( size_t iP = 0; iP < _polygons.size(); ++iP )
     {
       _Face& polygon = _polygons[ iP ];
@@ -1439,58 +1844,157 @@ namespace
         if ( polygon._links[ iL ]._link->_faces.size() < 2 )
           freeLinks.push_back( & polygon._links[ iL ]);
     }
-    // make closed chains of free links
     int nbFreeLinks = freeLinks.size();
     if ( 0 < nbFreeLinks && nbFreeLinks < 3 ) return;
+
+    set<TGeomID> usedFaceIDs;
+
+    // make closed chains of free links
     while ( nbFreeLinks > 0 )
     {
-      nodes.clear();
       _polygons.resize( _polygons.size() + 1 );
       _Face& polygon = _polygons.back();
-      polygon._links.clear();
+      polygon._polyLinks.reserve( 20 );
+      polygon._links.reserve( 20 );
 
-      // get a remaining link to start from
       _OrientedLink* curLink = 0;
-      for ( size_t iL = 0; iL < freeLinks.size() && !curLink; ++iL )
-        if (( curLink = freeLinks[ iL ] ))
-          freeLinks[ iL ] = 0;
-      nodes.push_back( curLink->LastNode() );
-      polygon._links.push_back( *curLink );
-
-      // find all links connected to curLink
-      _Node* curNode = 0;
-      do
+      _Node*         curNode;
+      if ( !hasEdgeIntersections )
       {
-        curNode = curLink->FirstNode();
-        curLink = 0;
+        // get a remaining link to start from
         for ( size_t iL = 0; iL < freeLinks.size() && !curLink; ++iL )
-          if ( freeLinks[ iL ] && freeLinks[ iL ]->LastNode() == curNode )
-          {
-            curLink = freeLinks[ iL ];
+          if (( curLink = freeLinks[ iL ] ))
             freeLinks[ iL ] = 0;
-            nodes.push_back( curNode );
-            polygon._links.push_back( *curLink );
+        polygon._links.push_back( *curLink );
+        --nbFreeLinks;
+        do
+        {
+          // find all links connected to curLink
+          curNode = curLink->FirstNode();
+          curLink = 0;
+          for ( size_t iL = 0; iL < freeLinks.size() && !curLink; ++iL )
+            if ( freeLinks[ iL ] && freeLinks[ iL ]->LastNode() == curNode )
+            {
+              curLink = freeLinks[ iL ];
+              freeLinks[ iL ] = 0;
+              polygon._links.push_back( *curLink );
+              --nbFreeLinks;
+            }
+        } while ( curLink );
+      }
+      else // there are intersections with EDGEs
+      {
+        TGeomID curFace;
+        // get a remaining link to start from, one lying on minimal
+        // nb of FACEs
+        {
+          map< vector< TGeomID >, int > facesOfLink;
+          map< vector< TGeomID >, int >::iterator f2l;
+          for ( size_t iL = 0; iL < freeLinks.size(); ++iL )
+            if ( freeLinks[ iL ] )
+            {
+              f2l = facesOfLink.insert
+                ( make_pair( freeLinks[ iL ]->GetNotUsedFace( usedFaceIDs ), iL )).first;
+              if ( f2l->first.size() == 1 )
+                break;
+            }
+          f2l = facesOfLink.begin();
+          if ( f2l->first.empty() )
+            return;
+          curFace = f2l->first[0];
+          curLink = freeLinks[ f2l->second ];
+          freeLinks[ f2l->second ] = 0;
+        }
+        usedFaceIDs.insert( curFace );
+        polygon._links.push_back( *curLink );
+        --nbFreeLinks;
+
+        // find all links bounding a FACE of curLink
+        do
+        {
+          // go forward from curLink
+          curNode = curLink->LastNode();
+          curLink = 0;
+          for ( size_t iL = 0; iL < freeLinks.size() && !curLink; ++iL )
+            if ( freeLinks[ iL ] &&
+                 freeLinks[ iL ]->FirstNode() == curNode &&
+                 freeLinks[ iL ]->LastNode()->IsOnFace( curFace ))
+            {
+              curLink = freeLinks[ iL ];
+              freeLinks[ iL ] = 0;
+              polygon._links.push_back( *curLink );
+              --nbFreeLinks;
+            }
+        } while ( curLink );
+
+        std::reverse( polygon._links.begin(), polygon._links.end() );
+
+        curLink = & polygon._links.back();
+        do
+        {
+          // go backward from curLink
+          curNode = curLink->FirstNode();
+          curLink = 0;
+          for ( size_t iL = 0; iL < freeLinks.size() && !curLink; ++iL )
+            if ( freeLinks[ iL ] &&
+                 freeLinks[ iL ]->LastNode() == curNode &&
+                 freeLinks[ iL ]->FirstNode()->IsOnFace( curFace ))
+            {
+              curLink = freeLinks[ iL ];
+              freeLinks[ iL ] = 0;
+              polygon._links.push_back( *curLink );
+              --nbFreeLinks;
+            }
+        } while ( curLink );
+
+        curNode = polygon._links.back().FirstNode();
+
+        if ( polygon._links[0].LastNode() != curNode )
+        {
+          if ( !_vertexNodes.empty() )
+          {
+            // add links with _vertexNodes if not already used
+            for ( size_t iN = 0; iN < _vertexNodes.size(); ++iN )
+              if ( _vertexNodes[ iN ].IsOnFace( curFace ))
+              {
+                bool used = ( curNode == &_vertexNodes[ iN ] );
+                for ( size_t iL = 0; iL < polygon._links.size() && !used; ++iL )
+                  used = ( &_vertexNodes[ iN ] ==  polygon._links[ iL ].LastNode() );
+                if ( !used )
+                {
+                  polyLink._nodes[0] = &_vertexNodes[ iN ];
+                  polyLink._nodes[1] = curNode;
+                  polygon._polyLinks.push_back( polyLink );
+                  polygon._links.push_back( _OrientedLink( &polygon._polyLinks.back() ));
+                  freeLinks.push_back( &polygon._links.back() );
+                  ++nbFreeLinks;
+                  curNode = &_vertexNodes[ iN ];
+                }
+                // TODO: to reorder _vertexNodes within polygon, if there are several ones
+              }
           }
-      } while ( curLink );
+          polyLink._nodes[0] = polygon._links[0].LastNode();
+          polyLink._nodes[1] = curNode;
+          polygon._polyLinks.push_back( polyLink );
+          polygon._links.push_back( _OrientedLink( &polygon._polyLinks.back() ));
+          freeLinks.push_back( &polygon._links.back() );
+          ++nbFreeLinks;
+        }
 
-      nbFreeLinks -= polygon._links.size();
+      } // if there are intersections with EDGEs
 
-      if ( curNode != nodes.front() || polygon._links.size() < 3 )
+      if ( polygon._links.size() < 3 ||
+           polygon._links[0].LastNode() != polygon._links.back().FirstNode() )
         return; // closed polygon not found -> invalid polyhedron
 
-      quantities.push_back( nodes.size() );
-      for ( size_t i = 0; i < nodes.size(); ++i )
-        polyhedraNodes.push_back( nodes[i]->Node() );
-
-      // add polygon to its links and reverse links
-      for ( size_t i = 0; i < polygon._links.size(); ++i )
+      // add polygon to its links
+      for ( size_t iL = 0; iL < polygon._links.size(); ++iL )
       {
-        polygon._links[i].Reverse();
-        polygon._links[i]._link->_faces.push_back( &polygon );
+        polygon._links[ iL ]._link->_faces.reserve( 2 );
+        polygon._links[ iL ]._link->_faces.push_back( &polygon );
+        polygon._links[ iL ].Reverse();
       }
-
-      //const size_t firstPoly = _polygons.size();
-    }
+    } // while ( nbFreeLinks > 0 )
 
     if ( ! checkPolyhedronSize() )
     {
@@ -1505,20 +2009,32 @@ namespace
     else if ( nbNodes == 6 && _polygons.size() == 5 ) isClassicElem = addPenta();
     else if ( nbNodes == 5 && _polygons.size() == 5 ) isClassicElem = addPyra ();
     if ( !isClassicElem )
-      _volumeDefs.set( polyhedraNodes, quantities );
+    {
+      _volumeDefs._nodes.clear();
+      _volumeDefs._quantities.clear();
+
+      for ( size_t iF = 0; iF < _polygons.size(); ++iF )
+      {
+        const size_t nbLinks = _polygons[ iF ]._links.size();
+        _volumeDefs._quantities.push_back( nbLinks );
+        for ( size_t iL = 0; iL < nbLinks; ++iL )
+          _volumeDefs._nodes.push_back( _polygons[ iF ]._links[ iL ].FirstNode() );
+      }
+    }
   }
   //================================================================================
   /*!
    * \brief Create elements in the mesh
    */
-  int Hexahedron::MakeElements(SMESH_MesherHelper& helper)
+  int Hexahedron::MakeElements(SMESH_MesherHelper&                      helper,
+                               const map< TGeomID, vector< TGeomID > >& edge2faceIDsMap)
   {
     SMESHDS_Mesh* mesh = helper.GetMeshDS();
 
     size_t nbCells[3] = { _grid->_coords[0].size() - 1,
                           _grid->_coords[1].size() - 1,
                           _grid->_coords[2].size() - 1 };
-    const size_t nbGridCells = nbCells[0] *nbCells [1] * nbCells[2];
+    const size_t nbGridCells = nbCells[0] * nbCells[1] * nbCells[2];
     vector< Hexahedron* > intersectedHex( nbGridCells, 0 );
     int nbIntHex = 0;
 
@@ -1535,10 +2051,10 @@ namespace
       for ( ; lineInd.More(); ++lineInd )
       {
         GridLine& line = _grid->_lines[ iDir ][ lineInd.LineIndex() ];
-        multiset< IntersectionPoint >::const_iterator ip = line._intPoints.begin();
+        multiset< F_IntersectPoint >::const_iterator ip = line._intPoints.begin();
         for ( ; ip != line._intPoints.end(); ++ip )
         {
-          if ( !ip->_node ) continue;
+          //if ( !ip->_node ) continue;
           lineInd.SetIndexOnLine( ip->_indexOnLine );
           for ( int iL = 0; iL < 4; ++iL ) // loop on 4 cells sharing a link
           {
@@ -1561,12 +2077,15 @@ namespace
             }
             const int iLink = iL + iDir * 4;
             hex->_hexLinks[iLink]._intNodes.push_back( _Node( 0, &(*ip) ));
-            hex->_nbIntNodes++;
+            hex->_nbIntNodes += bool( ip->_node );
           }
         }
       }
     }
 
+    // implement geom edges into the mesh
+    addEdges( helper, intersectedHex, edge2faceIDsMap );
+
     // add not split hexadrons to the mesh
     int nbAdded = 0;
     vector<int> intHexInd( nbIntHex );
@@ -1638,6 +2157,366 @@ namespace
 
   //================================================================================
   /*!
+   * \brief Implements geom edges into the mesh
+   */
+  void Hexahedron::addEdges(SMESH_MesherHelper&                      helper,
+                            vector< Hexahedron* >&                   hexes,
+                            const map< TGeomID, vector< TGeomID > >& edge2faceIDsMap)
+  {
+    if ( edge2faceIDsMap.empty() ) return;
+
+    // Prepare planes for intersecting with EDGEs
+    GridPlanes pln[3];
+    {
+      gp_XYZ origPnt = ( _grid->_coords[0][0] * _grid->_axes[0] +
+                         _grid->_coords[1][0] * _grid->_axes[1] +
+                         _grid->_coords[2][0] * _grid->_axes[2] );
+      for ( int iDirZ = 0; iDirZ < 3; ++iDirZ ) // iDirZ gives normal direction to planes
+      {
+        GridPlanes& planes = pln[ iDirZ ];
+        int iDirX = ( iDirZ + 1 ) % 3;
+        int iDirY = ( iDirZ + 2 ) % 3;
+        planes._uNorm  = ( _grid->_axes[ iDirY ] ^ _grid->_axes[ iDirZ ] ).Normalized();
+        planes._vNorm  = ( _grid->_axes[ iDirZ ] ^ _grid->_axes[ iDirX ] ).Normalized();
+        planes._zNorm  = ( _grid->_axes[ iDirX ] ^ _grid->_axes[ iDirY ] ).Normalized();
+        double   uvDot = planes._uNorm * planes._vNorm;
+        planes._factor = sqrt( 1. - uvDot * uvDot );
+        planes._origins.resize( _grid->_coords[ iDirZ ].size() );
+        planes._zProjs.resize ( _grid->_coords[ iDirZ ].size() );
+        planes._origins[0] = origPnt;
+        planes._zProjs [0] = 0;
+        const double       zFactor = _grid->_axes[ iDirZ ] * planes._zNorm;
+        const vector< double > & u = _grid->_coords[ iDirZ ];
+        for ( int i = 1; i < planes._origins.size(); ++i )
+        {
+          planes._origins[i] = origPnt + _grid->_axes[ iDirZ ] * ( u[i] - u[0] );
+          planes._zProjs [i] = zFactor * ( u[i] - u[0] );
+        }
+      }
+    }
+    const double deflection = _grid->_minCellSize / 20.;
+    const double tol        = _grid->_tol;
+    // int facets[6] = { SMESH_Block::ID_F0yz, SMESH_Block::ID_F1yz,
+    //                   SMESH_Block::ID_Fx0z, SMESH_Block::ID_Fx1z,
+    //                   SMESH_Block::ID_Fxy0, SMESH_Block::ID_Fxy1 };
+    E_IntersectPoint ip;
+    //ip._faceIDs.reserve(2);
+
+    // Intersect EDGEs with the planes
+    map< TGeomID, vector< TGeomID > >::const_iterator e2fIt = edge2faceIDsMap.begin();
+    for ( ; e2fIt != edge2faceIDsMap.end(); ++e2fIt )
+    {
+      const TGeomID  edgeID = e2fIt->first;
+      const TopoDS_Edge & E = TopoDS::Edge( _grid->_shapes( edgeID ));
+      BRepAdaptor_Curve curve( E );
+
+      ip._faceIDs = e2fIt->second;
+      ip._shapeID = edgeID;
+
+      // discretize the EGDE
+      GCPnts_UniformDeflection discret( curve, deflection, true );
+      if ( !discret.IsDone() || discret.NbPoints() < 2 )
+        continue;
+
+      // perform intersection
+      for ( int iDirZ = 0; iDirZ < 3; ++iDirZ )
+      {
+        GridPlanes& planes = pln[ iDirZ ];
+        int      iDirX = ( iDirZ + 1 ) % 3;
+        int      iDirY = ( iDirZ + 2 ) % 3;
+        double    xLen = _grid->_coords[ iDirX ].back() - _grid->_coords[ iDirX ][0];
+        double    yLen = _grid->_coords[ iDirY ].back() - _grid->_coords[ iDirY ][0];
+        double zFactor = _grid->_axes[ iDirZ ] * planes._zNorm;
+        int dIJK[3], d000[3] = { 0,0,0 };
+
+        // locate the 1st point of a segment within the grid
+        gp_XYZ p1     = discret.Value( 1 ).XYZ();
+        double u1     = discret.Parameter( 1 );
+        double zProj1 = planes._zNorm * ( p1 - planes._origins[0] );
+        gp_Pnt orig   = planes._origins[0] + planes._zNorm * zProj1;
+        gp_XY uv      = planes.GetUV( p1, orig );
+        int iX1       = int( uv.X() / xLen * ( _grid->_coords[ iDirX ].size() - 1. ));
+        int iY1       = int( uv.Y() / yLen * ( _grid->_coords[ iDirY ].size() - 1. ));
+        int iZ1       = int( zProj1 / planes._zProjs.back() * ( planes._zProjs.size() - 1. ));
+        locateValue( iX1, uv.X(), _grid->_coords[ iDirX ], dIJK[ iDirX ], tol );
+        locateValue( iY1, uv.Y(), _grid->_coords[ iDirY ], dIJK[ iDirY ], tol );
+        locateValue( iZ1, zProj1, planes._zProjs         , dIJK[ iDirZ ], tol );
+
+        int ijk[3]; // grid index where a segment intersect a plane
+        ijk[ iDirX ] = iX1;
+        ijk[ iDirY ] = iY1;
+        ijk[ iDirZ ] = iZ1;
+        ip._uvw[ iDirX ] = uv.X()           + _grid->_coords[ iDirX ][0];
+        ip._uvw[ iDirY ] = uv.Y()           + _grid->_coords[ iDirY ][0];
+        ip._uvw[ iDirZ ] = zProj1 / zFactor + _grid->_coords[ iDirZ ][0];
+
+        // add the 1st vertex point to a hexahedron
+        if ( iDirZ == 0 )
+        {
+          //ip._shapeID = _grid->_shapes.Add( helper.IthVertex( 0, curve.Edge(),/*CumOri=*/false));
+          ip._point   = p1;
+          _grid->_edgeIntP.push_back( ip );
+          if ( !addIntersection( _grid->_edgeIntP.back(), hexes, ijk, d000 ))
+            _grid->_edgeIntP.pop_back();
+        }
+        for ( int iP = 2; iP <= discret.NbPoints(); ++iP )
+        {
+          // locate the 2nd point of a segment within the grid
+          gp_XYZ p2     = discret.Value( iP ).XYZ();
+          double u2     = discret.Parameter( iP );
+          double zProj2 = planes._zNorm * ( p2 - planes._origins[0] );
+          int iZ2       = iZ1;
+          locateValue( iZ2, zProj2, planes._zProjs, dIJK[ iDirZ ], tol );
+
+          // treat intersections with planes between 2 end points of a segment
+          int dZ = ( iZ1 <= iZ2 ) ? +1 : -1;
+          int iZ = iZ1 + ( iZ1 < iZ2 );
+          for ( int i = 0, nb = Abs( iZ1 - iZ2 ); i < nb; ++i, iZ += dZ )
+          {
+            ip._point = findIntPoint( u1, zProj1, u2, zProj2,
+                                      planes._zProjs[ iZ ],
+                                      curve, planes._zNorm, planes._origins[0] );
+            gp_XY uv  = planes.GetUV( ip._point, planes._origins[ iZ ]);
+            locateValue( ijk[ iDirX ], uv.X(), _grid->_coords[ iDirX ], dIJK[ iDirX ], tol );
+            locateValue( ijk[ iDirY ], uv.Y(), _grid->_coords[ iDirY ], dIJK[ iDirY ], tol );
+            ijk[ iDirZ ] = iZ;
+            ip._uvw[ iDirX ] = uv.X()                         + _grid->_coords[ iDirX ][0];
+            ip._uvw[ iDirY ] = uv.Y()                         + _grid->_coords[ iDirY ][0];
+            ip._uvw[ iDirZ ] = planes._zProjs[ iZ ] / zFactor + _grid->_coords[ iDirZ ][0];
+
+            // add ip to hex "above" the plane
+            _grid->_edgeIntP.push_back( ip );
+            dIJK[ iDirZ ] = 0;
+            bool added = addIntersection(_grid->_edgeIntP.back(), hexes, ijk, dIJK);
+
+            // add ip to hex "below" the plane
+            ijk[ iDirZ ] = iZ-1;
+            if ( !addIntersection( _grid->_edgeIntP.back(), hexes, ijk, dIJK ) &&
+                 !added)
+              _grid->_edgeIntP.pop_back();
+          }
+          iZ1    = iZ2;
+          p1     = p2;
+          u1     = u2;
+          zProj1 = zProj2;
+        }
+        // add the 2nd vertex point to a hexahedron
+        if ( iDirZ == 0 )
+        {
+          orig = planes._origins[0] + planes._zNorm * zProj1;
+          uv   = planes.GetUV( p1, orig );
+          locateValue( ijk[ iDirX ], uv.X(), _grid->_coords[ iDirX ], dIJK[ iDirX ], tol );
+          locateValue( ijk[ iDirY ], uv.Y(), _grid->_coords[ iDirY ], dIJK[ iDirY ], tol );
+          ijk[ iDirZ ] = iZ1;
+          ip._uvw[ iDirX ] = uv.X()           + _grid->_coords[ iDirX ][0];
+          ip._uvw[ iDirY ] = uv.Y()           + _grid->_coords[ iDirY ][0];
+          ip._uvw[ iDirZ ] = zProj1 / zFactor + _grid->_coords[ iDirZ ][0];
+          ip._point = p1;
+          _grid->_edgeIntP.push_back( ip );
+          if ( !addIntersection( _grid->_edgeIntP.back(), hexes, ijk, d000 ))
+            _grid->_edgeIntP.pop_back();
+        }
+      } // loop on 3 grid directions
+    } // loop on EDGEs
+
+    // Create nodes at found intersections
+    // const E_IntersectPoint* eip;
+    // for ( size_t i = 0; i < hexes.size(); ++i )
+    // {
+    //   Hexahedron* h = hexes[i];
+    //   if ( !h ) continue;
+    //   for ( int iF = 0; iF < 6; ++iF )
+    //   {
+    //     _Face& quad = h->_hexQuads[ iF ];
+    //     for ( size_t iP = 0; iP < quad._edgeNodes.size(); ++iP )
+    //       if ( !quad._edgeNodes[ iP ]._node )
+    //         if (( eip = quad._edgeNodes[ iP ].EdgeIntPnt() ))
+    //           quad._edgeNodes[ iP ]._intPoint->_node = helper.AddNode( eip->_point.X(),
+    //                                                                    eip->_point.Y(),
+    //                                                                    eip->_point.Z() );
+    //   }
+    //   for ( size_t iP = 0; iP < hexes[i]->_vertexNodes.size(); ++iP )
+    //     if (( eip = h->_vertexNodes[ iP ].EdgeIntPnt() ))
+    //       h->_vertexNodes[ iP ]._intPoint->_node = helper.AddNode( eip->_point.X(),
+    //                                                                eip->_point.Y(),
+    //                                                                eip->_point.Z() );
+    // }
+  }
+
+  //================================================================================
+  /*!
+   * \brief Finds intersection of a curve with a plane
+   *  \param [in] u1 - parameter of one curve point
+   *  \param [in] proj1 - projection of the curve point to the plane normal
+   *  \param [in] u2 - parameter of another curve point
+   *  \param [in] proj2 - projection of the other curve point to the plane normal
+   *  \param [in] proj - projection of a point where the curve intersects the plane
+   *  \param [in] curve - the curve
+   *  \param [in] axis - the plane normal
+   *  \param [in] origin - the plane origin
+   *  \return gp_Pnt - the found intersection point
+   */
+  //================================================================================
+
+  gp_Pnt Hexahedron::findIntPoint( double u1, double proj1,
+                                   double u2, double proj2,
+                                   double proj,
+                                   BRepAdaptor_Curve& curve,
+                                   const gp_XYZ& axis,
+                                   const gp_XYZ& origin)
+  {
+    double r = (( proj - proj1 ) / ( proj2 - proj1 ));
+    double u = u1 * ( 1 - r ) + u2 * r;
+    gp_Pnt p = curve.Value( u );
+    double newProj =  axis * ( p.XYZ() - origin );
+    if ( Abs( proj - newProj ) > _grid->_tol / 10. )
+    {
+      if ( r > 0.5 )
+        return findIntPoint( u2, proj2, u, newProj, proj, curve, axis, origin );
+      else
+        return findIntPoint( u1, proj2, u, newProj, proj, curve, axis, origin );
+    }
+    return p;
+  }
+
+  //================================================================================
+  /*!
+   * \brief Returns index of a hexahedron sub-entities holding a point
+   *  \param [in] ip - intersection point
+   *  \param [out] facets - 0-3 facets holding a point
+   *  \param [out] sub - index of a vertex or an edge holding a point
+   *  \return int - number of facets holding a point
+   */
+  int Hexahedron::getEntity( const E_IntersectPoint* ip, int* facets, int& sub )
+  {
+    enum { X = 1, Y = 2, Z = 4 }; // == 001, 010, 100
+    int nbFacets = 0;
+    int vertex = 0, egdeMask = 0;
+
+    if ( Abs( _grid->_coords[0][ _i   ] - ip->_uvw[0] ) < _grid->_tol ) {
+      facets[ nbFacets++ ] = SMESH_Block::ID_F0yz;
+      egdeMask |= X;
+    }
+    else if ( Abs( _grid->_coords[0][ _i+1 ] - ip->_uvw[0] ) < _grid->_tol ) {
+      facets[ nbFacets++ ] = SMESH_Block::ID_F1yz;
+      vertex   |= X;
+      egdeMask |= X;
+    }
+    if ( Abs( _grid->_coords[1][ _j   ] - ip->_uvw[1] ) < _grid->_tol ) {
+      facets[ nbFacets++ ] = SMESH_Block::ID_Fx0z;
+      egdeMask |= Y;
+    }
+    else if ( Abs( _grid->_coords[1][ _j+1 ] - ip->_uvw[1] ) < _grid->_tol ) {
+      facets[ nbFacets++ ] = SMESH_Block::ID_Fx1z;
+      vertex   |= Y;
+      egdeMask |= Y;
+    }
+    if ( Abs( _grid->_coords[2][ _k   ] - ip->_uvw[2] ) < _grid->_tol ) {
+      facets[ nbFacets++ ] = SMESH_Block::ID_Fxy0;
+      egdeMask |= Z;
+    }
+    else if ( Abs( _grid->_coords[2][ _k+1 ] - ip->_uvw[2] ) < _grid->_tol ) {
+      facets[ nbFacets++ ] = SMESH_Block::ID_Fxy1;
+      vertex   |= Z;
+      egdeMask |= Z;
+    }
+
+    switch ( nbFacets )
+    {
+    case 0: sub = 0;         break;
+    case 1: sub = facets[0]; break;
+    case 2: {
+      const int edge [3][8] = {
+        { SMESH_Block::ID_E00z, SMESH_Block::ID_E10z,
+          SMESH_Block::ID_E01z, SMESH_Block::ID_E11z },
+        { SMESH_Block::ID_E0y0, SMESH_Block::ID_E1y0, 0, 0,
+          SMESH_Block::ID_E0y1, SMESH_Block::ID_E1y1 },
+        { SMESH_Block::ID_Ex00, 0, SMESH_Block::ID_Ex10, 0,
+          SMESH_Block::ID_Ex01, 0, SMESH_Block::ID_Ex11 }
+      };
+      switch ( egdeMask ) {
+      case X | Y: sub = edge[ 0 ][ vertex ]; break;
+      case X | Z: sub = edge[ 1 ][ vertex ]; break;
+      default:    sub = edge[ 2 ][ vertex ];
+      }
+      break;
+    }
+    //case 3:
+    default:
+      sub = vertex + SMESH_Block::ID_FirstV;
+    }
+
+    return nbFacets;
+  }
+  //================================================================================
+  /*!
+   * \brief Adds intersection with an EDGE
+   */
+  bool Hexahedron::addIntersection( const E_IntersectPoint& ip,
+                                    vector< Hexahedron* >&  hexes,
+                                    int ijk[], int dIJK[] )
+  {
+    bool added = false;
+
+    size_t hexIndex[4] = {
+      _grid->CellIndex( ijk[0], ijk[1], ijk[2] ),
+      dIJK[0] ? _grid->CellIndex( ijk[0]+dIJK[0], ijk[1], ijk[2] ) : -1,
+      dIJK[1] ? _grid->CellIndex( ijk[0], ijk[1]+dIJK[1], ijk[2] ) : -1,
+      dIJK[2] ? _grid->CellIndex( ijk[0], ijk[1], ijk[2]+dIJK[2] ) : -1
+    };
+    for ( int i = 0; i < 4; ++i )
+    {
+      if ( 0 <= hexIndex[i] && hexIndex[i] < hexes.size() && hexes[ hexIndex[i] ] )
+      {
+        Hexahedron* h = hexes[ hexIndex[i] ];
+        // check if ip is really inside the hex
+#ifdef _DEBUG_
+        if (( _grid->_coords[0][ h->_i   ] - _grid->_tol > ip._uvw[0] ) ||
+            ( _grid->_coords[0][ h->_i+1 ] + _grid->_tol < ip._uvw[0] ) ||
+            ( _grid->_coords[1][ h->_j   ] - _grid->_tol > ip._uvw[1] ) ||
+            ( _grid->_coords[1][ h->_j+1 ] + _grid->_tol < ip._uvw[1] ) ||
+            ( _grid->_coords[2][ h->_k   ] - _grid->_tol > ip._uvw[2] ) ||
+            ( _grid->_coords[2][ h->_k+1 ] + _grid->_tol < ip._uvw[2] ))
+          throw SALOME_Exception("ip outside a hex");
+#endif
+        h->_edgeIntPnts.push_back( & ip );
+        added = true;
+      }
+    }
+    return added;
+  }
+  //================================================================================
+  /*!
+   * \brief Finds nodes at a path from one node to another via intersections with EDGEs
+   */
+  bool Hexahedron::findChain( _Node*          n1,
+                              _Node*          n2,
+                              _Face&          quad,
+                              vector<_Node*>& chn )
+  {
+    chn.clear();
+    chn.push_back( n1 );
+    bool found = false;
+    do
+    {
+      found = false;
+      for ( size_t iP = 0; iP < quad._edgeNodes.size(); ++iP )
+        if (( std::find( ++chn.begin(), chn.end(), & quad._edgeNodes[iP]) == chn.end()) &&
+            chn.back()->IsLinked( quad._edgeNodes[ iP ]._intPoint ))
+        {
+          chn.push_back( & quad._edgeNodes[ iP ]);
+          found = true;
+          break;
+        }
+    } while ( found && chn.back() != n2 );
+
+    if ( chn.back() != n2 )
+      chn.push_back( n2 );
+
+    return chn.size() > 2;
+  }
+  //================================================================================
+  /*!
    * \brief Adds computed elements to the mesh
    */
   int Hexahedron::addElements(SMESH_MesherHelper& helper)
@@ -1646,8 +2525,19 @@ namespace
     // add elements resulted from hexahedron intersection
     //for ( size_t i = 0; i < _volumeDefs.size(); ++i )
     {
-      vector< const SMDS_MeshNode* >& nodes = _volumeDefs._nodes;
-      
+      vector< const SMDS_MeshNode* > nodes( _volumeDefs._nodes.size() );
+      for ( size_t iN = 0; iN < nodes.size(); ++iN )
+        if ( !( nodes[iN] = _volumeDefs._nodes[iN]->Node() ))
+        {
+          if ( const E_IntersectPoint* eip = _volumeDefs._nodes[iN]->EdgeIntPnt() )
+            nodes[iN] = _volumeDefs._nodes[iN]->_intPoint->_node =
+              helper.AddNode( eip->_point.X(),
+                              eip->_point.Y(),
+                              eip->_point.Z() );
+          else
+            throw SALOME_Exception("Bug: no node at intersection point");
+        }
+
       if ( !_volumeDefs._quantities.empty() )
       {
         helper.AddPolyhedralVolume( nodes, _volumeDefs._quantities );
@@ -1679,8 +2569,11 @@ namespace
    */
   bool Hexahedron::isInHole() const
   {
+    if ( !_vertexNodes.empty() )
+      return false;
+
     const int ijk[3] = { _i, _j, _k };
-    IntersectionPoint curIntPnt;
+    F_IntersectPoint curIntPnt;
 
     // consider a cell to be in a hole if all links in any direction
     // comes OUT of geometry
@@ -1698,19 +2591,19 @@ namespace
       {
         const _Link& link = _hexLinks[ iL + 4*iDir ];
         // check transition of the first node of a link
-        const IntersectionPoint* firstIntPnt = 0;
+        const F_IntersectPoint* firstIntPnt = 0;
         if ( link._nodes[0]->Node() ) // 1st node is a hexa corner
         {
           curIntPnt._paramOnLine = coords[ ijk[ iDir ]] - coords[0];
           const GridLine& line = _grid->_lines[ iDir ][ lineIndex[ iL ]];
-          multiset< IntersectionPoint >::const_iterator ip =
+          multiset< F_IntersectPoint >::const_iterator ip =
             line._intPoints.upper_bound( curIntPnt );
           --ip;
           firstIntPnt = &(*ip);
         }
         else if ( !link._intNodes.empty() )
         {
-          firstIntPnt = link._intNodes[0]._intPoint;
+          firstIntPnt = link._intNodes[0].FaceIntPnt();
         }
 
         if ( firstIntPnt )
@@ -1736,10 +2629,10 @@ namespace
     {
       const _Face& polygon = _polygons[iP];
       gp_XYZ area (0,0,0);
-      SMESH_TNodeXYZ p1 ( polygon._links[ 0 ].FirstNode()->Node() );
+      gp_XYZ p1 = polygon._links[ 0 ].FirstNode()->Point().XYZ();
       for ( size_t iL = 0; iL < polygon._links.size(); ++iL )
       {
-        SMESH_TNodeXYZ p2 ( polygon._links[ iL ].LastNode()->Node() );
+        gp_XYZ p2 = polygon._links[ iL ].LastNode()->Point().XYZ();
         area += p1 ^ p2;
         p1 = p2;
       }
@@ -1764,12 +2657,12 @@ namespace
          _polygons[4]._links.size() != 4 ||
          _polygons[5]._links.size() != 4   )
       return false;
-    const SMDS_MeshNode* nodes[8];
+    _Node* nodes[8];
     int nbN = 0;
     for ( int iL = 0; iL < 4; ++iL )
     {
       // a base node
-      nodes[iL] = _polygons[0]._links[iL].FirstNode()->Node();
+      nodes[iL] = _polygons[0]._links[iL].FirstNode();
       ++nbN;
 
       // find a top node above the base node
@@ -1781,13 +2674,13 @@ namespace
         if ( quad->_links[i]._link == link )
         {
           // 1st node of a link opposite to <link> in <quad>
-          nodes[iL+4] = quad->_links[(i+2)%4].FirstNode()->Node();
+          nodes[iL+4] = quad->_links[(i+2)%4].FirstNode();
           ++nbN;
           break;
         }
     }
     if ( nbN == 8 )
-      _volumeDefs.set( vector< const SMDS_MeshNode* >( nodes, nodes+8 ));
+      _volumeDefs.set( vector< _Node* >( nodes, nodes+8 ));
 
     return nbN == 8;
   }
@@ -1797,10 +2690,10 @@ namespace
    */
   bool Hexahedron::addTetra()
   {
-    const SMDS_MeshNode* nodes[4];
-    nodes[0] = _polygons[0]._links[0].FirstNode()->Node();
-    nodes[1] = _polygons[0]._links[1].FirstNode()->Node();
-    nodes[2] = _polygons[0]._links[2].FirstNode()->Node();
+    _Node* nodes[4];
+    nodes[0] = _polygons[0]._links[0].FirstNode();
+    nodes[1] = _polygons[0]._links[1].FirstNode();
+    nodes[2] = _polygons[0]._links[2].FirstNode();
 
     _Link* link = _polygons[0]._links[0]._link;
     ASSERT( link->_faces.size() > 1 );
@@ -1810,8 +2703,8 @@ namespace
     for ( int i = 0; i < 3; ++i )
       if ( tria->_links[i]._link == link )
       {
-        nodes[3] = tria->_links[(i+1)%3].LastNode()->Node();
-        _volumeDefs.set( vector< const SMDS_MeshNode* >( nodes, nodes+4 ));
+        nodes[3] = tria->_links[(i+1)%3].LastNode();
+        _volumeDefs.set( vector< _Node* >( nodes, nodes+4 ));
         return true;
       }
 
@@ -1831,12 +2724,12 @@ namespace
     if ( iTri < 0 ) return false;
 
     // find nodes
-    const SMDS_MeshNode* nodes[6];
+    _Node* nodes[6];
     int nbN = 0;
     for ( int iL = 0; iL < 3; ++iL )
     {
       // a base node
-      nodes[iL] = _polygons[ iTri ]._links[iL].FirstNode()->Node();
+      nodes[iL] = _polygons[ iTri ]._links[iL].FirstNode();
       ++nbN;
 
       // find a top node above the base node
@@ -1849,13 +2742,13 @@ namespace
         if ( quad->_links[i]._link == link )
         {
           // 1st node of a link opposite to <link> in <quad>
-          nodes[iL+3] = quad->_links[(i+2)%4].FirstNode()->Node();
+          nodes[iL+3] = quad->_links[(i+2)%4].FirstNode();
           ++nbN;
           break;
         }
     }
     if ( nbN == 6 )
-      _volumeDefs.set( vector< const SMDS_MeshNode* >( nodes, nodes+6 ));
+      _volumeDefs.set( vector< _Node* >( nodes, nodes+6 ));
 
     return ( nbN == 6 );
   }
@@ -1873,11 +2766,11 @@ namespace
     if ( iQuad < 0 ) return false;
 
     // find nodes
-    const SMDS_MeshNode* nodes[5];
-    nodes[0] = _polygons[iQuad]._links[0].FirstNode()->Node();
-    nodes[1] = _polygons[iQuad]._links[1].FirstNode()->Node();
-    nodes[2] = _polygons[iQuad]._links[2].FirstNode()->Node();
-    nodes[3] = _polygons[iQuad]._links[3].FirstNode()->Node();
+    _Node* nodes[5];
+    nodes[0] = _polygons[iQuad]._links[0].FirstNode();
+    nodes[1] = _polygons[iQuad]._links[1].FirstNode();
+    nodes[2] = _polygons[iQuad]._links[2].FirstNode();
+    nodes[3] = _polygons[iQuad]._links[3].FirstNode();
 
     _Link* link = _polygons[iQuad]._links[0]._link;
     ASSERT( link->_faces.size() > 1 );
@@ -1888,8 +2781,8 @@ namespace
     for ( int i = 0; i < 3; ++i )
       if ( tria->_links[i]._link == link )
       {
-        nodes[4] = tria->_links[(i+1)%3].LastNode()->Node();
-        _volumeDefs.set( vector< const SMDS_MeshNode* >( nodes, nodes+5 ));
+        nodes[4] = tria->_links[(i+1)%3].LastNode();
+        _volumeDefs.set( vector< _Node* >( nodes, nodes+5 ));
         return true;
       }
 
@@ -1930,25 +2823,37 @@ bool StdMeshers_Cartesian_3D::Compute(SMESH_Mesh &         theMesh,
   {
     Grid grid;
 
-    TopTools_MapOfShape faceMap;
-    for ( TopExp_Explorer fExp( theShape, TopAbs_FACE ); fExp.More(); fExp.Next() )
-      if ( !faceMap.Add( fExp.Current() ))
-        faceMap.Remove( fExp.Current() ); // remove a face shared by two solids
-
+    vector< TopoDS_Shape > faceVec;
+    {
+      TopTools_MapOfShape faceMap;
+      for ( TopExp_Explorer fExp( theShape, TopAbs_FACE ); fExp.More(); fExp.Next() )
+        if ( faceMap.Add( fExp.Current() )) // skip a face shared by two solids
+          faceVec.push_back( fExp.Current() );
+    }
     Bnd_Box shapeBox;
-    vector<FaceGridIntersector> facesItersectors( faceMap.Extent() );
-    TopTools_MapIteratorOfMapOfShape faceMppIt( faceMap );
-    for ( int i = 0; faceMppIt.More(); faceMppIt.Next(), ++i )
+    vector<FaceGridIntersector> facesItersectors( faceVec.size() );
+    map< TGeomID, vector< TGeomID > > edge2faceIDsMap;
+    TopExp_Explorer eExp;
+    for ( int i = 0; i < faceVec.size(); ++i )
     {
-      facesItersectors[i]._face = TopoDS::Face( faceMppIt.Key() );
-      facesItersectors[i]._grid = &grid;
+      facesItersectors[i]._face   = TopoDS::Face    ( faceVec[i] );
+      facesItersectors[i]._faceID = grid._shapes.Add( faceVec[i] );
+      facesItersectors[i]._grid   = &grid;
       shapeBox.Add( facesItersectors[i].GetFaceBndBox() );
+
+      if ( _hyp->GetToAddEdges() )
+        for ( eExp.Init( faceVec[i], TopAbs_EDGE ); eExp.More(); eExp.Next() )
+        {
+          const TopoDS_Edge& edge = TopoDS::Edge( eExp.Current() );
+          if ( !SMESH_Algo::isDegenerated( edge ))
+            edge2faceIDsMap[ grid._shapes.Add( edge )].push_back( facesItersectors[i]._faceID );
+        }
     }
 
     vector<double> xCoords, yCoords, zCoords;
     _hyp->GetCoordinates( xCoords, yCoords, zCoords, shapeBox );
 
-    grid.SetCoordinates( xCoords, yCoords, zCoords, theShape );
+    grid.SetCoordinates( xCoords, yCoords, zCoords, _hyp->GetAxisDirs(), theShape );
 
     // check if the grid encloses the shape
     if ( !_hyp->IsGridBySpacing(0) ||
@@ -2015,12 +2920,12 @@ bool StdMeshers_Cartesian_3D::Compute(SMESH_Mesh &         theMesh,
 
     // create volume elements
     Hexahedron hex( _hyp->GetSizeThreshold(), &grid );
-    int nbAdded = hex.MakeElements( helper );
+    int nbAdded = hex.MakeElements( helper, edge2faceIDsMap );
 
     SMESHDS_Mesh* meshDS = theMesh.GetMeshDS();
     if ( nbAdded > 0 )
     {
-      // make all SOLIDS computed
+      // make all SOLIDs computed
       if ( SMESHDS_SubMesh* sm1 = meshDS->MeshElements( solidExp.Current()) )
       {
         SMDS_ElemIteratorPtr volIt = sm1->GetElements();
@@ -2038,23 +2943,28 @@ bool StdMeshers_Cartesian_3D::Compute(SMESH_Mesh &         theMesh,
     // remove free nodes
     if ( SMESHDS_SubMesh * smDS = meshDS->MeshElements( helper.GetSubShapeID() ))
     {
-      // intersection nodes
+      TIDSortedNodeSet nodesToRemove;
+      // get intersection nodes
       for ( int iDir = 0; iDir < 3; ++iDir )
       {
         vector< GridLine >& lines = grid._lines[ iDir ];
         for ( size_t i = 0; i < lines.size(); ++i )
         {
-          multiset< IntersectionPoint >::iterator ip = lines[i]._intPoints.begin();
+          multiset< F_IntersectPoint >::iterator ip = lines[i]._intPoints.begin();
           for ( ; ip != lines[i]._intPoints.end(); ++ip )
             if ( ip->_node && ip->_node->NbInverseElements() == 0 )
-              meshDS->RemoveFreeNode( ip->_node, smDS, /*fromGroups=*/false );
+              nodesToRemove.insert( nodesToRemove.end(), ip->_node );
         }
       }
-      // grid nodes
+      // get grid nodes
       for ( size_t i = 0; i < grid._nodes.size(); ++i )
-        if ( !grid._isBndNode[i] ) // nodes on boundary are already removed
-          if ( grid._nodes[i] && grid._nodes[i]->NbInverseElements() == 0 )
-            meshDS->RemoveFreeNode( grid._nodes[i], smDS, /*fromGroups=*/false );
+        if ( grid._nodes[i] && grid._nodes[i]->NbInverseElements() == 0 )
+          nodesToRemove.insert( nodesToRemove.end(), grid._nodes[i] );
+
+      // do remove
+      TIDSortedNodeSet::iterator n = nodesToRemove.begin();
+      for ( ; n != nodesToRemove.end(); ++n )
+        meshDS->RemoveFreeNode( *n, smDS, /*fromGroups=*/false );
     }
 
     return nbAdded;
index 4e5e047..5382b7e 100644 (file)
@@ -251,7 +251,8 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(const StdMeshers_FaceSide*  theSide,
  */
 //================================================================================
 
-StdMeshers_FaceSide::StdMeshers_FaceSide(UVPtStructVec& theSideNodes)
+StdMeshers_FaceSide::StdMeshers_FaceSide(UVPtStructVec&     theSideNodes,
+                                         const TopoDS_Face& theFace)
 {
   myEdge.resize( 1 );
   myEdgeID.resize( 1, -1 );
@@ -272,12 +273,42 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(UVPtStructVec& theSideNodes)
   if ( !myPoints.empty() )
   {
     myPoints[0].normParam = 0;
-    gp_Pnt pPrev = SMESH_TNodeXYZ( myPoints[0].node );
-    for ( size_t i = 1; i < myPoints.size(); ++i )
+    if ( myPoints[0].node &&
+         myPoints.back().node &&
+         myPoints[ myNbPonits/2 ].node )
     {
-      gp_Pnt p = SMESH_TNodeXYZ( myPoints[i].node );
-      myLength += ( myPoints[i].normParam = p.Distance( pPrev ));
-      pPrev = p;
+      gp_Pnt pPrev = SMESH_TNodeXYZ( myPoints[0].node );
+      for ( size_t i = 1; i < myPoints.size(); ++i )
+      {
+        gp_Pnt p = SMESH_TNodeXYZ( myPoints[i].node );
+        myLength += p.Distance( pPrev );
+        myPoints[i].normParam = myLength;
+        pPrev = p;
+      }
+    }
+    else if ( !theFace.IsNull() )
+    {
+      TopLoc_Location loc;
+      Handle(Geom_Surface) surf = BRep_Tool::Surface( theFace, loc );
+      gp_Pnt pPrev = surf->Value( myPoints[0].u, myPoints[0].v );
+      for ( size_t i = 1; i < myPoints.size(); ++i )
+      {
+        gp_Pnt p = surf->Value( myPoints[i].u, myPoints[i].v );
+        myLength += p.Distance( pPrev );
+        myPoints[i].normParam = myLength;
+        pPrev = p;
+      }
+    }
+    else
+    {
+      gp_Pnt2d pPrev = myPoints[0].UV();
+      for ( size_t i = 1; i < myPoints.size(); ++i )
+      {
+        gp_Pnt2d p = myPoints[i].UV();
+        myLength += p.Distance( pPrev );
+        myPoints[i].normParam = myLength;
+        pPrev = p;
+      }
     }
     if ( myLength > std::numeric_limits<double>::min() )
       for ( size_t i = 1; i < myPoints.size(); ++i )
@@ -926,6 +957,18 @@ gp_Pnt2d StdMeshers_FaceSide::Value2d(double U) const
     return myC2d[ i ]->Value(par);
 
   }
+  else if ( !myPoints.empty() )
+  {
+    int i = U * double( myPoints.size()-1 );
+    while ( i > 0 && myPoints[ i ].normParam > U )
+      --i;
+    while ( i+1 < myPoints.size() && myPoints[ i+1 ].normParam < U )
+      ++i;
+    double r = (( U - myPoints[ i ].normParam ) /
+                ( myPoints[ i+1 ].normParam - myPoints[ i ].normParam ));
+    return ( myPoints[ i   ].UV() * ( 1 - r ) +
+             myPoints[ i+1 ].UV() * r );
+  }
   return myDefaultPnt2d;
 }
 
index 4f714d4..c08d75c 100644 (file)
@@ -35,6 +35,7 @@
 #include <Geom2d_Curve.hxx>
 #include <GeomAdaptor_Curve.hxx>
 #include <TopoDS_Edge.hxx>
+#include <TopoDS_Face.hxx>
 #include <TopoDS_Vertex.hxx>
 #include <gp_Pnt2d.hxx>
 
@@ -47,7 +48,6 @@ class SMESH_Mesh;
 class Adaptor2d_Curve2d;
 class Adaptor3d_Curve;
 class BRepAdaptor_CompCurve;
-class TopoDS_Face;
 struct SMESH_ComputeError;
 class StdMeshers_FaceSide;
 
@@ -96,7 +96,43 @@ public:
   /*!
    * \brief Create a side from an UVPtStructVec
    */
-  StdMeshers_FaceSide(UVPtStructVec& theSideNodes);
+  StdMeshers_FaceSide(UVPtStructVec&     theSideNodes,
+                      const TopoDS_Face& theFace = TopoDS_Face());
+
+  // static "consrtuctors"
+  static StdMeshers_FaceSidePtr New(const TopoDS_Face&   Face,
+                                    const TopoDS_Edge&   Edge,
+                                    SMESH_Mesh*          Mesh,
+                                    const bool           IsForward,
+                                    const bool           IgnoreMediumNodes,
+                                    SMESH_ProxyMesh::Ptr ProxyMesh = SMESH_ProxyMesh::Ptr())
+  { return StdMeshers_FaceSidePtr
+      ( new StdMeshers_FaceSide( Face,Edge,Mesh,IsForward,IgnoreMediumNodes,ProxyMesh ));
+  }
+  static StdMeshers_FaceSidePtr New (const TopoDS_Face&      Face,
+                                     std::list<TopoDS_Edge>& Edges,
+                                     SMESH_Mesh*             Mesh,
+                                     const bool              IsForward,
+                                     const bool              IgnoreMediumNodes,
+                                     SMESH_ProxyMesh::Ptr    ProxyMesh = SMESH_ProxyMesh::Ptr())
+  { return StdMeshers_FaceSidePtr
+      ( new StdMeshers_FaceSide( Face,Edges,Mesh,IsForward,IgnoreMediumNodes,ProxyMesh ));
+  }
+  static StdMeshers_FaceSidePtr New (const StdMeshers_FaceSide*  Side,
+                                     const SMDS_MeshNode*        Node,
+                                     const gp_Pnt2d*             Pnt2d1,
+                                     const gp_Pnt2d*             Pnt2d2=NULL,
+                                     const Handle(Geom2d_Curve)& C2d=NULL,
+                                     const double                UFirst=0.,
+                                     const double                ULast=1.)
+  { return StdMeshers_FaceSidePtr
+      ( new StdMeshers_FaceSide( Side,Node,Pnt2d1,Pnt2d2,C2d,UFirst,ULast ));
+  }
+  static StdMeshers_FaceSidePtr New (UVPtStructVec&     theSideNodes,
+                                     const TopoDS_Face& theFace = TopoDS_Face())
+  {
+    return StdMeshers_FaceSidePtr( new StdMeshers_FaceSide( theSideNodes, theFace ));
+  }
 
   /*!
    * \brief Return wires of a face as StdMeshers_FaceSide's
diff --git a/src/StdMeshers/StdMeshers_Geometric1D.cxx b/src/StdMeshers/StdMeshers_Geometric1D.cxx
new file mode 100644 (file)
index 0000000..8d2be97
--- /dev/null
@@ -0,0 +1,204 @@
+// Copyright (C) 2007-2013  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 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
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+
+//  SMESH SMESH : implementaion of SMESH idl descriptions
+//  File   : StdMeshers_Geometric1D.cxx
+//  Module : SMESH
+//
+#include "StdMeshers_Geometric1D.hxx"
+
+#include "SMESH_Mesh.hxx"
+
+#include <BRepAdaptor_Curve.hxx>
+#include <GCPnts_AbscissaPoint.hxx>
+#include <SMESH_Algo.hxx>
+#include <TopExp.hxx>
+#include <TopTools_IndexedMapOfShape.hxx>
+#include <TopoDS.hxx>
+#include <TopoDS_Edge.hxx>
+
+//=============================================================================
+/*!
+ * Constructor
+ */
+//=============================================================================
+
+StdMeshers_Geometric1D::StdMeshers_Geometric1D(int hypId, int studyId, SMESH_Gen * gen)
+  :StdMeshers_Reversible1D(hypId, studyId, gen)
+{
+  _begLength = 1.;
+  _ratio = 1.;
+  _name = "GeometricProgression";
+}
+
+//=============================================================================
+/*!
+ * Sets length of the first segment
+ */
+//=============================================================================
+
+void StdMeshers_Geometric1D::SetStartLength(double length)
+  throw(SALOME_Exception)
+{
+  if ( _begLength != length )
+  {
+    if (length <= 0)
+      throw SALOME_Exception(LOCALIZED("length must be positive"));
+    _begLength = length;
+    NotifySubMeshesHypothesisModification();
+  }
+}
+
+//=============================================================================
+/*!
+ * Sets value of Common Ratio
+ */
+//=============================================================================
+
+void StdMeshers_Geometric1D::SetCommonRatio(double factor)
+  throw(SALOME_Exception)
+{
+  if ( _ratio != factor )
+  {
+    if (factor == 0)
+      throw SALOME_Exception(LOCALIZED("Zero factor is not allowed"));
+    _ratio = factor;
+    NotifySubMeshesHypothesisModification();
+  }
+}
+
+//=============================================================================
+/*!
+ * Returns length of the first segment 
+ */
+//=============================================================================
+
+double StdMeshers_Geometric1D::GetStartLength() const
+{
+  return _begLength;
+}
+
+//=============================================================================
+/*!
+ * Returns value of Common Ratio
+ */
+//=============================================================================
+
+double StdMeshers_Geometric1D::GetCommonRatio() const
+{
+  return _ratio;
+}
+
+//=============================================================================
+/*!
+ *  
+ */
+//=============================================================================
+
+ostream & StdMeshers_Geometric1D::SaveTo(ostream & save)
+{
+  save << _begLength << " " << _ratio << " ";
+
+  StdMeshers_Reversible1D::SaveTo( save );
+
+  return save;
+}
+
+//=============================================================================
+/*!
+ *  
+ */
+//=============================================================================
+
+istream & StdMeshers_Geometric1D::LoadFrom(istream & load)
+{
+  bool isOK = true;
+  isOK = (load >> _begLength);
+  isOK = (load >> _ratio);
+
+  if (isOK)
+    StdMeshers_Reversible1D::LoadFrom( load );
+
+  return load;
+}
+
+//================================================================================
+/*!
+ * \brief Initialize start and end length by the mesh built on the geometry
+ * \param theMesh - the built mesh
+ * \param theShape - the geometry of interest
+ * \retval bool - true if parameter values have been successfully defined
+ */
+//================================================================================
+
+bool StdMeshers_Geometric1D::SetParametersByMesh(const SMESH_Mesh*   theMesh,
+                                                 const TopoDS_Shape& theShape)
+{
+  if ( !theMesh || theShape.IsNull() )
+    return false;
+
+  _begLength = _ratio = 0.;
+
+  int nbEdges = 0;
+  TopTools_IndexedMapOfShape edgeMap;
+  TopExp::MapShapes( theShape, TopAbs_EDGE, edgeMap );
+  for ( int i = 1; i <= edgeMap.Extent(); ++i )
+  {
+    const TopoDS_Edge& edge = TopoDS::Edge( edgeMap( i ));
+    BRepAdaptor_Curve C( edge );
+
+    vector< double > params;
+    if ( SMESH_Algo::GetNodeParamOnEdge( theMesh->GetMeshDS(), edge, params ))
+    {
+      nbEdges++;
+      double l1 = GCPnts_AbscissaPoint::Length( C, params[0], params[1] );
+      _begLength += l1;
+      if ( params.size() > 2 && l1 > 1e-100 )
+        _ratio += GCPnts_AbscissaPoint::Length( C, params[1], params[2]) / l1;
+      else
+        _ratio += 1;
+    }
+  }
+  if ( nbEdges ) {
+    _begLength /= nbEdges;
+    _ratio     /= nbEdges;
+  }
+  else {
+    _begLength = 1;
+    _ratio     = 1;
+  }
+  return nbEdges;
+}
+
+//================================================================================
+/*!
+ * \brief Initialize my parameter values by default parameters.
+ *  \retval bool - true if parameter values have been successfully defined
+ */
+//================================================================================
+
+bool StdMeshers_Geometric1D::SetParametersByDefaults(const TDefaults&  dflts,
+                                                     const SMESH_Mesh* /*mesh*/)
+{
+  return ( _begLength = dflts._elemLength );
+}
+
diff --git a/src/StdMeshers/StdMeshers_Geometric1D.hxx b/src/StdMeshers/StdMeshers_Geometric1D.hxx
new file mode 100644 (file)
index 0000000..436c317
--- /dev/null
@@ -0,0 +1,67 @@
+// Copyright (C) 2007-2013  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 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
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+
+//  SMESH SMESH : implementaion of SMESH idl descriptions
+//  File   : StdMeshers_Geometric1D.hxx
+//  Module : SMESH
+//
+#ifndef _SMESH_Geometric1D_HXX_
+#define _SMESH_Geometric1D_HXX_
+
+#include "SMESH_StdMeshers.hxx"
+
+#include "StdMeshers_Reversible1D.hxx"
+#include "Utils_SALOME_Exception.hxx"
+
+class STDMESHERS_EXPORT StdMeshers_Geometric1D: public StdMeshers_Reversible1D
+{
+public:
+  StdMeshers_Geometric1D(int hypId, int studyId, SMESH_Gen* gen);
+
+  void SetStartLength(double length) throw(SALOME_Exception);
+  void SetCommonRatio(double factor) throw(SALOME_Exception);
+
+  double GetStartLength() const;
+  double GetCommonRatio() const;
+
+  virtual std::ostream & SaveTo(std::ostream & save);
+  virtual std::istream & LoadFrom(std::istream & load);
+
+  /*!
+   * \brief Initialize start and end length by the mesh built on the geometry
+    * \param theMesh - the built mesh
+    * \param theShape - the geometry of interest
+    * \retval bool - true if parameter values have been successfully defined
+   */
+  virtual bool SetParametersByMesh(const SMESH_Mesh* theMesh, const TopoDS_Shape& theShape);
+
+  /*!
+   * \brief Initialize my parameter values by default parameters.
+   *  \retval bool - true if parameter values have been successfully defined
+   */
+  virtual bool SetParametersByDefaults(const TDefaults& dflts, const SMESH_Mesh* theMesh=0);
+
+protected:
+  double _begLength, _ratio;
+};
+
+#endif
index 949ba67..a0d205b 100644 (file)
@@ -222,8 +222,8 @@ namespace
    * \brief Finds FaceQuadStruct having a side equal to a given one and rearranges
    *  the found FaceQuadStruct::side to have the given side at a Q_BOTTOM place
    */
-  FaceQuadStructPtr getQuadWithBottom( StdMeshers_FaceSide* side,
-                                       FaceQuadStructPtr    quad[ 6 ])
+  FaceQuadStructPtr getQuadWithBottom( StdMeshers_FaceSidePtr side,
+                                       FaceQuadStructPtr      quad[ 6 ])
   {
     FaceQuadStructPtr foundQuad;
     for ( int i = 1; i < 6; ++i )
@@ -231,7 +231,7 @@ namespace
       if ( !quad[i] ) continue;
       for ( unsigned iS = 0; iS < quad[i]->side.size(); ++iS )
       {
-        const StdMeshers_FaceSide* side2 = quad[i]->side[iS];
+        const StdMeshers_FaceSidePtr side2 = quad[i]->side[iS];
         if (( side->FirstVertex().IsSame( side2->FirstVertex() ) ||
               side->FirstVertex().IsSame( side2->LastVertex() ))
             &&
@@ -241,7 +241,7 @@ namespace
         {
           if ( iS != Q_BOTTOM )
           {
-            vector< StdMeshers_FaceSide*> newSides;
+            vector< FaceQuadStruct::Side > newSides;
             for ( unsigned j = iS; j < quad[i]->side.size(); ++j )
               newSides.push_back( quad[i]->side[j] );
             for ( unsigned j = 0; j < iS; ++j )
@@ -391,7 +391,7 @@ bool StdMeshers_Hexa_3D::Compute(SMESH_Mesh &         aMesh,
   for ( int i = 0; i < 6; ++i )
   {
     const TopoDS_Face& F = aCubeSide[i]._quad->face;
-    StdMeshers_FaceSide* baseQuadSide = aCubeSide[i]._quad->side[ Q_BOTTOM ];
+    StdMeshers_FaceSidePtr baseQuadSide = aCubeSide[i]._quad->side[ Q_BOTTOM ];
     list<TopoDS_Edge> baseEdges( baseQuadSide->Edges().begin(), baseQuadSide->Edges().end() );
 
     // assure correctness of node positions on baseE:
index 07e9a16..33bd5e1 100644 (file)
@@ -107,7 +107,7 @@ namespace {
            algo->myProxyMesh->GetMesh() != helper->GetMesh() )
         algo->myProxyMesh.reset( new SMESH_ProxyMesh( *helper->GetMesh() ));
 
-      algo->myQuadStruct.reset();
+      algo->myQuadList.clear();
 
       if ( helper )
         algo->_quadraticMesh = helper->GetIsQuadratic();
@@ -166,15 +166,15 @@ namespace {
   //================================================================================
 
   bool setBottomEdge( const TopoDS_Edge&   botE,
-                      faceQuadStruct::Ptr& quad,
+                      FaceQuadStruct::Ptr& quad,
                       const TopoDS_Shape&  face)
   {
-    quad->side[ QUAD_TOP_SIDE  ]->Reverse();
-    quad->side[ QUAD_LEFT_SIDE ]->Reverse();
+    quad->side[ QUAD_TOP_SIDE  ].grid->Reverse();
+    quad->side[ QUAD_LEFT_SIDE ].grid->Reverse();
     int edgeIndex = 0;
     for ( size_t i = 0; i < quad->side.size(); ++i )
     {
-      StdMeshers_FaceSide* quadSide = quad->side[i];
+      StdMeshers_FaceSidePtr quadSide = quad->side[i];
       for ( int iE = 0; iE < quadSide->NbEdges(); ++iE )
         if ( botE.IsSame( quadSide->Edge( iE )))
         {
@@ -681,7 +681,7 @@ bool StdMeshers_Prism_3D::Compute(SMESH_Mesh& theMesh, const TopoDS_Shape& theSh
               continue; // already computed prism
             }
             // find a source FACE of the SOLID: it's a FACE sharing a bottom EDGE with wFace
-            const TopoDS_Edge& wEdge = (*wQuad)->side[ QUAD_TOP_SIDE ]->Edge(0);
+            const TopoDS_Edge& wEdge = (*wQuad)->side[ QUAD_TOP_SIDE ].grid->Edge(0);
             PShapeIteratorPtr faceIt = myHelper->GetAncestors( wEdge, *myHelper->GetMesh(),
                                                                TopAbs_FACE);
             while ( const TopoDS_Shape* f = faceIt->next() )
@@ -879,7 +879,7 @@ bool StdMeshers_Prism_3D::getWallFaces( Prism_3D::TPrismTopo & thePrism,
     int nbKnownFaces;
     do {
       nbKnownFaces = faceMap.Extent();
-      StdMeshers_FaceSide *rightSide, *topSide; // sides of the quad
+      StdMeshers_FaceSidePtr rightSide, topSide; // sides of the quad
       for ( size_t i = 0; i < thePrism.myWallQuads.size(); ++i )
       {
         rightSide = thePrism.myWallQuads[i].back()->side[ QUAD_RIGHT_SIDE ];
@@ -911,8 +911,8 @@ bool StdMeshers_Prism_3D::getWallFaces( Prism_3D::TPrismTopo & thePrism,
     {
       for ( size_t i = 0; i < thePrism.myWallQuads.size(); ++i )
       {
-        StdMeshers_FaceSide* topSide = thePrism.myWallQuads[i].back()->side[ QUAD_TOP_SIDE ];
-        const TopoDS_Edge &     topE = topSide->Edge( 0 );
+        StdMeshers_FaceSidePtr topSide = thePrism.myWallQuads[i].back()->side[ QUAD_TOP_SIDE ];
+        const TopoDS_Edge &       topE = topSide->Edge( 0 );
         if ( topSide->NbEdges() > 1 )
           return toSM( error(COMPERR_BAD_SHAPE, TCom("Side face #") <<
                              shapeID( thePrism.myWallQuads[i].back()->face )
@@ -958,8 +958,8 @@ bool StdMeshers_Prism_3D::getWallFaces( Prism_3D::TPrismTopo & thePrism,
   // Check that the top FACE shares all the top EDGEs
   for ( size_t i = 0; i < thePrism.myWallQuads.size(); ++i )
   {
-    StdMeshers_FaceSide* topSide = thePrism.myWallQuads[i].back()->side[ QUAD_TOP_SIDE ];
-    const TopoDS_Edge &     topE = topSide->Edge( 0 );
+    StdMeshers_FaceSidePtr topSide = thePrism.myWallQuads[i].back()->side[ QUAD_TOP_SIDE ];
+    const TopoDS_Edge &       topE = topSide->Edge( 0 );
     if ( !myHelper->IsSubShape( topE, thePrism.myTop ))
       return toSM( error( TCom("Wrong source face (#") << shapeID( thePrism.myBottom )));
   }
@@ -1205,7 +1205,7 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism)
     int wgt = 0; // "weight"
     for ( ; quad != thePrism.myWallQuads[iW].end(); ++quad )
     {
-      StdMeshers_FaceSide* lftSide = (*quad)->side[ QUAD_LEFT_SIDE ];
+      StdMeshers_FaceSidePtr lftSide = (*quad)->side[ QUAD_LEFT_SIDE ];
       for ( int i = 0; i < lftSide->NbEdges(); ++i )
       {
         ++wgt;
@@ -1224,7 +1224,7 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism)
       quad = thePrism.myWallQuads[iW].begin();
       for ( ; quad != thePrism.myWallQuads[iW].end(); ++quad )
         for ( int i = 0; i < NB_QUAD_SIDES; ++i )
-          (*quad)->side[ i ]->SetIgnoreMediumNodes( true );
+          (*quad)->side[ i ].grid->SetIgnoreMediumNodes( true );
     }
   }
 
@@ -1237,8 +1237,8 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism)
     Prism_3D::TQuadList::const_iterator quad = quads.begin();
     for ( ; quad != quads.end(); ++quad )
     {
-      StdMeshers_FaceSide* rgtSide = (*quad)->side[ QUAD_RIGHT_SIDE ]; // tgt
-      StdMeshers_FaceSide* lftSide = (*quad)->side[ QUAD_LEFT_SIDE ];  // src
+      StdMeshers_FaceSidePtr rgtSide = (*quad)->side[ QUAD_RIGHT_SIDE ]; // tgt
+      StdMeshers_FaceSidePtr lftSide = (*quad)->side[ QUAD_LEFT_SIDE ];  // src
       bool swapLeftRight = ( lftSide->NbSegments( /*update=*/true ) == 0 &&
                              rgtSide->NbSegments( /*update=*/true )  > 0 );
       if ( swapLeftRight )
@@ -1373,8 +1373,8 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism)
       // to compute stuctured quad mesh on wall FACEs
       // ---------------------------------------------------
       {
-        const TopoDS_Edge& botE = (*quad)->side[ QUAD_BOTTOM_SIDE ]->Edge(0);
-        const TopoDS_Edge& topE = (*quad)->side[ QUAD_TOP_SIDE    ]->Edge(0);
+        const TopoDS_Edge& botE = (*quad)->side[ QUAD_BOTTOM_SIDE ].grid->Edge(0);
+        const TopoDS_Edge& topE = (*quad)->side[ QUAD_TOP_SIDE    ].grid->Edge(0);
         SMESH_subMesh*    botSM = mesh->GetSubMesh( botE );
         SMESH_subMesh*    topSM = mesh->GetSubMesh( topE );
         SMESH_subMesh*    srcSM = botSM;
@@ -2352,7 +2352,7 @@ bool StdMeshers_PrismAsBlock::Init(SMESH_MesherHelper*         helper,
     Prism_3D::TQuadList::const_iterator quad = thePrism.myWallQuads[ iE ].begin();
     for ( ; quad != thePrism.myWallQuads[ iE ].end(); ++quad )
     {
-      const TopoDS_Edge& quadBot = (*quad)->side[ QUAD_BOTTOM_SIDE ]->Edge( 0 );
+      const TopoDS_Edge& quadBot = (*quad)->side[ QUAD_BOTTOM_SIDE ].grid->Edge( 0 );
       if ( !myHelper->LoadNodeColumns( faceColumns, (*quad)->face, quadBot, meshDS ))
         return error(COMPERR_BAD_INPUT_MESH, TCom("Can't find regular quadrangle mesh ")
                      << "on a side face #" << MeshDS()->ShapeToIndex( (*quad)->face ));
@@ -2373,7 +2373,7 @@ bool StdMeshers_PrismAsBlock::Init(SMESH_MesherHelper*         helper,
     Prism_3D::TQuadList::const_iterator quad = thePrism.myWallQuads[ iE ].begin();
     for ( ; quad != thePrism.myWallQuads[ iE ].end(); ++quad )
     {
-      const TopoDS_Edge& quadBot = (*quad)->side[ QUAD_BOTTOM_SIDE ]->Edge( 0 );
+      const TopoDS_Edge& quadBot = (*quad)->side[ QUAD_BOTTOM_SIDE ].grid->Edge( 0 );
       if ( !myHelper->LoadNodeColumns( faceColumns, (*quad)->face, quadBot, meshDS ))
         return error(COMPERR_BAD_INPUT_MESH, TCom("Can't find regular quadrangle mesh ")
                      << "on a side face #" << MeshDS()->ShapeToIndex( (*quad)->face ));
index 6bb6adc..96570ec 100644 (file)
@@ -1018,9 +1018,9 @@ bool StdMeshers_Projection_2D::Compute(SMESH_Mesh& theMesh, const TopoDS_Shape&
       map< double, const SMDS_MeshNode* >::iterator u_oldNode, u_newNode, u_newOnSeam, newEnd;
       set< const SMDS_MeshNode* > seamNodes;
 
-      // mapper puts on a seam edge nodes from 2 edges
+      // mapper changed, no more "mapper puts on a seam edge nodes from 2 edges"
       if ( isSeam && ! getBoundaryNodes ( sm, tgtFace, u2nodesOnSeam, seamNodes ))
-        RETURN_BAD_RESULT("getBoundaryNodes() failed");
+        ;//RETURN_BAD_RESULT("getBoundaryNodes() failed");
 
       SMDS_NodeIteratorPtr nIt = smDS->GetNodes();
       while ( nIt->more() )
index bfe999d..15bd194 100644 (file)
@@ -64,7 +64,8 @@ namespace {
     /*!
      * \brief Return an edge from which hypotheses are propagated from
      */
-    static TopoDS_Edge GetSource(SMESH_subMesh * submesh);
+    static TopoDS_Edge GetSource(SMESH_subMesh * submesh,
+                                 bool&           isPropagOfDistribution);
     /*!
      * \brief Does it's main job
      */
@@ -90,23 +91,28 @@ StdMeshers_Propagation::StdMeshers_Propagation (int hypId, int studyId, SMESH_Ge
   _name = GetName();
   _param_algo_dim = -1; // 1D auxiliary
 }
+StdMeshers_PropagOfDistribution::StdMeshers_PropagOfDistribution (int hypId,
+                                                                  int studyId,
+                                                                  SMESH_Gen * gen)
+  : StdMeshers_Propagation(hypId, studyId, gen)                        { _name = GetName(); }
 StdMeshers_Propagation::~StdMeshers_Propagation()                      {}
 string StdMeshers_Propagation::GetName ()                              { return "Propagation"; }
+string StdMeshers_PropagOfDistribution::GetName ()                     { return "PropagOfDistribution"; }
 ostream & StdMeshers_Propagation::SaveTo (ostream & save)              { return save; }
 istream & StdMeshers_Propagation::LoadFrom (istream & load)            { return load; }
-ostream & operator << (ostream & save, StdMeshers_Propagation & hyp)   { return hyp.SaveTo(save); }
-istream & operator >> (istream & load, StdMeshers_Propagation & hyp)   { return hyp.LoadFrom(load); }
 bool StdMeshers_Propagation::SetParametersByMesh(const SMESH_Mesh*,
                                                  const TopoDS_Shape& ) { return false; }
 bool StdMeshers_Propagation::SetParametersByDefaults(const TDefaults&,const SMESH_Mesh*) { return false; }
 void StdMeshers_Propagation::SetPropagationMgr(SMESH_subMesh* subMesh) { PropagationMgr::Set( subMesh ); }
 /*!
- * \brief Return an edge from which hypotheses are propagated from
+ * \brief Return an edge from which hypotheses are propagated
  */
-TopoDS_Edge StdMeshers_Propagation::GetPropagationSource(SMESH_Mesh& theMesh,
-                                                         const TopoDS_Shape& theEdge)
+TopoDS_Edge StdMeshers_Propagation::GetPropagationSource(SMESH_Mesh&         theMesh,
+                                                         const TopoDS_Shape& theEdge,
+                                                         bool&               isPropagOfDistribution)
 {
-  return PropagationMgr::GetSource(theMesh.GetSubMeshContaining( theEdge ));
+  return PropagationMgr::GetSource( theMesh.GetSubMeshContaining( theEdge ),
+                                    isPropagOfDistribution);
 }
 
 //=============================================================================
@@ -126,11 +132,13 @@ namespace {
   struct PropagationMgrData : public EventListenerData
   {
     bool myForward; //!< true if a curve of edge in chain is codirected with one of source edge
+    bool myIsPropagOfDistribution; //!< type of Propagation hyp
     PropagationMgrData( SubMeshState state=WAIT_PROPAG_HYP ): EventListenerData(true) {
-      myType = state; myForward = true;
+      myType = state; myForward = true; myIsPropagOfDistribution = false;
     }
     void Init() {
       myType = WAIT_PROPAG_HYP;  mySubMeshes.clear(); myForward = true;
+      myIsPropagOfDistribution = false;
     }
     SubMeshState State() const {
       return (SubMeshState) myType;
@@ -217,8 +225,13 @@ namespace {
   const SMESH_Hypothesis* getProagationHyp (SMESH_Mesh&         theMesh,
                                             const TopoDS_Shape& theEdge)
   {
-    static SMESH_HypoFilter propagHypFilter
-      ( SMESH_HypoFilter::HasName( StdMeshers_Propagation::GetName ()));
+    static SMESH_HypoFilter propagHypFilter;
+    if ( propagHypFilter.IsEmpty() )
+    {
+      propagHypFilter.
+        Init( SMESH_HypoFilter::HasName( StdMeshers_Propagation::GetName ())).
+        Or  ( SMESH_HypoFilter::HasName( StdMeshers_PropagOfDistribution::GetName ()));
+    }
     return theMesh.GetHypothesis( theEdge, propagHypFilter, true );
   }
   //================================================================================
@@ -248,6 +261,10 @@ namespace {
     PropagationMgrData* chainData = getData( theMainSubMesh );
     chainData->SetState( HAS_PROPAG_HYP );
 
+    if ( const SMESH_Hypothesis * propagHyp = getProagationHyp( *mesh, theMainEdge ))
+      chainData->myIsPropagOfDistribution =
+        ( StdMeshers_PropagOfDistribution::GetName() == propagHyp->GetName() );
+
     // Edge submeshes, to which the 1D hypothesis will be propagated from theMainEdge
     list<SMESH_subMesh*> & chain = chainData->mySubMeshes;
     chain.clear();
@@ -462,17 +479,21 @@ namespace {
   {
     if ( findData( submesh )) return;
     DBGMSG( "PropagationMgr::Set() on  " << submesh->GetId() );
-    EventListenerData* data = new PropagationMgrData();
+    PropagationMgrData* data = new PropagationMgrData();
     submesh->SetEventListener( getListener(), data, submesh );
 
     const SMESH_Hypothesis * propagHyp =
       getProagationHyp( *submesh->GetFather(), submesh->GetSubShape() );
     if ( propagHyp )
+    {
+      data->myIsPropagOfDistribution =
+        ( StdMeshers_PropagOfDistribution::GetName() == propagHyp->GetName() );
       getListener()->ProcessEvent( SMESH_subMesh::ADD_HYP,
                                    SMESH_subMesh::ALGO_EVENT,
                                    submesh,
                                    data,
                                    propagHyp);
+    }
   }
   //================================================================================
   /*!
@@ -480,7 +501,8 @@ namespace {
    */
   //================================================================================
 
-  TopoDS_Edge PropagationMgr::GetSource(SMESH_subMesh * submesh)
+  TopoDS_Edge PropagationMgr::GetSource(SMESH_subMesh * submesh,
+                                        bool&           isPropagOfDistribution)
   {
     if ( PropagationMgrData* data = findData( submesh )) {
       if ( data->State() == IN_CHAIN ) {
@@ -489,6 +511,9 @@ namespace {
           TopoDS_Shape edge = sm->GetSubShape();
           edge = edge.Oriented( data->myForward ? TopAbs_FORWARD : TopAbs_REVERSED );
           DBGMSG( " GetSource() = edge " << sm->GetId() << " REV = " << (!data->myForward));
+          isPropagOfDistribution = false;
+          if ( PropagationMgrData* data = findData( sm ))
+            isPropagOfDistribution = data->myIsPropagOfDistribution;
           if ( edge.ShapeType() == TopAbs_EDGE )
             return TopoDS::Edge( edge );
         }
@@ -502,9 +527,9 @@ namespace {
    */
   //================================================================================
 
-  void PropagationMgr::ProcessEvent(const int          event,
-                                    const int          eventType,
-                                    SMESH_subMesh*     subMesh,
+  void PropagationMgr::ProcessEvent(const int                       event,
+                                    const int                       eventType,
+                                    SMESH_subMesh*                  subMesh,
                                     SMESH_subMeshEventListenerData* listenerData,
                                     const SMESH_Hypothesis*         hyp)
   {
@@ -516,7 +541,8 @@ namespace {
       return;
     DBGMSG( "PropagationMgr::ProcessEvent() on  " << subMesh->GetId() );
 
-    bool isPropagHyp = ( StdMeshers_Propagation::GetName() == hyp->GetName() );
+    bool isPropagHyp = ( StdMeshers_Propagation::GetName()          == hyp->GetName() ||
+                         StdMeshers_PropagOfDistribution::GetName() == hyp->GetName() );
 
     PropagationMgrData* data = static_cast<PropagationMgrData*>( listenerData );
     switch ( data->State() ) {
index c5a1616..e8a2b6e 100644 (file)
@@ -50,8 +50,6 @@ class STDMESHERS_EXPORT StdMeshers_Propagation:public SMESH_Hypothesis
 
   virtual std::ostream & SaveTo(std::ostream & save);
   virtual std::istream & LoadFrom(std::istream & load);
-  friend std::ostream & operator <<(std::ostream & save, StdMeshers_Propagation & hyp);
-  friend std::istream & operator >>(std::istream & load, StdMeshers_Propagation & hyp);
 
   static std::string GetName ();
 
@@ -69,7 +67,9 @@ class STDMESHERS_EXPORT StdMeshers_Propagation:public SMESH_Hypothesis
     * \param theEdge - edge to which hypotheses are propagated
     * \retval TopoDS_Edge - source edge, also passing orientation
    */
-  static TopoDS_Edge GetPropagationSource(SMESH_Mesh& theMesh, const TopoDS_Shape& theEdge);
+  static TopoDS_Edge GetPropagationSource(SMESH_Mesh&         theMesh,
+                                          const TopoDS_Shape& theEdge,
+                                          bool&               isPropagOfDistribution );
 
   /*!
    * \brief Initialize my parameter values by the mesh built on the geometry
@@ -88,4 +88,19 @@ class STDMESHERS_EXPORT StdMeshers_Propagation:public SMESH_Hypothesis
   virtual bool SetParametersByDefaults(const TDefaults& dflts, const SMESH_Mesh* theMesh=0);
 
 };
+
+// =======================================================================
+/*!
+ * \brief Propagation Of Distribution hypothesis
+ */
+// =======================================================================
+
+class STDMESHERS_EXPORT StdMeshers_PropagOfDistribution: public StdMeshers_Propagation
+{
+ public:
+  StdMeshers_PropagOfDistribution(int hypId, int studyId, SMESH_Gen * gen);
+
+  static std::string GetName();
+};
+
 #endif
index 5369e59..e8ab5ad 100644 (file)
@@ -87,6 +87,43 @@ void StdMeshers_QuadrangleParams::SetQuadType (StdMeshers_QuadType type)
   }
 }
 
+//================================================================================
+/*!
+ * \brief Set positions of enforced nodes
+ */
+//================================================================================
+
+void StdMeshers_QuadrangleParams::
+SetEnforcedNodes( const std::vector< TopoDS_Shape >& shapes,
+                  const std::vector< gp_Pnt >&       points )
+{
+  bool isChanged = ( shapes        != _enforcedVertices ||
+                     points.size() != _enforcedPoints.size() );
+  for ( size_t i = 0; i < points.size() && !isChanged; ++i )
+    isChanged = ( _enforcedPoints[ i ].SquareDistance( points[i] ) > 1e-100 );
+      
+  if ( isChanged )
+  {
+    _enforcedVertices = shapes;
+    _enforcedPoints   = points;
+    NotifySubMeshesHypothesisModification();
+  }
+}
+
+//================================================================================
+/*!
+ * \brief Returns positions of enforced nodes
+ */
+//================================================================================
+
+void StdMeshers_QuadrangleParams::
+GetEnforcedNodes( std::vector< TopoDS_Shape >& shapes,
+                  std::vector< gp_Pnt >&       points ) const
+{
+  shapes = _enforcedVertices;
+  points = _enforcedPoints;
+}
+
 //=============================================================================
 /*!
  *
@@ -98,6 +135,13 @@ ostream & StdMeshers_QuadrangleParams::SaveTo(ostream & save)
     save << _triaVertexID << " UNDEFINED " << int(_quadType);
   else
     save << _triaVertexID << " " << _objEntry << " " << int(_quadType);
+
+  save << " " << _enforcedPoints.size();
+  for ( size_t i = 0; i < _enforcedPoints.size(); ++i )
+    save << " " << _enforcedPoints[i].X()
+         << " " << _enforcedPoints[i].Y()
+         << " " << _enforcedPoints[i].Z();
+
   return save;
 }
 
@@ -122,29 +166,25 @@ istream & StdMeshers_QuadrangleParams::LoadFrom(istream & load)
   if (isOK)
     _quadType = StdMeshers_QuadType(type);
 
+  // _enforcedVertices are loaded at StdMeshers_I level
+  // because GEOM objects are referred by study entry.
+
+  int nbP = 0;
+  double x,y,z;
+  if ( load >> nbP && nbP > 0 )
+  {
+    _enforcedPoints.reserve( nbP );
+    while ( _enforcedPoints.size() < _enforcedPoints.capacity() )
+      if ( load >> x &&
+           load >> y &&
+           load >> z )
+        _enforcedPoints.push_back( gp_Pnt( x,y,z ));
+      else
+        break;
+  }
   return load;
 }
 
-//=============================================================================
-/*!
- *
- */
-//=============================================================================
-ostream & operator <<(ostream & save, StdMeshers_QuadrangleParams & hyp)
-{
-  return hyp.SaveTo( save );
-}
-
-//=============================================================================
-/*!
- *
- */
-//=============================================================================
-istream & operator >>(istream & load, StdMeshers_QuadrangleParams & hyp)
-{
-  return hyp.LoadFrom( load );
-}
-
 //================================================================================
 /*!
  * \brief Redifined method
index 9d8374d..afc6a85 100644 (file)
 #define _SMESH_QUADRANGLEPARAMS_HXX_
 
 #include "SMESH_StdMeshers.hxx"
-
 #include "SMESH_Hypothesis.hxx"
-#include "Utils_SALOME_Exception.hxx"
+
+#include <gp_Pnt.hxx>
+
+#include <vector>
+#include <string>
 
 enum StdMeshers_QuadType
   {
@@ -38,8 +41,7 @@ enum StdMeshers_QuadType
     QUAD_NB_TYPES
   };
 
-class STDMESHERS_EXPORT StdMeshers_QuadrangleParams:
-  public SMESH_Hypothesis
+class STDMESHERS_EXPORT StdMeshers_QuadrangleParams: public SMESH_Hypothesis
 {
 public:
   StdMeshers_QuadrangleParams(int hypId, int studyId, SMESH_Gen* gen);
@@ -54,12 +56,13 @@ public:
   void SetQuadType (StdMeshers_QuadType type);
   StdMeshers_QuadType GetQuadType() const { return _quadType; }
 
+  void SetEnforcedNodes( const std::vector< TopoDS_Shape >& shapes,
+                         const std::vector< gp_Pnt >&       points );
+  void GetEnforcedNodes( std::vector< TopoDS_Shape >& shapes,
+                         std::vector< gp_Pnt >&       points ) const;
+
   virtual std::ostream & SaveTo(std::ostream & save);
   virtual std::istream & LoadFrom(std::istream & load);
-  friend std::ostream& operator << (std::ostream & save,
-                                    StdMeshers_QuadrangleParams & hyp);
-  friend std::istream& operator >> (std::istream & load,
-                                    StdMeshers_QuadrangleParams & hyp);
 
   /*!
    * \brief Initialize start and end length by the mesh built on the geometry
@@ -78,9 +81,11 @@ public:
                                        const SMESH_Mesh* theMesh=0);
 
 protected:
-  int                 _triaVertexID;
-  std::string         _objEntry;
-  StdMeshers_QuadType _quadType;
+  int                         _triaVertexID;
+  std::string                 _objEntry;
+  StdMeshers_QuadType         _quadType;
+  std::vector< TopoDS_Shape > _enforcedVertices;
+  std::vector< gp_Pnt >       _enforcedPoints;
 };
 
 #endif
index 62fadc5..1275988 100644 (file)
 #include "SMESH_Block.hxx"
 #include "SMESH_Comment.hxx"
 #include "SMESH_Gen.hxx"
+#include "SMESH_HypoFilter.hxx"
 #include "SMESH_Mesh.hxx"
+#include "SMESH_MeshAlgos.hxx"
 #include "SMESH_MesherHelper.hxx"
 #include "SMESH_subMesh.hxx"
 #include "StdMeshers_FaceSide.hxx"
 #include "StdMeshers_QuadrangleParams.hxx"
 #include "StdMeshers_ViscousLayers2D.hxx"
 
+#include <BRepClass_FaceClassifier.hxx>
 #include <BRep_Tool.hxx>
 #include <GeomAPI_ProjectPointOnSurf.hxx>
 #include <Geom_Surface.hxx>
@@ -85,8 +88,9 @@ StdMeshers_Quadrangle_2D::StdMeshers_Quadrangle_2D (int hypId, int studyId,
     myTrianglePreference(false),
     myTriaVertexID(-1),
     myNeedSmooth(false),
+    myParams( NULL ),
     myQuadType(QUAD_STANDARD),
-    myHelper( 0 )
+    myHelper( NULL )
 {
   MESSAGE("StdMeshers_Quadrangle_2D::StdMeshers_Quadrangle_2D");
   _name = "Quadrangle_2D";
@@ -119,15 +123,16 @@ bool StdMeshers_Quadrangle_2D::CheckHypothesis
                           const TopoDS_Shape&                  aShape,
                           SMESH_Hypothesis::Hypothesis_Status& aStatus)
 {
-  myTriaVertexID = -1;
-  myQuadType = QUAD_STANDARD;
+  myTriaVertexID         = -1;
+  myQuadType             = QUAD_STANDARD;
   myQuadranglePreference = false;
-  myTrianglePreference = false;
-  myQuadStruct.reset();
-  myHelper = NULL;
+  myTrianglePreference   = false;
+  myHelper               = (SMESH_MesherHelper*)NULL;
+  myParams               = NULL;
+  myQuadList.clear();
 
   bool isOk = true;
-  aStatus = SMESH_Hypothesis::HYP_OK;
+  aStatus   = SMESH_Hypothesis::HYP_OK;
 
   const list <const SMESHDS_Hypothesis * >& hyps =
     GetUsedHypothesis(aMesh, aShape, false);
@@ -138,11 +143,11 @@ bool StdMeshers_Quadrangle_2D::CheckHypothesis
   // First assigned hypothesis (if any) is processed now
   if (hyps.size() > 0) {
     aHyp = hyps.front();
-    if (strcmp("QuadrangleParams", aHyp->GetName()) == 0) {
-      const StdMeshers_QuadrangleParams* aHyp1 = 
-        (const StdMeshers_QuadrangleParams*)aHyp;
-      myTriaVertexID = aHyp1->GetTriaVertex();
-      myQuadType = aHyp1->GetQuadType();
+    if (strcmp("QuadrangleParams", aHyp->GetName()) == 0)
+    {
+      myParams = (const StdMeshers_QuadrangleParams*)aHyp;
+      myTriaVertexID = myParams->GetTriaVertex();
+      myQuadType     = myParams->GetQuadType();
       if (myQuadType == QUAD_QUADRANGLE_PREF ||
           myQuadType == QUAD_QUADRANGLE_PREF_REVERSED)
         myQuadranglePreference = true;
@@ -221,22 +226,25 @@ bool StdMeshers_Quadrangle_2D::Compute (SMESH_Mesh&         aMesh,
   FaceQuadStruct::Ptr quad = CheckNbEdges( aMesh, F, /*considerMesh=*/true );
   if (!quad)
     return false;
-  myQuadStruct = quad;
+  myQuadList.clear();
+  myQuadList.push_back( quad );
+
+  if ( !getEnforcedUV() )
+    return false;
 
   updateDegenUV( quad );
 
+  int n1 = quad->side[0].NbPoints();
+  int n2 = quad->side[1].NbPoints();
+  int n3 = quad->side[2].NbPoints();
+  int n4 = quad->side[3].NbPoints();
+
   enum { NOT_COMPUTED = -1, COMPUTE_FAILED = 0, COMPUTE_OK = 1 };
   int res = NOT_COMPUTED;
   if (myQuadranglePreference)
   {
-    int n1    = quad->side[0]->NbPoints();
-    int n2    = quad->side[1]->NbPoints();
-    int n3    = quad->side[2]->NbPoints();
-    int n4    = quad->side[3]->NbPoints();
     int nfull = n1+n2+n3+n4;
-    int ntmp  = nfull/2;
-    ntmp = ntmp*2;
-    if (nfull == ntmp && ((n1 != n3) || (n2 != n4)))
+    if ((nfull % 2) == 0 && ((n1 != n3) || (n2 != n4)))
     {
       // special path genarating only quandrangle faces
       res = computeQuadPref( aMesh, F, quad );
@@ -244,10 +252,6 @@ bool StdMeshers_Quadrangle_2D::Compute (SMESH_Mesh&         aMesh,
   }
   else if (myQuadType == QUAD_REDUCED)
   {
-    int n1     = quad->side[0]->NbPoints();
-    int n2     = quad->side[1]->NbPoints();
-    int n3     = quad->side[2]->NbPoints();
-    int n4     = quad->side[3]->NbPoints();
     int n13    = n1 - n3;
     int n24    = n2 - n4;
     int n13tmp = n13/2; n13tmp = n13tmp*2;
@@ -275,7 +279,10 @@ bool StdMeshers_Quadrangle_2D::Compute (SMESH_Mesh&         aMesh,
 
   if ( res == NOT_COMPUTED )
   {
-    res = computeQuadDominant( aMesh, F, quad );
+    if ( n1 != n3 || n2 != n4 )
+      res = computeTriangles( aMesh, F, quad );
+    else
+      res = computeQuadDominant( aMesh, F );
   }
 
   if ( res == COMPUTE_OK && myNeedSmooth )
@@ -286,6 +293,83 @@ bool StdMeshers_Quadrangle_2D::Compute (SMESH_Mesh&         aMesh,
 
 //================================================================================
 /*!
+ * \brief Compute quadrangles and triangles on the quad
+ */
+//================================================================================
+
+bool StdMeshers_Quadrangle_2D::computeTriangles(SMESH_Mesh&         aMesh,
+                                                const TopoDS_Face&  aFace,
+                                                FaceQuadStruct::Ptr quad)
+{
+  int nb = quad->side[0].grid->NbPoints();
+  int nr = quad->side[1].grid->NbPoints();
+  int nt = quad->side[2].grid->NbPoints();
+  int nl = quad->side[3].grid->NbPoints();
+
+  // rotate the quad to have nbNodeOut sides on TOP [and LEFT]
+  if ( nb > nt )
+    quad->shift( nl > nr ? 3 : 2, true );
+  else if ( nr > nl )
+    quad->shift( 1, true );
+  else if ( nl > nr )
+    quad->shift( nt > nb ? 0 : 3, true );
+
+  if ( !setNormalizedGrid( quad ))
+    return false;
+
+  if ( quad->nbNodeOut( QUAD_BOTTOM_SIDE ))
+  {
+    splitQuad( quad, 0, 1 );
+  }
+  if ( quad->nbNodeOut( QUAD_TOP_SIDE    ))
+  {
+    splitQuad( quad, 0, quad->jSize-2 );
+  }
+  FaceQuadStruct::Ptr newQuad = myQuadList.back();
+  if ( quad != newQuad ) // split done
+  {
+    // make quad be a greatest one
+    if ( quad->side[ QUAD_LEFT_SIDE ].NbPoints() == 2 ||
+         quad->side[ QUAD_RIGHT_SIDE ].NbPoints() == 2  )
+      quad = newQuad;
+    if ( !setNormalizedGrid( quad ))
+      return false;
+  }
+
+  if ( quad->nbNodeOut( QUAD_RIGHT_SIDE ))
+  {
+    splitQuad( quad, quad->iSize-2, 0 );
+  }
+  if ( quad->nbNodeOut( QUAD_LEFT_SIDE    ))
+  {
+    splitQuad( quad, 1, 0 );
+  }
+
+  return computeQuadDominant( aMesh, aFace );
+}
+
+//================================================================================
+/*!
+ * \brief Compute quadrangles and possibly triangles on all quads of myQuadList
+ */
+//================================================================================
+
+bool StdMeshers_Quadrangle_2D::computeQuadDominant(SMESH_Mesh&         aMesh,
+                                                   const TopoDS_Face&  aFace)
+{
+  if ( !addEnforcedNodes() )
+    return false;
+
+  std::list< FaceQuadStruct::Ptr >::iterator quad = myQuadList.begin();
+  for ( ; quad != myQuadList.end(); ++quad )
+    if ( !computeQuadDominant( aMesh, aFace, *quad ))
+      return false;
+
+  return true;
+}
+
+//================================================================================
+/*!
  * \brief Compute quadrangles and possibly triangles
  */
 //================================================================================
@@ -294,37 +378,28 @@ bool StdMeshers_Quadrangle_2D::computeQuadDominant(SMESH_Mesh&         aMesh,
                                                    const TopoDS_Face&  aFace,
                                                    FaceQuadStruct::Ptr quad)
 {
-  // set normalized grid on unit square in parametric domain
+  // --- set normalized grid on unit square in parametric domain
 
-  if (!setNormalizedGrid(aMesh, aFace, quad))
+  if ( !setNormalizedGrid( quad ))
     return false;
 
-  // --- compute 3D values on points, store points & quadrangles
-
-  int nbdown  = quad->side[0]->NbPoints();
-  int nbup    = quad->side[2]->NbPoints();
-
-  int nbright = quad->side[1]->NbPoints();
-  int nbleft  = quad->side[3]->NbPoints();
+  // --- create nodes on points, and create quadrangles
 
-  int nbhoriz  = Min(nbdown, nbup);
-  int nbvertic = Min(nbright, nbleft);
+  int nbhoriz  = quad->iSize;
+  int nbvertic = quad->jSize;
 
   // internal mesh nodes
   SMESHDS_Mesh *  meshDS = aMesh.GetMeshDS();
   Handle(Geom_Surface) S = BRep_Tool::Surface(aFace);
-  int i, j, geomFaceID = meshDS->ShapeToIndex(aFace);
-  for (i = 1; i < nbhoriz - 1; i++) {
-    for (j = 1; j < nbvertic - 1; j++) {
-      int ij = j * nbhoriz + i;
-      double u = quad->uv_grid[ij].u;
-      double v = quad->uv_grid[ij].v;
-      gp_Pnt P = S->Value(u, v);
-      SMDS_MeshNode * node = meshDS->AddNode(P.X(), P.Y(), P.Z());
-      meshDS->SetNodeOnFace(node, geomFaceID, u, v);
-      quad->uv_grid[ij].node = node;
+  int i,j,    geomFaceID = meshDS->ShapeToIndex(aFace);
+  for (i = 1; i < nbhoriz - 1; i++)
+    for (j = 1; j < nbvertic - 1; j++)
+    {
+      UVPtStruct& uvPnt = quad->UVPt( i, j );
+      gp_Pnt P          = S->Value( uvPnt.u, uvPnt.v );
+      uvPnt.node        = meshDS->AddNode(P.X(), P.Y(), P.Z());
+      meshDS->SetNodeOnFace( uvPnt.node, geomFaceID, uvPnt.u, uvPnt.v );
     }
-  }
   
   // mesh faces
 
@@ -340,21 +415,20 @@ bool StdMeshers_Quadrangle_2D::computeQuadDominant(SMESH_Mesh&         aMesh,
   //              i
   //             [0]
   
-  i = 0;
   int ilow = 0;
   int iup = nbhoriz - 1;
-  if (quad->isEdgeOut[3]) { ilow++; } else { if (quad->isEdgeOut[1]) iup--; }
+  if (quad->nbNodeOut(3)) { ilow++; } else { if (quad->nbNodeOut(1)) iup--; }
   
   int jlow = 0;
   int jup = nbvertic - 1;
-  if (quad->isEdgeOut[0]) { jlow++; } else { if (quad->isEdgeOut[2]) jup--; }
+  if (quad->nbNodeOut(0)) { jlow++; } else { if (quad->nbNodeOut(2)) jup--; }
   
   // regular quadrangles
   for (i = ilow; i < iup; i++) {
     for (j = jlow; j < jup; j++) {
       const SMDS_MeshNode *a, *b, *c, *d;
-      a = quad->uv_grid[      * nbhoriz + i    ].node;
-      b = quad->uv_grid[      * nbhoriz + i + 1].node;
+      a = quad->uv_grid[ j      * nbhoriz + i    ].node;
+      b = quad->uv_grid[ j      * nbhoriz + i + 1].node;
       c = quad->uv_grid[(j + 1) * nbhoriz + i + 1].node;
       d = quad->uv_grid[(j + 1) * nbhoriz + i    ].node;
       SMDS_MeshFace* face = myHelper->AddFace(a, b, c, d);
@@ -364,19 +438,25 @@ bool StdMeshers_Quadrangle_2D::computeQuadDominant(SMESH_Mesh&         aMesh,
     }
   }
 
-  const vector<UVPtStruct>& uv_e0 = quad->side[0]->GetUVPtStruct(true,0);
-  const vector<UVPtStruct>& uv_e1 = quad->side[1]->GetUVPtStruct(false,1);
-  const vector<UVPtStruct>& uv_e2 = quad->side[2]->GetUVPtStruct(true,1);
-  const vector<UVPtStruct>& uv_e3 = quad->side[3]->GetUVPtStruct(false,0);
+  // Boundary elements (must always be on an outer boundary of the FACE)
+  
+  const vector<UVPtStruct>& uv_e0 = quad->side[0].grid->GetUVPtStruct();
+  const vector<UVPtStruct>& uv_e1 = quad->side[1].grid->GetUVPtStruct();
+  const vector<UVPtStruct>& uv_e2 = quad->side[2].grid->GetUVPtStruct();
+  const vector<UVPtStruct>& uv_e3 = quad->side[3].grid->GetUVPtStruct();
 
   if (uv_e0.empty() || uv_e1.empty() || uv_e2.empty() || uv_e3.empty())
     return error(COMPERR_BAD_INPUT_MESH);
 
   double eps = Precision::Confusion();
 
-  // Boundary quadrangles
-  
-  if (quad->isEdgeOut[0]) {
+  int nbdown  = (int) uv_e0.size();
+  int nbup    = (int) uv_e2.size();
+  int nbright = (int) uv_e1.size();
+  int nbleft  = (int) uv_e3.size();
+
+  if (quad->nbNodeOut(0) && nbvertic == 2)
+  {
     // Down edge is out
     // 
     // |___|___|___|___|___|___|
@@ -393,8 +473,12 @@ bool StdMeshers_Quadrangle_2D::computeQuadDominant(SMESH_Mesh&         aMesh,
     // number of last node of the down edge to be processed
     int stop = nbdown - 1;
     // if right edge is out, we will stop at a node, previous to the last one
-    if (quad->isEdgeOut[1]) stop--;
-    
+    //if (quad->nbNodeOut(1)) stop--;
+    if ( quad->nbNodeOut( QUAD_RIGHT_SIDE ))
+      quad->UVPt( nbhoriz-1, 1 ).node = uv_e1[1].node;
+    if ( quad->nbNodeOut( QUAD_LEFT_SIDE ))
+      quad->UVPt( 0, 1 ).node = uv_e3[1].node;
+
     // for each node of the down edge find nearest node
     // in the first row of the regular grid and link them
     for (i = 0; i < stop; i++) {
@@ -449,7 +533,7 @@ bool StdMeshers_Quadrangle_2D::computeQuadDominant(SMESH_Mesh&         aMesh,
           if (face) meshDS->SetMeshElementOnShape(face, geomFaceID);
         }
         else {
-          splitQuad(meshDS, geomFaceID, a, b, c, d);
+          splitQuadFace(meshDS, geomFaceID, a, b, c, d);
         }
 
         // if node d is not at position g - make additional triangles
@@ -468,7 +552,8 @@ bool StdMeshers_Quadrangle_2D::computeQuadDominant(SMESH_Mesh&         aMesh,
       }
     }
   } else {
-    if (quad->isEdgeOut[2]) {
+    if (quad->nbNodeOut(2) && nbvertic == 2)
+    {
       // Up edge is out
       // 
       // <-<-<-<-<-<-<-<-<-<-<-<-< -- direction of processing
@@ -483,9 +568,16 @@ bool StdMeshers_Quadrangle_2D::computeQuadDominant(SMESH_Mesh&         aMesh,
 
       int g = nbhoriz - 1; // last processed node in the regular grid
 
+      ilow = 0;
+      iup = nbhoriz - 1;
+
       int stop = 0;
       // if left edge is out, we will stop at a second node
-      if (quad->isEdgeOut[3]) stop++;
+      //if (quad->nbNodeOut(3)) stop++;
+      if ( quad->nbNodeOut( QUAD_RIGHT_SIDE ))
+        quad->UVPt( nbhoriz-1, 0 ).node = uv_e1[ nbright-2 ].node;
+      if ( quad->nbNodeOut( QUAD_LEFT_SIDE ))
+        quad->UVPt( 0, 0 ).node = uv_e3[ nbleft-2 ].node;
 
       // for each node of the up edge find nearest node
       // in the first row of the regular grid and link them
@@ -536,10 +628,10 @@ bool StdMeshers_Quadrangle_2D::computeQuadDominant(SMESH_Mesh&         aMesh,
             if (face) meshDS->SetMeshElementOnShape(face, geomFaceID);
           }
           else {
-            splitQuad(meshDS, geomFaceID, a, b, c, d);
+            splitQuadFace(meshDS, geomFaceID, a, b, c, d);
           }
 
-          if (near + 1 < g) { // if d not is at g - make additional triangles
+          if (near + 1 < g) { // if d is not at g - make additional triangles
             for (int k = near + 1; k < g; k++) {
               c = quad->uv_grid[nbhoriz*(nbvertic - 2) + k].node;
               if (k + 1 > iup)
@@ -557,12 +649,14 @@ bool StdMeshers_Quadrangle_2D::computeQuadDominant(SMESH_Mesh&         aMesh,
   }
 
   // right or left boundary quadrangles
-  if (quad->isEdgeOut[1]) {
-//    MESSAGE("right edge is out");
+  if (quad->nbNodeOut( QUAD_RIGHT_SIDE ) && nbhoriz == 2)
+  {
     int g = 0; // last processed node in the grid
     int stop = nbright - 1;
-    if (quad->isEdgeOut[2]) stop--;
-    for (i = 0; i < stop; i++) {
+    i = 0;
+    if (quad->side[ QUAD_RIGHT_SIDE ].from != i    ) i++;
+    if (quad->side[ QUAD_RIGHT_SIDE ].to   != stop ) stop--;
+    for ( ; i < stop; i++) {
       const SMDS_MeshNode *a, *b, *c, *d;
       a = uv_e1[i].node;
       b = uv_e1[i + 1].node;
@@ -609,7 +703,7 @@ bool StdMeshers_Quadrangle_2D::computeQuadDominant(SMESH_Mesh&         aMesh,
           if (face) meshDS->SetMeshElementOnShape(face, geomFaceID);
         }
         else {
-          splitQuad(meshDS, geomFaceID, a, b, c, d);
+          splitQuadFace(meshDS, geomFaceID, a, b, c, d);
         }
 
         if (near - 1 > g) { // if d not is at g - make additional triangles
@@ -627,12 +721,14 @@ bool StdMeshers_Quadrangle_2D::computeQuadDominant(SMESH_Mesh&         aMesh,
       }
     }
   } else {
-    if (quad->isEdgeOut[3]) {
+    if (quad->nbNodeOut(3) && nbhoriz == 2) {
 //      MESSAGE("left edge is out");
       int g = nbvertic - 1; // last processed node in the grid
       int stop = 0;
-      if (quad->isEdgeOut[0]) stop++;
-      for (i = nbleft - 1; i > stop; i--) {
+      i = nbleft - 1;
+      if (quad->side[3].from != stop ) stop++;
+      if (quad->side[3].to   != i    ) i--;
+      for (; i > stop; i--) {
         const SMDS_MeshNode *a, *b, *c, *d;
         a = uv_e3[i].node;
         b = uv_e3[i - 1].node;
@@ -678,7 +774,7 @@ bool StdMeshers_Quadrangle_2D::computeQuadDominant(SMESH_Mesh&         aMesh,
             if (face) meshDS->SetMeshElementOnShape(face, geomFaceID);
           }
           else {
-            splitQuad(meshDS, geomFaceID, a, b, c, d);
+            splitQuadFace(meshDS, geomFaceID, a, b, c, d);
           }
 
           if (near + 1 < g) { // if d not is at g - make additional triangles
@@ -828,8 +924,8 @@ FaceQuadStruct::Ptr StdMeshers_Quadrangle_2D::CheckNbEdges(SMESH_Mesh &
                                                            const TopoDS_Shape & aShape,
                                                            const bool           considerMesh)
 {
-  if ( myQuadStruct && myQuadStruct->face.IsSame( aShape ))
-    return myQuadStruct;
+  if ( !myQuadList.empty() && myQuadList.front()->face.IsSame( aShape ))
+    return myQuadList.front();
 
   TopoDS_Face F = TopoDS::Face(aShape);
   if ( F.Orientation() >= TopAbs_INTERNAL ) F.Orientation( TopAbs_FORWARD );
@@ -852,7 +948,6 @@ FaceQuadStruct::Ptr StdMeshers_Quadrangle_2D::CheckNbEdges(SMESH_Mesh &
     return FaceQuadStruct::Ptr();
   }
   FaceQuadStruct::Ptr quad( new FaceQuadStruct );
-  quad->uv_grid = 0;
   quad->side.reserve(nbEdgesInWire.front());
   quad->face = F;
 
@@ -870,17 +965,17 @@ FaceQuadStruct::Ptr StdMeshers_Quadrangle_2D::CheckNbEdges(SMESH_Mesh &
         else
           sideEdges.push_back( *edgeIt++ );
       if ( !sideEdges.empty() )
-        quad->side.push_back(new StdMeshers_FaceSide(F, sideEdges, &aMesh, iSide < QUAD_TOP_SIDE,
-                                                     ignoreMediumNodes, myProxyMesh));
+        quad->side.push_back( StdMeshers_FaceSide::New(F, sideEdges, &aMesh, iSide < QUAD_TOP_SIDE,
+                                                       ignoreMediumNodes, myProxyMesh));
       else
         --iSide;
     }
-    const vector<UVPtStruct>& UVPSleft  = quad->side[0]->GetUVPtStruct(true,0);
-    /*  vector<UVPtStruct>& UVPStop   = */quad->side[1]->GetUVPtStruct(false,1);
-    /*  vector<UVPtStruct>& UVPSright = */quad->side[2]->GetUVPtStruct(true,1);
+    const vector<UVPtStruct>& UVPSleft  = quad->side[0].GetUVPtStruct(true,0);
+    /*  vector<UVPtStruct>& UVPStop   = */quad->side[1].GetUVPtStruct(false,1);
+    /*  vector<UVPtStruct>& UVPSright = */quad->side[2].GetUVPtStruct(true,1);
     const SMDS_MeshNode* aNode = UVPSleft[0].node;
-    gp_Pnt2d aPnt2d(UVPSleft[0].u, UVPSleft[0].v);
-    quad->side.push_back(new StdMeshers_FaceSide(quad->side[1], aNode, &aPnt2d));
+    gp_Pnt2d aPnt2d = UVPSleft[0].UV();
+    quad->side.push_back( StdMeshers_FaceSide::New( quad->side[1].grid.get(), aNode, &aPnt2d ));
     myNeedSmooth = ( nbDegenEdges > 0 );
     return quad;
   }
@@ -922,15 +1017,15 @@ FaceQuadStruct::Ptr StdMeshers_Quadrangle_2D::CheckNbEdges(SMESH_Mesh &
       }
       if ( !sideEdges.empty() )
       {
-        quad->side.push_back(new StdMeshers_FaceSide(F, sideEdges, &aMesh, iSide < QUAD_TOP_SIDE,
-                                                     ignoreMediumNodes, myProxyMesh));
+        quad->side.push_back( StdMeshers_FaceSide::New( F, sideEdges, &aMesh, iSide < QUAD_TOP_SIDE,
+                                                        ignoreMediumNodes, myProxyMesh ));
         ++iSide;
       }
       else if ( !SMESH_Algo::isDegenerated( *edgeIt ) && // closed EDGE
                 myHelper->IthVertex( 0, *edgeIt ).IsSame( myHelper->IthVertex( 1, *edgeIt )))
       {
-        quad->side.push_back(new StdMeshers_FaceSide(F, *edgeIt++, &aMesh, iSide < QUAD_TOP_SIDE,
-                                                     ignoreMediumNodes, myProxyMesh));
+        quad->side.push_back( StdMeshers_FaceSide::New( F, *edgeIt++, &aMesh, iSide < QUAD_TOP_SIDE,
+                                                        ignoreMediumNodes, myProxyMesh));
         ++iSide;
       }
       if ( quad->side.size() == 4 )
@@ -1155,36 +1250,12 @@ StdMeshers_Quadrangle_2D::CheckAnd2Dcompute (SMESH_Mesh &         aMesh,
   if ( quad )
   {
     // set normalized grid on unit square in parametric domain
-    if ( ! setNormalizedGrid( aMesh, TopoDS::Face( aShape ), quad))
+    if ( ! setNormalizedGrid( quad ))
       quad.reset();
   }
   return quad;
 }
 
-//=============================================================================
-/*!
- *
- */
-//=============================================================================
-
-faceQuadStruct::~faceQuadStruct()
-{
-  for (size_t i = 0; i < side.size(); i++) {
-    if (side[i]) {
-      delete side[i];
-      for (size_t j = i+1; j < side.size(); j++)
-        if ( side[i] == side[j] )
-          side[j] = 0;
-    }
-  }
-  side.clear();
-
-  if (uv_grid) {
-    delete [] uv_grid;
-    uv_grid = 0;
-  }
-}
-
 namespace
 {
   inline const vector<UVPtStruct>& getUVPtStructIn(FaceQuadStruct::Ptr& quad, int i, int nbSeg)
@@ -1192,9 +1263,9 @@ namespace
     bool   isXConst   = (i == QUAD_BOTTOM_SIDE || i == QUAD_TOP_SIDE);
     double constValue = (i == QUAD_BOTTOM_SIDE || i == QUAD_LEFT_SIDE) ? 0 : 1;
     return
-      quad->isEdgeOut[i] ?
-      quad->side[i]->SimulateUVPtStruct(nbSeg,isXConst,constValue) :
-      quad->side[i]->GetUVPtStruct(isXConst,constValue);
+      quad->nbNodeOut(i) ?
+      quad->side[i].grid->SimulateUVPtStruct(nbSeg,isXConst,constValue) :
+      quad->side[i].grid->GetUVPtStruct     (isXConst,constValue);
   }
   inline gp_UV calcUV(double x, double y,
                       const gp_UV& a0,const gp_UV& a1,const gp_UV& a2,const gp_UV& a3,
@@ -1212,10 +1283,11 @@ namespace
  */
 //=============================================================================
 
-bool StdMeshers_Quadrangle_2D::setNormalizedGrid (SMESH_Mesh &          aMesh,
-                                                  const TopoDS_Face&    aFace,
-                                                  FaceQuadStruct::Ptr & quad)
+bool StdMeshers_Quadrangle_2D::setNormalizedGrid (FaceQuadStruct::Ptr quad)
 {
+  if ( !quad->uv_grid.empty() )
+    return true;
+
   // Algorithme décrit dans "Génération automatique de maillages"
   // P.L. GEORGE, MASSON, § 6.4.1 p. 84-85
   // traitement dans le domaine paramétrique 2d u,v
@@ -1233,81 +1305,133 @@ bool StdMeshers_Quadrangle_2D::setNormalizedGrid (SMESH_Mesh &          aMesh,
   //             =down
   //
 
-  int nbhoriz  = Min(quad->side[0]->NbPoints(), quad->side[2]->NbPoints());
-  int nbvertic = Min(quad->side[1]->NbPoints(), quad->side[3]->NbPoints());
-
-  quad->isEdgeOut[0] = (quad->side[0]->NbPoints() > quad->side[2]->NbPoints());
-  quad->isEdgeOut[1] = (quad->side[1]->NbPoints() > quad->side[3]->NbPoints());
-  quad->isEdgeOut[2] = (quad->side[2]->NbPoints() > quad->side[0]->NbPoints());
-  quad->isEdgeOut[3] = (quad->side[3]->NbPoints() > quad->side[1]->NbPoints());
-
-  UVPtStruct *uv_grid = quad->uv_grid = new UVPtStruct[nbvertic * nbhoriz];
+  int nbhoriz  = Min(quad->side[0].NbPoints(), quad->side[2].NbPoints());
+  int nbvertic = Min(quad->side[1].NbPoints(), quad->side[3].NbPoints());
 
-  const vector<UVPtStruct>& uv_e0 = getUVPtStructIn(quad, 0, nbhoriz - 1);
-  const vector<UVPtStruct>& uv_e1 = getUVPtStructIn(quad, 1, nbvertic - 1);
-  const vector<UVPtStruct>& uv_e2 = getUVPtStructIn(quad, 2, nbhoriz - 1);
-  const vector<UVPtStruct>& uv_e3 = getUVPtStructIn(quad, 3, nbvertic - 1);
+  if ( myQuadList.size() == 1 )
+  {
+    // all sub-quads must have NO sides with nbNodeOut > 0
+    quad->nbNodeOut(0) = Max( 0, quad->side[0].grid->NbPoints() - quad->side[2].grid->NbPoints());
+    quad->nbNodeOut(1) = Max( 0, quad->side[1].grid->NbPoints() - quad->side[3].grid->NbPoints());
+    quad->nbNodeOut(2) = Max( 0, quad->side[2].grid->NbPoints() - quad->side[0].grid->NbPoints());
+    quad->nbNodeOut(3) = Max( 0, quad->side[3].grid->NbPoints() - quad->side[1].grid->NbPoints());
+  }
+  int from[4] = {
+    quad->side[0].from,
+    quad->side[1].from,
+    quad->side[2].from,
+    quad->side[3].from
+  };
+  const vector<UVPtStruct>& uv_e0_vec = quad->side[ 0 ].GetUVPtStruct();
+  const vector<UVPtStruct>& uv_e1_vec = quad->side[ 1 ].GetUVPtStruct();
+  const vector<UVPtStruct>& uv_e2_vec = quad->side[ 2 ].GetUVPtStruct();
+  const vector<UVPtStruct>& uv_e3_vec = quad->side[ 3 ].GetUVPtStruct();
 
-  if (uv_e0.empty() || uv_e1.empty() || uv_e2.empty() || uv_e3.empty())
+  if (uv_e0_vec.empty() || uv_e1_vec.empty() || uv_e2_vec.empty() || uv_e3_vec.empty())
     //return error("Can't find nodes on sides");
     return error(COMPERR_BAD_INPUT_MESH);
 
+  UVPtStruct* uv_e0 = (UVPtStruct*) & uv_e0_vec[0] + from[0];
+  UVPtStruct* uv_e1 = (UVPtStruct*) & uv_e1_vec[0] + from[1];
+  UVPtStruct* uv_e2 = (UVPtStruct*) & uv_e2_vec[0] + from[2];
+  UVPtStruct* uv_e3 = (UVPtStruct*) & uv_e3_vec[0] + from[3];
+
+  quad->uv_grid.resize( nbvertic * nbhoriz );
+  quad->iSize = nbhoriz;
+  quad->jSize = nbvertic;
+  UVPtStruct *uv_grid = & quad->uv_grid[0];
+
+  quad->uv_box.Clear();
+
   // copy data of face boundary
-  {
+
+  { // BOTTOM
     const int j = 0;
-    for (int i = 0; i < nbhoriz; i++)       // down
+    const double x0 = uv_e0[ 0 ].normParam;
+    const double dx = uv_e0[ nbhoriz-1 ].normParam - uv_e0[ 0 ].normParam;
+    for (int i = 0; i < nbhoriz; i++) {      // down
+      uv_e0[i].x = ( uv_e0[i].normParam - x0 ) / dx;
+      uv_e0[i].y = 0.;
       uv_grid[ j * nbhoriz + i ] = uv_e0[i];
+      quad->uv_box.Add( uv_e0[i].UV() );
+    }
   }
-  {
+  { // RIGHT
     const int i = nbhoriz - 1;
-    for (int j = 0; j < nbvertic; j++)      // right
+    const double y0 = uv_e1[ 0 ].normParam;
+    const double dy = uv_e1[ nbvertic-1 ].normParam - uv_e1[ 0 ].normParam;
+    int j = 0, nb = nbvertic;
+    if ( quad->UVPt( i, j ).node ) ++j; // avoid copying from a split emulated side
+    for ( ; j < nb; j++) {      // right
+      uv_e1[j].x = 1.;
+      uv_e1[j].y = ( uv_e1[j].normParam - y0 ) / dy;
       uv_grid[ j * nbhoriz + i ] = uv_e1[j];
+      quad->uv_box.Add( uv_e1[j].UV() );
+    }
   }
-  {
+  { // TOP
     const int j = nbvertic - 1;
-    for (int i = 0; i < nbhoriz; i++)       // up
+    const double x0 = uv_e2[ 0 ].normParam;
+    const double dx = uv_e2[ nbhoriz-1 ].normParam - uv_e2[ 0 ].normParam;
+    int i = 0, nb = nbhoriz;
+    if ( quad->UVPt( nb-1, j ).node ) --nb; // avoid copying from a split emulated side
+    for (; i < nb; i++) {      // up
+      uv_e2[i].x = ( uv_e2[i].normParam - x0 ) / dx;
+      uv_e2[i].y = 1.;
       uv_grid[ j * nbhoriz + i ] = uv_e2[i];
+      quad->uv_box.Add( uv_e2[i].UV() );
+    }
   }
-  {
+  { // LEFT
     const int i = 0;
-    for (int j = 0; j < nbvertic; j++)      // left
+    const double y0 = uv_e3[ 0 ].normParam;
+    const double dy = uv_e3[ nbvertic-1 ].normParam - uv_e3[ 0 ].normParam;
+    int j = 0, nb = nbvertic;
+    if ( quad->UVPt( i, j    ).node ) ++j; // avoid copying from a split emulated side
+    if ( quad->UVPt( i, nb-1 ).node ) --nb;
+    for ( ; j < nb; j++) {     // left
+      uv_e3[j].x = 0.;
+      uv_e3[j].y = ( uv_e3[j].normParam - y0 ) / dy;
       uv_grid[ j * nbhoriz + i ] = uv_e3[j];
+      quad->uv_box.Add( uv_e3[j].UV() );
+    }
   }
 
   // normalized 2d parameters on grid
 
-  for (int i = 0; i < nbhoriz; i++) {
-    for (int j = 0; j < nbvertic; j++) {
-      int ij = j * nbhoriz + i;
-      // --- droite i cste : x = x0 + y(x1-x0)
-      double x0 = uv_e0[i].normParam;   // bas  - sud
-      double x1 = uv_e2[i].normParam;   // haut - nord
-      // --- droite j cste : y = y0 + x(y1-y0)
-      double y0 = uv_e3[j].normParam;   // gauche - ouest
-      double y1 = uv_e1[j].normParam;   // droite - est
+  for (int i = 1; i < nbhoriz-1; i++)
+  {
+    const double x0 = uv_e0[i].x;
+    const double x1 = uv_e2[i].x;
+    for (int j = 1; j < nbvertic-1; j++)
+    {
+      const double y0 = uv_e3[j].y;
+      const double y1 = uv_e1[j].y;
       // --- intersection : x=x0+(y0+x(y1-y0))(x1-x0)
       double x = (x0 + y0 * (x1 - x0)) / (1 - (y1 - y0) * (x1 - x0));
       double y = y0 + x * (y1 - y0);
+      int   ij = j * nbhoriz + i;
       uv_grid[ij].x = x;
       uv_grid[ij].y = y;
+      uv_grid[ij].node = NULL;
     }
   }
 
   // projection on 2d domain (u,v)
 
-  gp_UV a0 (uv_e0.front().u, uv_e0.front().v);
-  gp_UV a1 (uv_e0.back().u,  uv_e0.back().v );
-  gp_UV a2 (uv_e2.back().u,  uv_e2.back().v );
-  gp_UV a3 (uv_e2.front().u, uv_e2.front().v);
+  gp_UV a0 = uv_e0[0        ].UV();
+  gp_UV a1 = uv_e0[nbhoriz-1].UV();
+  gp_UV a2 = uv_e2[nbhoriz-1].UV();
+  gp_UV a3 = uv_e2[0        ].UV();
 
-  for (int i = 0; i < nbhoriz; i++)
+  for (int i = 1; i < nbhoriz-1; i++)
   {
-    gp_UV p0( uv_e0[i].u, uv_e0[i].v );
-    gp_UV p2( uv_e2[i].u, uv_e2[i].v );
-    for (int j = 0; j < nbvertic; j++)
+    gp_UV p0 = uv_e0[i].UV();
+    gp_UV p2 = uv_e2[i].UV();
+    for (int j = 1; j < nbvertic-1; j++)
     {
-      gp_UV p1( uv_e1[j].u, uv_e1[j].v );
-      gp_UV p3( uv_e3[j].u, uv_e3[j].v );
+      gp_UV p1 = uv_e1[j].UV();
+      gp_UV p3 = uv_e3[j].UV();
 
       int ij = j * nbhoriz + i;
       double x = uv_grid[ij].x;
@@ -1343,7 +1467,7 @@ static void shiftQuad(FaceQuadStruct::Ptr& quad, const int num)
 void FaceQuadStruct::shift( size_t nb, bool ori )
 {
   if ( nb == 0 ) return;
-  StdMeshers_FaceSide* sideArr[4] = { side[0], side[1], side[2], side[3] };
+  StdMeshers_FaceSidePtr sideArr[4] = { side[0], side[1], side[2], side[3] };
   for (int i = QUAD_BOTTOM_SIDE; i < NB_QUAD_SIDES; ++i) {
     int id = (i + nb) % NB_QUAD_SIDES;
     bool wasForward = (i  < QUAD_TOP_SIDE);
@@ -1367,10 +1491,10 @@ static gp_UV calcUV(double x0, double x1, double y0, double y1,
   double x = (x0 + y0 * (x1 - x0)) / (1 - (y1 - y0) * (x1 - x0));
   double y = y0 + x * (y1 - y0);
 
-  gp_UV p0 = quad->side[QUAD_BOTTOM_SIDE]->Value2d(x).XY();
-  gp_UV p1 = quad->side[QUAD_RIGHT_SIDE ]->Value2d(y).XY();
-  gp_UV p2 = quad->side[QUAD_TOP_SIDE   ]->Value2d(x).XY();
-  gp_UV p3 = quad->side[QUAD_LEFT_SIDE  ]->Value2d(y).XY();
+  gp_UV p0 = quad->side[QUAD_BOTTOM_SIDE].grid->Value2d(x).XY();
+  gp_UV p1 = quad->side[QUAD_RIGHT_SIDE ].grid->Value2d(y).XY();
+  gp_UV p2 = quad->side[QUAD_TOP_SIDE   ].grid->Value2d(x).XY();
+  gp_UV p3 = quad->side[QUAD_LEFT_SIDE  ].grid->Value2d(y).XY();
 
   gp_UV uv = calcUV(x,y, a0,a1,a2,a3, p0,p1,p2,p3);
 
@@ -1387,10 +1511,10 @@ static gp_UV calcUV2(double x, double y,
                      const gp_UV& a0, const gp_UV& a1,
                      const gp_UV& a2, const gp_UV& a3)
 {
-  gp_UV p0 = quad->side[QUAD_BOTTOM_SIDE]->Value2d(x).XY();
-  gp_UV p1 = quad->side[QUAD_RIGHT_SIDE ]->Value2d(y).XY();
-  gp_UV p2 = quad->side[QUAD_TOP_SIDE   ]->Value2d(x).XY();
-  gp_UV p3 = quad->side[QUAD_LEFT_SIDE  ]->Value2d(y).XY();
+  gp_UV p0 = quad->side[QUAD_BOTTOM_SIDE].grid->Value2d(x).XY();
+  gp_UV p1 = quad->side[QUAD_RIGHT_SIDE ].grid->Value2d(y).XY();
+  gp_UV p2 = quad->side[QUAD_TOP_SIDE   ].grid->Value2d(x).XY();
+  gp_UV p3 = quad->side[QUAD_LEFT_SIDE  ].grid->Value2d(y).XY();
 
   gp_UV uv = calcUV(x,y, a0,a1,a2,a3, p0,p1,p2,p3);
 
@@ -1418,24 +1542,37 @@ bool StdMeshers_Quadrangle_2D::computeQuadPref (SMESH_Mesh &        aMesh,
   bool WisF = true;
   int i,j,geomFaceID = meshDS->ShapeToIndex(aFace);
 
-  int nb = quad->side[0]->NbPoints();
-  int nr = quad->side[1]->NbPoints();
-  int nt = quad->side[2]->NbPoints();
-  int nl = quad->side[3]->NbPoints();
+  int nb = quad->side[0].grid->NbPoints();
+  int nr = quad->side[1].grid->NbPoints();
+  int nt = quad->side[2].grid->NbPoints();
+  int nl = quad->side[3].grid->NbPoints();
   int dh = abs(nb-nt);
   int dv = abs(nr-nl);
 
-  // rotate sides to be as in the picture below and to have
-  // dh >= dv and nt > nb
-  if ( dh >= dv )
-    shiftQuad( quad, ( nt > nb ) ? 0 : 2 );
+  if ( myForcedPnts.empty() )
+  {
+    // rotate sides to be as in the picture below and to have
+    // dh >= dv and nt > nb
+    if ( dh >= dv )
+      shiftQuad( quad, ( nt > nb ) ? 0 : 2 );
+    else
+      shiftQuad( quad, ( nr > nl ) ? 1 : 3 );
+  }
   else
-    shiftQuad( quad, ( nr > nl ) ? 1 : 3 );
+  {
+    // rotate the quad to have nt > nb [and nr > nl]
+    if ( nb > nt )
+      quad->shift( nr > nl ? 1 : 2, true );
+    else if ( nr > nl )
+      quad->shift( nb == nt ? 1 : 0, true );
+    else if ( nl > nr )
+      quad->shift( 3, true );
+  }
 
-  nb = quad->side[0]->NbPoints();
-  nr = quad->side[1]->NbPoints();
-  nt = quad->side[2]->NbPoints();
-  nl = quad->side[3]->NbPoints();
+  nb = quad->side[0].grid->NbPoints();
+  nr = quad->side[1].grid->NbPoints();
+  nt = quad->side[2].grid->NbPoints();
+  nl = quad->side[3].grid->NbPoints();
   dh = abs(nb-nt);
   dv = abs(nr-nl);
   int nbh  = Max(nb,nt);
@@ -1470,6 +1607,302 @@ bool StdMeshers_Quadrangle_2D::computeQuadPref (SMESH_Mesh &        aMesh,
   //      0------------0
   //       0  bottom  1
 
+  const vector<UVPtStruct>& uv_eb = quad->side[0].GetUVPtStruct(true,0);
+  const vector<UVPtStruct>& uv_er = quad->side[1].GetUVPtStruct(false,1);
+  const vector<UVPtStruct>& uv_et = quad->side[2].GetUVPtStruct(true,1);
+  const vector<UVPtStruct>& uv_el = quad->side[3].GetUVPtStruct(false,0);
+
+  if (uv_eb.size() != nb || uv_er.size() != nr || uv_et.size() != nt || uv_el.size() != nl)
+    return error(COMPERR_BAD_INPUT_MESH);
+
+  gp_UV a0,a1,a2,a3, p0,p1,p2,p3, uv;
+  double x,y;
+
+  a0 = uv_eb[ 0 ].UV();
+  a1 = uv_er[ 0 ].UV();
+  a2 = uv_er[ nr-1 ].UV();
+  a3 = uv_et[ 0 ].UV();
+
+  if ( !myForcedPnts.empty() )
+  {
+    if ( dv != 0 && dh != 0 )
+    {
+      const int dmin = Min( dv, dh );
+
+      // Make a side separating domains L and Cb
+      StdMeshers_FaceSidePtr sideLCb;
+      UVPtStruct p3dom; // a point where 3 domains meat
+      {                                                     //   dmin
+        vector<UVPtStruct> pointsLCb( dmin+1 );             // 1--------1
+        pointsLCb[0] = uv_eb[0];                            //  |   |  |
+        for ( int i = 1; i <= dmin; ++i )                   //  |   |Ct|
+        {                                                   //  | L |  |
+          x  = uv_et[ i ].normParam;                        //  |   |__|
+          y  = uv_er[ i ].normParam;                        //  |  /   |
+          p0 = quad->side[0].grid->Value2d( x ).XY();       //  | / Cb |dmin
+          p1 = uv_er[ i ].UV();                             //  |/     |
+          p2 = uv_et[ i ].UV();                             // 0--------0
+          p3 = quad->side[3].grid->Value2d( y ).XY();
+          uv = calcUV( x,y, a0,a1,a2,a3, p0,p1,p2,p3 );
+          pointsLCb[ i ].u = uv.X();
+          pointsLCb[ i ].v = uv.Y();
+        }
+        sideLCb = StdMeshers_FaceSide::New( pointsLCb, aFace );
+        p3dom   = pointsLCb.back();
+      }
+      // Make a side separating domains L and Ct
+      StdMeshers_FaceSidePtr sideLCt;
+      {
+        vector<UVPtStruct> pointsLCt( nl );
+        pointsLCt[0]     = p3dom;
+        pointsLCt.back() = uv_et[ dmin ];
+        x  = uv_et[ dmin ].normParam;
+        p0 = quad->side[0].grid->Value2d( x ).XY();
+        p2 = uv_et[ dmin ].UV();
+        for ( int i = 1; i < nl; ++i )
+        {
+          y  = uv_er[ i + dmin ].normParam;
+          p1 = uv_er[ i + dmin ].UV();
+          p3 = quad->side[3].grid->Value2d( y ).XY();
+          uv = calcUV( x,y, a0,a1,a2,a3, p0,p1,p2,p3 );
+          pointsLCt[ i ].u = uv.X();
+          pointsLCt[ i ].v = uv.Y();
+        }
+        sideLCt = StdMeshers_FaceSide::New( pointsLCt, aFace );
+      }
+      // Make a side separating domains Cb and Ct
+      StdMeshers_FaceSidePtr sideCbCt;
+      {
+        vector<UVPtStruct> pointsCbCt( nb );
+        pointsCbCt[0]     = p3dom;
+        pointsCbCt.back() = uv_er[ dmin ];
+        y  = uv_er[ dmin ].normParam;
+        p1 = uv_er[ dmin ].UV();
+        p3 = quad->side[3].grid->Value2d( y ).XY();
+        for ( int i = 1; i < nb-1; ++i )
+        {
+          x  = uv_et[ i + dmin ].normParam;
+          p2 = uv_et[ i + dmin ].UV();
+          p0 = quad->side[0].grid->Value2d( x ).XY();
+          uv = calcUV( x,y, a0,a1,a2,a3, p0,p1,p2,p3 );
+          pointsCbCt[ i ].u = uv.X();
+          pointsCbCt[ i ].v = uv.Y();
+        }
+        sideCbCt = StdMeshers_FaceSide::New( pointsCbCt, aFace );
+      }
+      // Make Cb quad
+      FaceQuadStruct* qCb = new FaceQuadStruct( quad->face );
+      myQuadList.push_back( FaceQuadStruct::Ptr( qCb ));
+      qCb->side.resize(4);
+      qCb->side[0] = quad->side[0];
+      qCb->side[1] = quad->side[1];
+      qCb->side[2] = sideCbCt;
+      qCb->side[3] = sideLCb;
+      qCb->side[1].to = dmin+1;
+      // Make L quad
+      FaceQuadStruct* qL = new FaceQuadStruct( quad->face );
+      myQuadList.push_back( FaceQuadStruct::Ptr( qL ));
+      qL->side.resize(4);
+      qL->side[0] = sideLCb;
+      qL->side[1] = sideLCt;
+      qL->side[2] = quad->side[2];
+      qL->side[3] = quad->side[3];
+      qL->side[2].to = dmin+1;
+      // Make Ct from the main quad
+      FaceQuadStruct::Ptr qCt = quad;
+      qCt->side[0] = sideCbCt;
+      qCt->side[3] = sideLCt;
+      qCt->side[1].from = dmin;
+      qCt->side[2].from = dmin;
+      qCt->uv_grid.clear();
+
+      // Connect sides
+      qCb->side[3].AddContact( dmin, & qCb->side[2], 0 );
+      qCb->side[3].AddContact( dmin, & qCt->side[3], 0 );
+      qCt->side[3].AddContact(    0, & qCt->side[0], 0 );
+      qCt->side[0].AddContact(    0, & qL ->side[0], dmin );
+      qL ->side[0].AddContact( dmin, & qL ->side[1], 0 );
+      qL ->side[0].AddContact( dmin, & qCb->side[2], 0 );
+
+      if ( dh == dv )
+        return computeQuadDominant( aMesh, aFace );
+      else
+        return computeQuadPref( aMesh, aFace, qCt );
+
+    } // if ( dv != 0 && dh != 0 )
+
+    // Case dv == 0
+    //
+    //     lw   nb  lw = dh/2
+    //    +------------+
+    //    |   |    |   |
+    //    |   | Ct |   |
+    //    | L |    | R |
+    //    |   |____|   |
+    //    |  /      \  |
+    //    | /   Cb   \ |
+    //    |/          \|
+    //    +------------+
+    const int lw = dh/2; // lateral width
+    const int bfrom = quad->side[0].from;
+    const int rfrom = quad->side[1].from;
+    const int tfrom = quad->side[2].from;
+    const int lfrom = quad->side[3].from;
+
+    const double   lL = quad->side[3].Length();
+    const double lLwL = quad->side[2].Length( tfrom, tfrom + lw + 1 );
+    const double yCbL = lLwL / ( lLwL + lL );
+
+    const double   lR = quad->side[1].Length();
+    const double lLwR = quad->side[2].Length( nt - lw - 1, nt );
+    const double yCbR = lLwR / ( lLwR + lR );
+
+    // Make sides separating domains Cb and L and R
+    StdMeshers_FaceSidePtr sideLCb, sideRCb;
+    UVPtStruct pTBL, pTBR; // points where 3 domains meat
+    {
+      vector<UVPtStruct> pointsLCb( lw+1 ), pointsRCb( lw+1 );
+      pointsLCb[0] = uv_eb[ 0  + bfrom ];
+      pointsRCb[0] = uv_eb[ nb + bfrom ];
+      for ( int i = 1, i2 = nt-2; i <= lw; ++i, --i2 )
+      {
+        x  = quad->side[2].Param( i );
+        y  = yCbL * i / lw;
+        p0 = quad->side[0].Value2d( x );
+        p1 = quad->side[1].Value2d( y );
+        p2 = uv_et[ i + tfrom ].UV();
+        p3 = quad->side[3].Value2d( y );
+        uv = calcUV( x,y, a0,a1,a2,a3, p0,p1,p2,p3 );
+        pointsLCb[ i ].u = uv.X();
+        pointsLCb[ i ].v = uv.Y();
+        pointsLCb[ i ].x = x;
+
+        x  = quad->side[2].Param( i2 );
+        y  = yCbR * i / lw;
+        p0 = quad->side[0].Value2d( x );
+        p2 = uv_et[ i2 + tfrom ].UV();
+        uv = calcUV( x,y, a0,a1,a2,a3, p0,p1,p2,p3 );
+        pointsRCb[ i ].u = uv.X();
+        pointsRCb[ i ].v = uv.Y();
+        pointsRCb[ i ].x = x;
+      }
+      sideLCb = StdMeshers_FaceSide::New( pointsLCb, aFace );
+      sideRCb = StdMeshers_FaceSide::New( pointsRCb, aFace );
+      pTBL    = pointsLCb.back();
+      pTBR    = pointsRCb.back();
+    }
+    // Make sides separating domains Ct and L and R
+    StdMeshers_FaceSidePtr sideLCt, sideRCt;
+    {
+      vector<UVPtStruct> pointsLCt( nl ), pointsRCt( nl );
+      pointsLCt[0]     = pTBL;
+      pointsLCt.back() = uv_et[ lw + tfrom ];
+      pointsRCt[0]     = pTBR;
+      pointsRCt.back() = uv_et[ lw + nb - 1 + tfrom ];
+      x  = pTBL.x;
+      p0 = quad->side[0].Value2d( x );
+      p2 = uv_et[ lw + tfrom ].UV();
+      int     iR = lw + nb - 1;
+      double  xR = pTBR.x;
+      gp_UV  p0R = quad->side[0].Value2d( xR );
+      gp_UV  p2R = uv_et[ iR + tfrom ].UV();
+      for ( int i = 1; i < nl; ++i )
+      {
+        y  = yCbL + ( 1. - yCbL ) * i / nl;
+        p1 = quad->side[1].Value2d( y );
+        p3 = quad->side[3].Value2d( y );
+        uv = calcUV( x,y, a0,a1,a2,a3, p0,p1,p2,p3 );
+        pointsLCt[ i ].u = uv.X();
+        pointsLCt[ i ].v = uv.Y();
+
+        y  = yCbR + ( 1. - yCbR ) * i / nl;
+        p1 = quad->side[1].Value2d( y );
+        p3 = quad->side[3].Value2d( y );
+        uv = calcUV( xR,y, a0,a1,a2,a3, p0,p1,p2,p3 );
+        pointsRCt[ i ].u = uv.X();
+        pointsRCt[ i ].v = uv.Y();
+      }
+      sideLCt = StdMeshers_FaceSide::New( pointsLCt, aFace );
+      sideRCt = StdMeshers_FaceSide::New( pointsRCt, aFace );
+    }
+    // Make a side separating domains Cb and Ct
+    StdMeshers_FaceSidePtr sideCbCt;
+    {
+      vector<UVPtStruct> pointsCbCt( nb );
+      pointsCbCt[0]     = pTBL;
+      pointsCbCt.back() = pTBR;
+      p1 = quad->side[1].Value2d( yCbR );
+      p3 = quad->side[3].Value2d( yCbL );
+      for ( int i = 1; i < nb-1; ++i )
+      {
+        x  = quad->side[2].Param( i + lw );
+        y  = yCbL + ( yCbR - yCbL ) * i / nb;
+        p2 = uv_et[ i + lw + tfrom ].UV();
+        p0 = quad->side[0].Value2d( x );
+        uv = calcUV( x,y, a0,a1,a2,a3, p0,p1,p2,p3 );
+        pointsCbCt[ i ].u = uv.X();
+        pointsCbCt[ i ].v = uv.Y();
+      }
+      sideCbCt = StdMeshers_FaceSide::New( pointsCbCt, aFace );
+    }
+    // Make Cb quad
+    FaceQuadStruct* qCb = new FaceQuadStruct( quad->face );
+    myQuadList.push_back( FaceQuadStruct::Ptr( qCb ));
+    qCb->side.resize(4);
+    qCb->side[0] = quad->side[0];
+    qCb->side[1] = sideRCb;
+    qCb->side[2] = sideCbCt;
+    qCb->side[3] = sideLCb;
+    // Make L quad
+    FaceQuadStruct* qL = new FaceQuadStruct( quad->face );
+    myQuadList.push_back( FaceQuadStruct::Ptr( qL ));
+    qL->side.resize(4);
+    qL->side[0] = sideLCb;
+    qL->side[1] = sideLCt;
+    qL->side[2] = quad->side[2];
+    qL->side[3] = quad->side[3];
+    qL->side[2].to = lw+1;
+    // Make R quad
+    FaceQuadStruct* qR = new FaceQuadStruct( quad->face );
+    myQuadList.push_back( FaceQuadStruct::Ptr( qR ));
+    qR->side.resize(4);
+    qR->side[0] = sideRCb;
+    qR->side[0].from = lw;
+    qR->side[0].to   = -1;
+    qR->side[1] = quad->side[1];
+    qR->side[2] = quad->side[2];
+    qR->side[2].from = nb + lw + tfrom;
+    qR->side[3] = sideRCt;
+    // Make Ct from the main quad
+    FaceQuadStruct::Ptr qCt = quad;
+    qCt->side[0] = sideCbCt;
+    qCt->side[1] = sideRCt;
+    qCt->side[2].from = lw + tfrom;
+    qCt->side[2].to   = nt - lw + tfrom;
+    qCt->side[3] = sideLCt;
+    qCt->uv_grid.clear();
+
+    // Connect sides
+    qCb->side[3].AddContact( lw, & qCb->side[2], 0 );
+    qCb->side[3].AddContact( lw, & qCt->side[3], 0 );
+    qCt->side[3].AddContact(  0, & qCt->side[0], 0 );
+    qCt->side[0].AddContact(  0, & qL ->side[0], lw );
+    qL ->side[0].AddContact( lw, & qL ->side[1], 0 );
+    qL ->side[0].AddContact( lw, & qCb->side[2], 0 );
+    //
+    qCb->side[1].AddContact( lw, & qCb->side[2], lw );
+    qCb->side[1].AddContact( lw, & qCt->side[1], 0 );
+    qCt->side[0].AddContact( lw, & qCt->side[1], 0 );
+    qCt->side[0].AddContact( lw, & qR ->side[0], lw );
+    qR ->side[3].AddContact( lw, & qR ->side[0], lw );
+    qR ->side[3].AddContact( lw, & qCb->side[2], lw );
+
+    if ( dh == dv )
+      return computeQuadDominant( aMesh, aFace );
+
+
+  }
+
   if ( dh > dv ) {
     addv = (dh-dv)/2;
     nbv  = nbv + addv;
@@ -1479,19 +1912,6 @@ bool StdMeshers_Quadrangle_2D::computeQuadPref (SMESH_Mesh &        aMesh,
     nbh  = nbh + addh;
   }
 
-  const vector<UVPtStruct>& uv_eb = quad->side[0]->GetUVPtStruct(true,0);
-  const vector<UVPtStruct>& uv_er = quad->side[1]->GetUVPtStruct(false,1);
-  const vector<UVPtStruct>& uv_et = quad->side[2]->GetUVPtStruct(true,1);
-  const vector<UVPtStruct>& uv_el = quad->side[3]->GetUVPtStruct(false,0);
-
-  if (uv_eb.size() != nb || uv_er.size() != nr || uv_et.size() != nt || uv_el.size() != nl)
-    return error(COMPERR_BAD_INPUT_MESH);
-
-  if ( !OldVersion )
-  {
-    // dh/2, Min(nb,nt), dh - dh/2, dv
-  }
-
   // arrays for normalized params
   TColStd_SequenceOfReal npb, npr, npt, npl;
   for (i=0; i<nb; i++) {
@@ -1523,11 +1943,6 @@ bool StdMeshers_Quadrangle_2D::computeQuadPref (SMESH_Mesh &        aMesh,
       npl.InsertAfter(1,npl.Value(2)-dpr);
     }
   }
-  
-  gp_XY a0(uv_eb.front().u, uv_eb.front().v);
-  gp_XY a1(uv_eb.back().u,  uv_eb.back().v);
-  gp_XY a2(uv_et.back().u,  uv_et.back().v);
-  gp_XY a3(uv_et.front().u, uv_et.front().v);
 
   int nnn = Min(nr,nl);
   // auxilary sequence of XY for creation nodes
@@ -1543,7 +1958,7 @@ bool StdMeshers_Quadrangle_2D::computeQuadPref (SMESH_Mesh &        aMesh,
       NodesL.SetValue(1,j,uv_el[j-1].node);
     if (dl>0) {
       // add top nodes
-      for (i=1; i<=dl; i++) 
+      for (i=1; i<=dl; i++)
         NodesL.SetValue(i+1,nl,uv_et[i].node);
       // create and add needed nodes
       TColgp_SequenceOfXY UVtmp;
@@ -1580,13 +1995,13 @@ bool StdMeshers_Quadrangle_2D::computeQuadPref (SMESH_Mesh &        aMesh,
           if (WisF) {