src/SMDS/SMDS_VolumeTool.cxx

   1 // Copyright (C) 2007-2011  CEA/DEN, EDF R&D, OPEN CASCADE
   2 //
   3 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
   4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
   5 //
   6 // This library is free software; you can redistribute it and/or
   7 // modify it under the terms of the GNU Lesser General Public
   8 // License as published by the Free Software Foundation; either
   9 // version 2.1 of the License.
  10 //
  11 // This library is distributed in the hope that it will be useful,
  12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  14 // Lesser General Public License for more details.
  15 //
  16 // You should have received a copy of the GNU Lesser General Public
  17 // License along with this library; if not, write to the Free Software
  18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  19 //
  20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  21 //
  22
  23 // File      : SMDS_VolumeTool.cxx
  24 // Created   : Tue Jul 13 12:22:13 2004
  25 // Author    : Edward AGAPOV (eap)
  26 //
  27 #ifdef _MSC_VER
  28 #pragma warning(disable:4786)
  29 #endif
  30
  31 #include "SMDS_VolumeTool.hxx"
  32
  33 #include "SMDS_MeshElement.hxx"
  34 #include "SMDS_MeshNode.hxx"
  35 #include "SMDS_VtkVolume.hxx"
  36 #include "SMDS_Mesh.hxx"
  37
  38 #include "utilities.h"
  39
  40 #include <map>
  41 #include <float.h>
  42 #include <math.h>
  43
  44 using namespace std;
  45
  46 // ======================================================
  47 // Node indices in faces depending on volume orientation
  48 // making most faces normals external
  49 // ======================================================
  50
  51 /*
  52 //           N3
  53 //           +
  54 //          /|\
  55 //         / | \
  56 //        /  |  \
  57 //    N0 +---|---+ N1                TETRAHEDRON
  58 //       \   |   /
  59 //        \  |  /
  60 //         \ | /
  61 //          \|/
  62 //           +
  63 //           N2
  64 */
  65 static int Tetra_F [4][4] = { // FORWARD == EXTERNAL
  66   { 0, 1, 2, 0 },              // All faces have external normals
  67   { 0, 3, 1, 0 },
  68   { 1, 3, 2, 1 },
  69   { 0, 2, 3, 0 }};
  70 static int Tetra_R [4][4] = { // REVERSED
  71   { 0, 1, 2, 0 },             // All faces but a bottom have external normals
  72   { 0, 1, 3, 0 },
  73   { 1, 2, 3, 1 },
  74   { 0, 3, 2, 0 }};
  75 static int Tetra_RE [4][4] = { // REVERSED -> FORWARD (EXTERNAL)
  76   { 0, 2, 1, 0 },              // All faces have external normals
  77   { 0, 1, 3, 0 },
  78   { 1, 2, 3, 1 },
  79   { 0, 3, 2, 0 }};
  80 static int Tetra_nbN [] = { 3, 3, 3, 3 };
  81
  82 //
  83 //     PYRAMID
  84 //
  85 static int Pyramid_F [5][5] = { // FORWARD == EXTERNAL
  86   { 0, 1, 2, 3, 0 },            // All faces have external normals
  87   { 0, 4, 1, 0, 4 },
  88   { 1, 4, 2, 1, 4 },
  89   { 2, 4, 3, 2, 4 },
  90   { 3, 4, 0, 3, 4 }};
  91 static int Pyramid_R [5][5] = { // REVERSED
  92   { 0, 1, 2, 3, 0 },            // All faces but a bottom have external normals
  93   { 0, 1, 4, 0, 4 },
  94   { 1, 2, 4, 1, 4 },
  95   { 2, 3, 4, 2, 4 },
  96   { 3, 0, 4, 3, 4 }};
  97 static int Pyramid_RE [5][5] = { // REVERSED -> FORWARD (EXTERNAL)
  98   { 0, 3, 2, 1, 0 },             // All faces but a bottom have external normals
  99   { 0, 1, 4, 0, 4 },
 100   { 1, 2, 4, 1, 4 },
 101   { 2, 3, 4, 2, 4 },
 102   { 3, 0, 4, 3, 4 }};
 103 static int Pyramid_nbN [] = { 4, 3, 3, 3, 3 };
 104
 105 /*
 106 //            + N4
 107 //           /|\
 108 //          / | \
 109 //         /  |  \
 110 //        /   |   \
 111 //    N3 +---------+ N5
 112 //       |    |    |
 113 //       |    + N1 |
 114 //       |   / \   |                PENTAHEDRON
 115 //       |  /   \  |
 116 //       | /     \ |
 117 //       |/       \|
 118 //    N0 +---------+ N2
 119 */
 120 static int Penta_F [5][5] = { // FORWARD
 121   { 0, 1, 2, 0, 0 },          // Top face has an internal normal, other - external
 122   { 3, 4, 5, 3, 3 },          // 0 is bottom, 1 is top face
 123   { 0, 2, 5, 3, 0 },
 124   { 1, 4, 5, 2, 1 },
 125   { 0, 3, 4, 1, 0 }};
 126 static int Penta_R [5][5] = { // REVERSED
 127   { 0, 1, 2, 0, 0 },          // Bottom face has an internal normal, other - external
 128   { 3, 4, 5, 3, 3 },          // 0 is bottom, 1 is top face
 129   { 0, 3, 5, 2, 0 },
 130   { 1, 2, 5, 4, 1 },
 131   { 0, 1, 4, 3, 0 }};
 132 static int Penta_FE [5][5] = { // FORWARD -> EXTERNAL
 133   { 0, 1, 2, 0, 0 },
 134   { 3, 5, 4, 3, 3 },
 135   { 0, 2, 5, 3, 0 },
 136   { 1, 4, 5, 2, 1 },
 137   { 0, 3, 4, 1, 0 }};
 138 static int Penta_RE [5][5] = { // REVERSED -> EXTERNAL
 139   { 0, 2, 1, 0, 0 },
 140   { 3, 4, 5, 3, 3 },
 141   { 0, 3, 5, 2, 0 },
 142   { 1, 2, 5, 4, 1 },
 143   { 0, 1, 4, 3, 0 }};
 144 static int Penta_nbN [] = { 3, 3, 4, 4, 4 };
 145
 146 /*
 147 //         N5+----------+N6
 148 //          /|         /|
 149 //         / |        / |
 150 //        /  |       /  |
 151 //     N4+----------+N7 |
 152 //       |   |      |   |           HEXAHEDRON
 153 //       |   |      |   |
 154 //       |   |      |   |
 155 //       | N1+------|---+N2
 156 //       |  /       |  /
 157 //       | /        | /
 158 //       |/         |/
 159 //     N0+----------+N3
 160 */
 161 static int Hexa_F [6][5] = { // FORWARD
 162   { 0, 1, 2, 3, 0 },         // opposite faces are neighbouring,
 163   { 4, 5, 6, 7, 4 },         // odd face(1,3,5) normal is internal, even(0,2,4) - external
 164   { 1, 0, 4, 5, 1 },         // same index nodes of opposite faces are linked
 165   { 2, 3, 7, 6, 2 },
 166   { 0, 3, 7, 4, 0 },
 167   { 1, 2, 6, 5, 1 }};
 168 // static int Hexa_R [6][5] = { // REVERSED
 169 //   { 0, 3, 2, 1, 0 },         // opposite faces are neighbouring,
 170 //   { 4, 7, 6, 5, 4 },         // odd face(1,3,5) normal is external, even(0,2,4) - internal
 171 //   { 1, 5, 4, 0, 1 },         // same index nodes of opposite faces are linked
 172 //   { 2, 6, 7, 3, 2 },
 173 //   { 0, 4, 7, 3, 0 },
 174 //   { 1, 5, 6, 2, 1 }};
 175 static int Hexa_FE [6][5] = { // FORWARD -> EXTERNAL
 176   { 0, 1, 2, 3, 0 } ,         // opposite faces are neighbouring,
 177   { 4, 7, 6, 5, 4 },          // all face normals are external,
 178   { 0, 4, 5, 1, 0 },          // links in opposite faces: 0-0, 1-3, 2-2, 3-1
 179   { 3, 2, 6, 7, 3 },
 180   { 0, 3, 7, 4, 0 },
 181   { 1, 5, 6, 2, 1 }};
 182 static int Hexa_RE [6][5] = { // REVERSED -> EXTERNAL
 183   { 0, 3, 2, 1, 0 },          // opposite faces are neighbouring,
 184   { 4, 5, 6, 7, 4 },          // all face normals are external,
 185   { 0, 1, 5, 4, 0 },          // links in opposite faces: 0-0, 1-3, 2-2, 3-1
 186   { 3, 7, 6, 2, 3 },
 187   { 0, 4, 7, 3, 0 },
 188   { 1, 2, 6, 5, 1 }};
 189 static int Hexa_nbN [] = { 4, 4, 4, 4, 4, 4 };
 190
 191
 192 /*
 193 //           N3
 194 //           +
 195 //          /|\
 196 //        7/ | \8
 197 //        /  |4 \                    QUADRATIC
 198 //    N0 +---|---+ N1                TETRAHEDRON
 199 //       \   +9  /
 200 //        \  |  /
 201 //        6\ | /5
 202 //          \|/
 203 //           +
 204 //           N2
 205 */
 206 static int QuadTetra_F [4][7] = { // FORWARD == EXTERNAL
 207   { 0, 4, 1, 5, 2, 6, 0 },        // All faces have external normals
 208   { 0, 7, 3, 8, 1, 4, 0 },
 209   { 1, 8, 3, 9, 2, 5, 1 },
 210   { 0, 6, 2, 9, 3, 7, 0 }};
 211 static int QuadTetra_R [4][7] = { // REVERSED
 212   { 0, 4, 1, 5, 2, 6, 0 },        // All faces but a bottom have external normals
 213   { 0, 4, 1, 8, 3, 7, 0 },
 214   { 1, 5, 2, 9, 3, 8, 1 },
 215   { 0, 7, 3, 9, 2, 6, 0 }};
 216 static int QuadTetra_RE [4][7] = { // REVERSED -> FORWARD (EXTERNAL)
 217   { 0, 6, 2, 5, 1, 4, 0 },              // All faces have external normals
 218   { 0, 4, 1, 8, 3, 7, 0 },
 219   { 1, 5, 2, 9, 3, 8, 1 },
 220   { 0, 7, 3, 9, 2, 6, 0 }};
 221 static int QuadTetra_nbN [] = { 6, 6, 6, 6 };
 222
 223 //
 224 //     QUADRATIC
 225 //     PYRAMID
 226 //
 227 //            +4
 228 //
 229 //
 230 //       10+-----+11
 231 //         |     |        9 - middle point for (0,4) etc.
 232 //         |     |
 233 //        9+-----+12
 234 //
 235 //            6
 236 //      1+----+----+2
 237 //       |         |
 238 //       |         |
 239 //      5+         +7
 240 //       |         |
 241 //       |         |
 242 //      0+----+----+3
 243 //            8
 244 static int QuadPyram_F [5][9] = {  // FORWARD == EXTERNAL
 245   { 0, 5, 1, 6, 2, 7, 3, 8, 0 },   // All faces have external normals
 246   { 0, 9, 4, 10,1, 5, 0, 4, 4 },
 247   { 1, 10,4, 11,2, 6, 1, 4, 4 },
 248   { 2, 11,4, 12,3, 7, 2, 4, 4 },
 249   { 3, 12,4, 9, 0, 8, 3, 4, 4 }};
 250 static int QuadPyram_R [5][9] = {  // REVERSED
 251   { 0, 5, 1, 6, 2, 7, 3, 8, 0 },   // All faces but a bottom have external normals
 252   { 0, 5, 1, 10,4, 9, 0, 4, 4 },
 253   { 1, 6, 2, 11,4, 10,1, 4, 4 },
 254   { 2, 7, 3, 12,4, 11,2, 4, 4 },
 255   { 3, 8, 0, 9, 4, 12,3, 4, 4 }};
 256 static int QuadPyram_RE [5][9] = { // REVERSED -> FORWARD (EXTERNAL)
 257   { 0, 8, 3, 7, 2, 6, 1, 5, 0 },   // All faces but a bottom have external normals
 258   { 0, 5, 1, 10,4, 9, 0, 4, 4 },
 259   { 1, 6, 2, 11,4, 10,1, 4, 4 },
 260   { 2, 7, 3, 12,4, 11,2, 4, 4 },
 261   { 3, 8, 0, 9, 4, 12,3, 4, 4 }};
 262 static int QuadPyram_nbN [] = { 8, 6, 6, 6, 6 };
 263
 264 /*
 265 //            + N4
 266 //           /|\
 267 //         9/ | \10
 268 //         /  |  \
 269 //        /   |   \
 270 //    N3 +----+----+ N5
 271 //       |    |11  |
 272 //       |    |    |
 273 //       |    +13  |                QUADRATIC
 274 //       |    |    |                PENTAHEDRON
 275 //       |    |    |
 276 //       |    |    |
 277 //       |    |    |
 278 //     12+    |    +14
 279 //       |    |    |
 280 //       |    |    |
 281 //       |    + N1 |
 282 //       |   / \   |
 283 //       | 6/   \7 |
 284 //       | /     \ |
 285 //       |/       \|
 286 //    N0 +---------+ N2
 287 //            8
 288 */
 289 static int QuadPenta_F [5][9] = {  // FORWARD
 290   { 0, 6, 1, 7, 2, 8, 0, 0, 0 },   // Top face has an internal normal, other - external
 291   { 3, 9, 4, 10,5, 11,3, 3, 3 },   // 0 is bottom, 1 is top face
 292   { 0, 8, 2, 14,5, 11,3, 12,0 },
 293   { 1, 13,4, 10,5, 14,2, 7, 1 },
 294   { 0, 12,3, 9, 4, 13,1, 6, 0 }};
 295 static int QuadPenta_R [5][9] = { // REVERSED
 296   { 0, 6, 1, 7, 2, 8, 0, 0, 0 },  // Bottom face has an internal normal, other - external
 297   { 3, 9, 4, 10,5, 11,3, 3, 3 },  // 0 is bottom, 1 is top face
 298   { 0, 12,3, 11,5, 14,2, 8, 0 },
 299   { 1, 7, 2, 14,5, 10,4, 13,1 },
 300   { 0, 6, 1, 13,4, 9, 3, 12,0 }};
 301 static int QuadPenta_FE [5][9] = { // FORWARD -> EXTERNAL
 302   { 0, 6, 1, 7, 2, 8, 0, 0, 0 },
 303   { 3,11, 5, 10,4, 9, 3, 3, 3 },
 304   { 0, 8, 2, 14,5, 11,3, 12,0 },
 305   { 1, 13,4, 10,5, 14,2, 7, 1 },
 306   { 0, 12,3, 9, 4, 13,1, 6, 0 }};
 307 static int QuadPenta_RE [5][9] = { // REVERSED -> EXTERNAL
 308   { 0, 8, 2, 7, 1, 6, 0, 0, 0 },
 309   { 3, 9, 4, 10,5, 11,3, 3, 3 },
 310   { 0, 12,3, 11,5, 14,2, 8, 0 },
 311   { 1, 7, 2, 14,5, 10,4, 13,1 },
 312   { 0, 6, 1, 13,4, 9, 3, 12,0 }};
 313 static int QuadPenta_nbN [] = { 6, 6, 8, 8, 8 };
 314
 315 /*
 316 //                 13
 317 //         N5+-----+-----+N6
 318 //          /|          /|
 319 //       12+ |       14+ |
 320 //        /  |        /  |
 321 //     N4+-----+-----+N7 |           QUADRATIC
 322 //       |   | 15    |   |           HEXAHEDRON
 323 //       |   |       |   |
 324 //       | 17+       |   +18
 325 //       |   |       |   |
 326 //       |   |       |   |
 327 //       |   |       |   |
 328 //     16+   |       +19 |
 329 //       |   |       |   |
 330 //       |   |     9 |   |
 331 //       | N1+-----+-|---+N2
 332 //       |  /        |  /
 333 //       | +8        | +10
 334 //       |/          |/
 335 //     N0+-----+-----+N3
 336 //             11
 337 */
 338 static int QuadHexa_F [6][9] = {  // FORWARD
 339   { 0, 8, 1, 9, 2, 10,3, 11,0 },  // opposite faces are neighbouring,
 340   { 4, 12,5, 13,6, 14,7, 15,4 },  // odd face(1,3,5) normal is internal, even(0,2,4) - external
 341   { 1, 8, 0, 16,4, 12,5, 17,1 },  // same index nodes of opposite faces are linked
 342   { 2, 10,3, 19,7, 14,6, 18,2 },
 343   { 0, 11,3, 19,7, 15,4, 16,0 },
 344   { 1, 9, 2, 18,6, 13,5, 17,1 }};
 345 // static int Hexa_R [6][5] = { // REVERSED
 346 //   { 0, 3, 2, 1, 0 },         // opposite faces are neighbouring,
 347 //   { 4, 7, 6, 5, 4 },         // odd face(1,3,5) normal is external, even(0,2,4) - internal
 348 //   { 1, 5, 4, 0, 1 },         // same index nodes of opposite faces are linked
 349 //   { 2, 6, 7, 3, 2 },
 350 //   { 0, 4, 7, 3, 0 },
 351 //   { 1, 5, 6, 2, 1 }};
 352 static int QuadHexa_FE [6][9] = {  // FORWARD -> EXTERNAL
 353   { 0, 8, 1, 9, 2, 10,3, 11,0 },   // opposite faces are neighbouring,
 354   { 4, 15,7, 14,6, 13,5, 12,4 },   // all face normals are external,
 355   { 0, 16,4, 12,5, 17,1, 8, 0 },   // links in opposite faces: 0-0, 1-3, 2-2, 3-1
 356   { 3, 10,2, 18,6, 14,7, 19,3 },
 357   { 0, 11,3, 19,7, 15,4, 16,0 },
 358   { 1, 17,5, 13,6, 18,2, 9, 1 }};
 359 static int QuadHexa_RE [6][9] = {  // REVERSED -> EXTERNAL
 360   { 0, 11,3, 10,2, 9, 1, 8, 0 },   // opposite faces are neighbouring,
 361   { 4, 12,5, 13,6, 14,7, 15,4 },   // all face normals are external,
 362   { 0, 8, 1, 17,5, 12,4, 16,0 },   // links in opposite faces: 0-0, 1-3, 2-2, 3-1
 363   { 3, 19,7, 14,6, 18,2, 10,3 },
 364   { 0, 16,4, 15,7, 19,3, 11,0 },
 365   { 1, 9, 2, 18,6, 13,5, 17,1 }};
 366 static int QuadHexa_nbN [] = { 8, 8, 8, 8, 8, 8 };
 367
 368
 369 // ========================================================
 370 // to perform some calculations without linkage to CASCADE
 371 // ========================================================
 372 namespace
 373 {
 374 struct XYZ {
 375   double x;
 376   double y;
 377   double z;
 378   XYZ()                               { x = 0; y = 0; z = 0; }
 379   XYZ( double X, double Y, double Z ) { x = X; y = Y; z = Z; }
 380   XYZ( const XYZ& other )             { x = other.x; y = other.y; z = other.z; }
 381   XYZ( const SMDS_MeshNode* n )       { x = n->X(); y = n->Y(); z = n->Z(); }
 382   inline XYZ operator-( const XYZ& other );
 383   inline XYZ Crossed( const XYZ& other );
 384   inline double Dot( const XYZ& other );
 385   inline double Magnitude();
 386 };
 387 inline XYZ XYZ::operator-( const XYZ& Right ) {
 388   return XYZ(x - Right.x, y - Right.y, z - Right.z);
 389 }
 390 inline XYZ XYZ::Crossed( const XYZ& Right ) {
 391   return XYZ (y * Right.z - z * Right.y,
 392               z * Right.x - x * Right.z,
 393               x * Right.y - y * Right.x);
 394 }
 395 inline double XYZ::Dot( const XYZ& Other ) {
 396   return(x * Other.x + y * Other.y + z * Other.z);
 397 }
 398 inline double XYZ::Magnitude() {
 399   return sqrt (x * x + y * y + z * z);
 400 }
 401 }
 402
 403 //=======================================================================
 404 //function : SMDS_VolumeTool
 405 //purpose  :
 406 //=======================================================================
 407
 408 SMDS_VolumeTool::SMDS_VolumeTool ()
 409      : myVolume( 0 ),
 410        myPolyedre( 0 ),
 411        myVolForward( true ),
 412        myNbFaces( 0 ),
 413        myVolumeNbNodes( 0 ),
 414        myVolumeNodes( NULL ),
 415        myExternalFaces( false ),
 416        myFaceNbNodes( 0 ),
 417        myCurFace( -1 ),
 418        myFaceNodeIndices( NULL ),
 419        myFaceNodes( NULL )
 420 {
 421   //MESSAGE("******************************************************** SMDS_VolumeToo");
 422 }
 423
 424 //=======================================================================
 425 //function : SMDS_VolumeTool
 426 //purpose  :
 427 //=======================================================================
 428
 429 SMDS_VolumeTool::SMDS_VolumeTool (const SMDS_MeshElement* theVolume)
 430      : myVolume( 0 ),
 431        myPolyedre( 0 ),
 432        myVolForward( true ),
 433        myNbFaces( 0 ),
 434        myVolumeNbNodes( 0 ),
 435        myVolumeNodes( NULL ),
 436        myExternalFaces( false ),
 437        myFaceNbNodes( 0 ),
 438        myCurFace( -1 ),
 439        myFaceNodeIndices( NULL ),
 440        myFaceNodes( NULL )
 441 {
 442   //MESSAGE("******************************************************** SMDS_VolumeToo");
 443   Set( theVolume );
 444 }
 445
 446 //=======================================================================
 447 //function : SMDS_VolumeTool
 448 //purpose  :
 449 //=======================================================================
 450
 451 SMDS_VolumeTool::~SMDS_VolumeTool()
 452 {
 453   if ( myVolumeNodes != NULL ) delete [] myVolumeNodes;
 454   if ( myFaceNodes != NULL   ) delete [] myFaceNodes;
 455
 456   myFaceNodeIndices = NULL;
 457   myVolumeNodes = myFaceNodes = NULL;
 458 }
 459
 460 //=======================================================================
 461 //function : SetVolume
 462 //purpose  : Set volume to iterate on
 463 //=======================================================================
 464
 465 bool SMDS_VolumeTool::Set (const SMDS_MeshElement* theVolume)
 466 {
 467   myVolume = 0;
 468   myPolyedre = 0;
 469
 470   myVolForward = true;
 471   myNbFaces = 0;
 472   myVolumeNbNodes = 0;
 473   if (myVolumeNodes != NULL) {
 474     delete [] myVolumeNodes;
 475     myVolumeNodes = NULL;
 476   }
 477
 478   myExternalFaces = false;
 479   myFaceNbNodes = 0;
 480
 481   myCurFace = -1;
 482   myFaceNodeIndices = NULL;
 483   if (myFaceNodes != NULL) {
 484     delete [] myFaceNodes;
 485     myFaceNodes = NULL;
 486   }
 487
 488   if ( theVolume && theVolume->GetType() == SMDSAbs_Volume )
 489   {
 490     myVolume = theVolume;
 491
 492     myNbFaces = theVolume->NbFaces();
 493     myVolumeNbNodes = theVolume->NbNodes();
 494
 495     // set volume nodes
 496     int iNode = 0;
 497     myVolumeNodes = new const SMDS_MeshNode* [myVolumeNbNodes];
 498     SMDS_ElemIteratorPtr nodeIt = myVolume->nodesIterator();
 499     while ( nodeIt->more() ) {
 500       myVolumeNodes[ iNode++ ] = static_cast<const SMDS_MeshNode*>( nodeIt->next() );
 501     }
 502
 503     if (myVolume->IsPoly()) {
 504       myPolyedre = dynamic_cast<const SMDS_VtkVolume*>( myVolume );
 505       if (!myPolyedre) {
 506         MESSAGE("Warning: bad volumic element");
 507         return false;
 508       }
 509     }
 510     else {
 511       switch ( myVolumeNbNodes ) {
 512       case 4:
 513       case 5:
 514       case 6:
 515       case 8:
 516       case 10:
 517       case 13:
 518       case 15:
 519       case 20: {
 520         // define volume orientation
 521         XYZ botNormal;
 522         GetFaceNormal( 0, botNormal.x, botNormal.y, botNormal.z );
 523         const SMDS_MeshNode* botNode = myVolumeNodes[ 0 ];
 524         int topNodeIndex = myVolume->NbCornerNodes() - 1;
 525         while ( !IsLinked( 0, topNodeIndex, /*ignoreMediumNodes=*/true )) --topNodeIndex;
 526         const SMDS_MeshNode* topNode = myVolumeNodes[ topNodeIndex ];
 527         XYZ upDir (topNode->X() - botNode->X(),
 528                    topNode->Y() - botNode->Y(),
 529                    topNode->Z() - botNode->Z() );
 530         myVolForward = ( botNormal.Dot( upDir ) < 0 );
 531         break;
 532       }
 533       default:
 534         break;
 535       }
 536     }
 537   }
 538   return ( myVolume != 0 );
 539 }
 540
 541 //=======================================================================
 542 //function : Inverse
 543 //purpose  : Inverse volume
 544 //=======================================================================
 545
 546 #define SWAP_NODES(nodes,i1,i2)           \
 547 {                                         \
 548   const SMDS_MeshNode* tmp = nodes[ i1 ]; \
 549   nodes[ i1 ] = nodes[ i2 ];              \
 550   nodes[ i2 ] = tmp;                      \
 551 }
 552 void SMDS_VolumeTool::Inverse ()
 553 {
 554   if ( !myVolume ) return;
 555
 556   if (myVolume->IsPoly()) {
 557     MESSAGE("Warning: attempt to inverse polyhedral volume");
 558     return;
 559   }
 560
 561   myVolForward = !myVolForward;
 562   myCurFace = -1;
 563
 564   // inverse top and bottom faces
 565   switch ( myVolumeNbNodes ) {
 566   case 4:
 567     SWAP_NODES( myVolumeNodes, 1, 2 );
 568     break;
 569   case 5:
 570     SWAP_NODES( myVolumeNodes, 1, 3 );
 571     break;
 572   case 6:
 573     SWAP_NODES( myVolumeNodes, 1, 2 );
 574     SWAP_NODES( myVolumeNodes, 4, 5 );
 575     break;
 576   case 8:
 577     SWAP_NODES( myVolumeNodes, 1, 3 );
 578     SWAP_NODES( myVolumeNodes, 5, 7 );
 579     break;
 580
 581   case 10:
 582     SWAP_NODES( myVolumeNodes, 1, 2 );
 583     SWAP_NODES( myVolumeNodes, 4, 6 );
 584     SWAP_NODES( myVolumeNodes, 8, 9 );
 585     break;
 586   case 13:
 587     SWAP_NODES( myVolumeNodes, 1, 3 );
 588     SWAP_NODES( myVolumeNodes, 5, 8 );
 589     SWAP_NODES( myVolumeNodes, 6, 7 );
 590     SWAP_NODES( myVolumeNodes, 10, 12 );
 591     break;
 592   case 15:
 593     SWAP_NODES( myVolumeNodes, 1, 2 );
 594     SWAP_NODES( myVolumeNodes, 4, 5 );
 595     SWAP_NODES( myVolumeNodes, 6, 8 );
 596     SWAP_NODES( myVolumeNodes, 9, 11 );
 597     SWAP_NODES( myVolumeNodes, 13, 14 );
 598     break;
 599   case 20:
 600     SWAP_NODES( myVolumeNodes, 1, 3 );
 601     SWAP_NODES( myVolumeNodes, 5, 7 );
 602     SWAP_NODES( myVolumeNodes, 8, 11 );
 603     SWAP_NODES( myVolumeNodes, 9, 10 );
 604     SWAP_NODES( myVolumeNodes, 12, 15 );
 605     SWAP_NODES( myVolumeNodes, 13, 14 );
 606     SWAP_NODES( myVolumeNodes, 17, 19 );
 607     break;
 608   default:;
 609   }
 610 }
 611
 612 //=======================================================================
 613 //function : GetVolumeType
 614 //purpose  :
 615 //=======================================================================
 616
 617 SMDS_VolumeTool::VolumeType SMDS_VolumeTool::GetVolumeType() const
 618 {
 619   if ( myPolyedre )
 620     return POLYHEDA;
 621
 622   if ( myVolume ) {
 623 //    static const VolumeType types[] = {
 624 //      TETRA,     // myVolumeNbNodes = 4
 625 //      PYRAM,     // myVolumeNbNodes = 5
 626 //      PENTA,     // myVolumeNbNodes = 6
 627 //      UNKNOWN,   // myVolumeNbNodes = 7
 628 //      HEXA       // myVolumeNbNodes = 8
 629 //    };
 630 //    return types[ myVolumeNbNodes - 4 ];
 631     switch(myVolumeNbNodes) {
 632     case 4: return TETRA; break;
 633     case 5: return PYRAM; break;
 634     case 6: return PENTA; break;
 635     case 8: return HEXA; break;
 636     case 10: return QUAD_TETRA; break;
 637     case 13: return QUAD_PYRAM; break;
 638     case 15: return QUAD_PENTA; break;
 639     case 20: return QUAD_HEXA; break;
 640     default: break;
 641     }
 642   }
 643
 644   return UNKNOWN;
 645 }
 646
 647 //=======================================================================
 648 //function : getTetraVolume
 649 //purpose  :
 650 //=======================================================================
 651
 652 static double getTetraVolume(const SMDS_MeshNode* n1,
 653                              const SMDS_MeshNode* n2,
 654                              const SMDS_MeshNode* n3,
 655                              const SMDS_MeshNode* n4)
 656 {
 657   double X1 = n1->X();
 658   double Y1 = n1->Y();
 659   double Z1 = n1->Z();
 660
 661   double X2 = n2->X();
 662   double Y2 = n2->Y();
 663   double Z2 = n2->Z();
 664
 665   double X3 = n3->X();
 666   double Y3 = n3->Y();
 667   double Z3 = n3->Z();
 668
 669   double X4 = n4->X();
 670   double Y4 = n4->Y();
 671   double Z4 = n4->Z();
 672
 673   double Q1 = -(X1-X2)*(Y3*Z4-Y4*Z3);
 674   double Q2 =  (X1-X3)*(Y2*Z4-Y4*Z2);
 675   double R1 = -(X1-X4)*(Y2*Z3-Y3*Z2);
 676   double R2 = -(X2-X3)*(Y1*Z4-Y4*Z1);
 677   double S1 =  (X2-X4)*(Y1*Z3-Y3*Z1);
 678   double S2 = -(X3-X4)*(Y1*Z2-Y2*Z1);
 679
 680   return (Q1+Q2+R1+R2+S1+S2)/6.0;
 681 }
 682
 683 //=======================================================================
 684 //function : GetSize
 685 //purpose  : Return element volume
 686 //=======================================================================
 687
 688 double SMDS_VolumeTool::GetSize() const
 689 {
 690   double V = 0.;
 691   if ( !myVolume )
 692     return 0.;
 693
 694   if ( myVolume->IsPoly() )
 695   {
 696     if ( !myPolyedre )
 697       return 0.;
 698
 699     SMDS_Mesh *mesh = SMDS_Mesh::_meshList[myPolyedre->getMeshId()];
 700     // split a polyhedron into tetrahedrons
 701
 702     SMDS_VolumeTool* me = const_cast< SMDS_VolumeTool* > ( this );
 703     XYZ baryCenter;
 704     me->GetBaryCenter(baryCenter.x, baryCenter.y, baryCenter.z);
 705     SMDS_MeshNode *bcNode = mesh->AddNode( baryCenter.x, baryCenter.y, baryCenter.z );
 706
 707     for ( int f = 0; f < NbFaces(); ++f )
 708     {
 709       bool externalFace = me->IsFaceExternal( f ); // it calls setFace()
 710       for ( int n = 2; n < myFaceNbNodes; ++n )
 711       {
 712         double Vn = getTetraVolume( myFaceNodes[ 0 ],
 713                                     myFaceNodes[ n-1 ],
 714                                     myFaceNodes[ n ],
 715                                     bcNode );
 716 ///         cout <<"++++   " << Vn << "   nodes " <<myFaceNodes[ 0 ]->GetID() << " " <<myFaceNodes[ n-1 ]->GetID() << " " <<myFaceNodes[ n ]->GetID() << "        < " << V << endl;
 717         V += externalFace ? -Vn : Vn;
 718       }
 719     }
 720     mesh->RemoveNode(bcNode);
 721   }
 722   else
 723   {
 724     const static int ind[] = {
 725       0, 1, 3, 6, 11, 19, 32, 46, 66};
 726     const static int vtab[][4] = {
 727       // tetrahedron
 728       { 0, 1, 2, 3 },
 729       // pyramid
 730       { 0, 1, 3, 4 },
 731       { 1, 2, 3, 4 },
 732       // pentahedron
 733       { 0, 1, 2, 3 },
 734       { 1, 5, 3, 4 },
 735       { 1, 5, 2, 3 },
 736       // hexahedron
 737       { 1, 4, 3, 0 },
 738       { 4, 1, 6, 5 },
 739       { 1, 3, 6, 2 },
 740       { 4, 6, 3, 7 },
 741       { 1, 4, 6, 3 },
 742
 743       // quadratic tetrahedron
 744       { 0, 4, 6, 7 },
 745       { 1, 5, 4, 8 },
 746       { 2, 6, 5, 9 },
 747       { 7, 8, 9, 3 },
 748       { 4, 6, 7, 9 },
 749       { 4, 5, 6, 9 },
 750       { 4, 7, 8, 9 },
 751       { 4, 5, 9, 8 },
 752
 753       // quadratic pyramid
 754       { 0, 5, 8, 9 },
 755       { 1, 5,10, 6 },
 756       { 2, 6,11, 7 },
 757       { 3, 7,12, 8 },
 758       { 4, 9,11,10 },
 759       { 4, 9,12,11 },
 760       { 10, 5, 9, 8 },
 761       { 10, 8, 9,12 },
 762       { 10, 8,12, 7 },
 763       { 10, 7,12,11 },
 764       { 10, 7,11, 6 },
 765       { 10, 5, 8, 6 },
 766       { 10, 6, 8, 7 },
 767
 768       // quadratic pentahedron
 769       { 12, 0, 8, 6 },
 770       { 12, 8, 7, 6 },
 771       { 12, 8, 2, 7 },
 772       { 12, 6, 7, 1 },
 773       { 12, 1, 7,13 },
 774       { 12, 7, 2,13 },
 775       { 12, 2,14,13 },
 776
 777       { 12, 3, 9,11 },
 778       { 12,11, 9,10 },
 779       { 12,11,10, 5 },
 780       { 12, 9, 4,10 },
 781       { 12,14, 5,10 },
 782       { 12,14,10, 4 },
 783       { 12,14, 4,13 },
 784
 785       // quadratic hexahedron
 786       { 16, 0,11, 8 },
 787       { 16,11, 9, 8 },
 788       { 16, 8, 9, 1 },
 789       { 16,11, 3,10 },
 790       { 16,11,10, 9 },
 791       { 16,10, 2, 9 },
 792       { 16, 3,19, 2 },
 793       { 16, 2,19,18 },
 794       { 16, 2,18,17 },
 795       { 16, 2,17, 1 },
 796
 797       { 16, 4,12,15 },
 798       { 16,12, 5,13 },
 799       { 16,12,13,15 },
 800       { 16,13, 6,14 },
 801       { 16,13,14,15 },
 802       { 16,14, 7,15 },
 803       { 16, 6, 5,17 },
 804       { 16,18, 6,17 },
 805       { 16,18, 7, 6 },
 806       { 16,18,19, 7 },
 807
 808     };
 809
 810     int type = GetVolumeType();
 811     int n1 = ind[type];
 812     int n2 = ind[type+1];
 813
 814     for (int i = n1; i <  n2; i++) {
 815       V -= getTetraVolume( myVolumeNodes[ vtab[i][0] ],
 816                            myVolumeNodes[ vtab[i][1] ],
 817                            myVolumeNodes[ vtab[i][2] ],
 818                            myVolumeNodes[ vtab[i][3] ]);
 819     }
 820   }
 821   return V;
 822 }
 823
 824 //=======================================================================
 825 //function : GetBaryCenter
 826 //purpose  :
 827 //=======================================================================
 828
 829 bool SMDS_VolumeTool::GetBaryCenter(double & X, double & Y, double & Z) const
 830 {
 831   X = Y = Z = 0.;
 832   if ( !myVolume )
 833     return false;
 834
 835   for ( int i = 0; i < myVolumeNbNodes; i++ ) {
 836     X += myVolumeNodes[ i ]->X();
 837     Y += myVolumeNodes[ i ]->Y();
 838     Z += myVolumeNodes[ i ]->Z();
 839   }
 840   X /= myVolumeNbNodes;
 841   Y /= myVolumeNbNodes;
 842   Z /= myVolumeNbNodes;
 843
 844   return true;
 845 }
 846
 847 //================================================================================
 848 /*!
 849  * \brief Classify a point
 850  *  \param tol - thickness of faces
 851  */
 852 //================================================================================
 853
 854 bool SMDS_VolumeTool::IsOut(double X, double Y, double Z, double tol)
 855 {
 856   // LIMITATION: for convex volumes only
 857   XYZ p( X,Y,Z );
 858   for ( int iF = 0; iF < myNbFaces; ++iF )
 859   {
 860     XYZ faceNormal;
 861     if ( !GetFaceNormal( iF, faceNormal.x, faceNormal.y, faceNormal.z ))
 862       continue;
 863     if ( !IsFaceExternal( iF ))
 864       faceNormal = XYZ() - faceNormal; // reverse
 865
 866     XYZ face2p( p - XYZ( myFaceNodes[0] ));
 867     if ( face2p.Dot( faceNormal ) > tol )
 868       return true;
 869   }
 870   return false;
 871 }
 872
 873 //=======================================================================
 874 //function : SetExternalNormal
 875 //purpose  : Node order will be so that faces normals are external
 876 //=======================================================================
 877
 878 void SMDS_VolumeTool::SetExternalNormal ()
 879 {
 880   myExternalFaces = true;
 881   myCurFace = -1;
 882 }
 883
 884 //=======================================================================
 885 //function : NbFaceNodes
 886 //purpose  : Return number of nodes in the array of face nodes
 887 //=======================================================================
 888
 889 int SMDS_VolumeTool::NbFaceNodes( int faceIndex )
 890 {
 891     if ( !setFace( faceIndex ))
 892       return 0;
 893     return myFaceNbNodes;
 894 }
 895
 896 //=======================================================================
 897 //function : GetFaceNodes
 898 //purpose  : Return pointer to the array of face nodes.
 899 //           To comfort link iteration, the array
 900 //           length == NbFaceNodes( faceIndex ) + 1 and
 901 //           the last node == the first one.
 902 //=======================================================================
 903
 904 const SMDS_MeshNode** SMDS_VolumeTool::GetFaceNodes( int faceIndex )
 905 {
 906   if ( !setFace( faceIndex ))
 907     return 0;
 908   return myFaceNodes;
 909 }
 910
 911 //=======================================================================
 912 //function : GetFaceNodesIndices
 913 //purpose  : Return pointer to the array of face nodes indices
 914 //           To comfort link iteration, the array
 915 //           length == NbFaceNodes( faceIndex ) + 1 and
 916 //           the last node index == the first one.
 917 //=======================================================================
 918
 919 const int* SMDS_VolumeTool::GetFaceNodesIndices( int faceIndex )
 920 {
 921   if ( !setFace( faceIndex ))
 922     return 0;
 923
 924   if (myVolume->IsPoly())
 925   {
 926     myPolyIndices.resize( myFaceNbNodes + 1 );
 927     myFaceNodeIndices = & myPolyIndices[0];
 928     for ( int i = 0; i <= myFaceNbNodes; ++i )
 929       myFaceNodeIndices[i] = myVolume->GetNodeIndex( myFaceNodes[i] );
 930   }
 931   return myFaceNodeIndices;
 932 }
 933
 934 //=======================================================================
 935 //function : GetFaceNodes
 936 //purpose  : Return a set of face nodes.
 937 //=======================================================================
 938
 939 bool SMDS_VolumeTool::GetFaceNodes (int                        faceIndex,
 940                                     set<const SMDS_MeshNode*>& theFaceNodes )
 941 {
 942   if ( !setFace( faceIndex ))
 943     return false;
 944
 945   theFaceNodes.clear();
 946   int iNode, nbNode = myFaceNbNodes;
 947   for ( iNode = 0; iNode < nbNode; iNode++ )
 948     theFaceNodes.insert( myFaceNodes[ iNode ]);
 949
 950   return true;
 951 }
 952
 953 //=======================================================================
 954 //function : IsFaceExternal
 955 //purpose  : Check normal orientation of a given face
 956 //=======================================================================
 957
 958 bool SMDS_VolumeTool::IsFaceExternal( int faceIndex )
 959 {
 960   if ( myExternalFaces || !myVolume )
 961     return true;
 962
 963   if (myVolume->IsPoly()) {
 964     XYZ aNormal, baryCenter, p0 (myPolyedre->GetFaceNode(faceIndex + 1, 1));
 965     GetFaceNormal(faceIndex, aNormal.x, aNormal.y, aNormal.z);
 966     GetBaryCenter(baryCenter.x, baryCenter.y, baryCenter.z);
 967     XYZ insideVec (baryCenter - p0);
 968     if (insideVec.Dot(aNormal) > 0)
 969       return false;
 970     return true;
 971   }
 972
 973   switch ( myVolumeNbNodes ) {
 974   case 4:
 975   case 5:
 976   case 10:
 977   case 13:
 978     // only the bottom of a reversed tetrahedron can be internal
 979     return ( myVolForward || faceIndex != 0 );
 980   case 6:
 981   case 15:
 982     // in a forward pentahedron, the top is internal, in a reversed one - bottom
 983     return ( myVolForward ? faceIndex != 1 : faceIndex != 0 );
 984   case 8:
 985   case 20: {
 986     // in a forward hexahedron, even face normal is external, odd - internal
 987     bool odd = faceIndex % 2;
 988     return ( myVolForward ? !odd : odd );
 989   }
 990   default:;
 991   }
 992   return false;
 993 }
 994
 995 //=======================================================================
 996 //function : GetFaceNormal
 997 //purpose  : Return a normal to a face
 998 //=======================================================================
 999
1000 bool SMDS_VolumeTool::GetFaceNormal (int faceIndex, double & X, double & Y, double & Z)
1001 {
1002   if ( !setFace( faceIndex ))
1003     return false;
1004
1005   XYZ p1 ( myFaceNodes[0] );
1006   XYZ p2 ( myFaceNodes[1] );
1007   XYZ p3 ( myFaceNodes[2] );
1008   XYZ aVec12( p2 - p1 );
1009   XYZ aVec13( p3 - p1 );
1010   XYZ cross = aVec12.Crossed( aVec13 );
1011
1012   //if ( myFaceNbNodes == 4 ) {
1013   if ( myFaceNbNodes >3 ) {
1014     XYZ p4 ( myFaceNodes[3] );
1015     XYZ aVec14( p4 - p1 );
1016     XYZ cross2 = aVec13.Crossed( aVec14 );
1017     cross.x += cross2.x;
1018     cross.y += cross2.y;
1019     cross.z += cross2.z;
1020   }
1021
1022   double size = cross.Magnitude();
1023   if ( size <= DBL_MIN )
1024     return false;
1025
1026   X = cross.x / size;
1027   Y = cross.y / size;
1028   Z = cross.z / size;
1029
1030   return true;
1031 }
1032
1033 //================================================================================
1034 /*!
1035  * \brief Return barycenter of a face
1036  */
1037 //================================================================================
1038
1039 bool SMDS_VolumeTool::GetFaceBaryCenter (int faceIndex, double & X, double & Y, double & Z)
1040 {
1041   if ( !setFace( faceIndex ))
1042     return false;
1043
1044   X = Y = Z = 0.0;
1045   for ( int i = 0; i < myFaceNbNodes; ++i )
1046   {
1047     X += myFaceNodes[i]->X() / myFaceNbNodes;
1048     Y += myFaceNodes[i]->Y() / myFaceNbNodes;
1049     Z += myFaceNodes[i]->Z() / myFaceNbNodes;
1050   }
1051   return true;
1052 }
1053
1054 //=======================================================================
1055 //function : GetFaceArea
1056 //purpose  : Return face area
1057 //=======================================================================
1058
1059 double SMDS_VolumeTool::GetFaceArea( int faceIndex )
1060 {
1061   if (myVolume->IsPoly()) {
1062     MESSAGE("Warning: attempt to obtain area of a face of polyhedral volume");
1063     return 0;
1064   }
1065
1066   if ( !setFace( faceIndex ))
1067     return 0;
1068
1069   XYZ p1 ( myFaceNodes[0] );
1070   XYZ p2 ( myFaceNodes[1] );
1071   XYZ p3 ( myFaceNodes[2] );
1072   XYZ aVec12( p2 - p1 );
1073   XYZ aVec13( p3 - p1 );
1074   double area = aVec12.Crossed( aVec13 ).Magnitude() * 0.5;
1075
1076   if ( myFaceNbNodes == 4 ) {
1077     XYZ p4 ( myFaceNodes[3] );
1078     XYZ aVec14( p4 - p1 );
1079     area += aVec14.Crossed( aVec13 ).Magnitude() * 0.5;
1080   }
1081   return area;
1082 }
1083
1084 //=======================================================================
1085 //function : GetOppFaceIndex
1086 //purpose  : Return index of the opposite face if it exists, else -1.
1087 //=======================================================================
1088
1089 int SMDS_VolumeTool::GetOppFaceIndex( int faceIndex ) const
1090 {
1091   int ind = -1;
1092   if (myVolume->IsPoly()) {
1093     MESSAGE("Warning: attempt to obtain opposite face on polyhedral volume");
1094     return ind;
1095   }
1096
1097   if ( faceIndex >= 0 && faceIndex < NbFaces() ) {
1098     switch ( myVolumeNbNodes ) {
1099     case 6:
1100       if ( faceIndex == 0 || faceIndex == 1 )
1101         ind = 1 - faceIndex;
1102         break;
1103     case 8:
1104       ind = faceIndex + ( faceIndex % 2 ? -1 : 1 );
1105       break;
1106     default:;
1107     }
1108   }
1109   return ind;
1110 }
1111
1112 //=======================================================================
1113 //function : IsLinked
1114 //purpose  : return true if theNode1 is linked with theNode2
1115 // If theIgnoreMediumNodes then corner nodes of quadratic cell are considered linked as well
1116 //=======================================================================
1117
1118 bool SMDS_VolumeTool::IsLinked (const SMDS_MeshNode* theNode1,
1119                                 const SMDS_MeshNode* theNode2,
1120                                 const bool           theIgnoreMediumNodes) const
1121 {
1122   if ( !myVolume )
1123     return false;
1124
1125   if (myVolume->IsPoly()) {
1126     if (!myPolyedre) {
1127       MESSAGE("Warning: bad volumic element");
1128       return false;
1129     }
1130     bool isLinked = false;
1131     int iface;
1132     for (iface = 1; iface <= myNbFaces && !isLinked; iface++) {
1133       int inode, nbFaceNodes = myPolyedre->NbFaceNodes(iface);
1134
1135       for (inode = 1; inode <= nbFaceNodes && !isLinked; inode++) {
1136         const SMDS_MeshNode* curNode = myPolyedre->GetFaceNode(iface, inode);
1137
1138         if (curNode == theNode1 || curNode == theNode2) {
1139           int inextnode = (inode == nbFaceNodes) ? 1 : inode + 1;
1140           const SMDS_MeshNode* nextNode = myPolyedre->GetFaceNode(iface, inextnode);
1141
1142           if ((curNode == theNode1 && nextNode == theNode2) ||
1143               (curNode == theNode2 && nextNode == theNode1)) {
1144             isLinked = true;
1145           }
1146         }
1147       }
1148     }
1149     return isLinked;
1150   }
1151
1152   // find nodes indices
1153   int i1 = -1, i2 = -1;
1154   for ( int i = 0; i < myVolumeNbNodes; i++ ) {
1155     if ( myVolumeNodes[ i ] == theNode1 )
1156       i1 = i;
1157     else if ( myVolumeNodes[ i ] == theNode2 )
1158       i2 = i;
1159   }
1160   return IsLinked( i1, i2 );
1161 }
1162
1163 //=======================================================================
1164 //function : IsLinked
1165 //purpose  : return true if the node with theNode1Index is linked
1166 //           with the node with theNode2Index
1167 // If theIgnoreMediumNodes then corner nodes of quadratic cell are considered linked as well
1168 //=======================================================================
1169
1170 bool SMDS_VolumeTool::IsLinked (const int theNode1Index,
1171                                 const int theNode2Index,
1172                                 bool      theIgnoreMediumNodes) const
1173 {
1174   if ( myVolume->IsPoly() ) {
1175     return IsLinked(myVolumeNodes[theNode1Index], myVolumeNodes[theNode2Index]);
1176   }
1177
1178   int minInd = min( theNode1Index, theNode2Index );
1179   int maxInd = max( theNode1Index, theNode2Index );
1180
1181   if ( minInd < 0 || maxInd > myVolumeNbNodes - 1 || maxInd == minInd )
1182     return false;
1183
1184   SMDSAbs_EntityType type = myVolume->GetEntityType();
1185   if ( myVolume->IsQuadratic() )
1186   {
1187     int firstMediumInd = myVolume->NbCornerNodes();
1188     if ( minInd >= firstMediumInd )
1189       return false; // medium nodes are not linked
1190     if ( maxInd < firstMediumInd ) // both nodes are corners
1191     {
1192       if ( theIgnoreMediumNodes )
1193         type = SMDSAbs_EntityType( int(type)-1 ); // check linkage of corner nodes
1194       else
1195         return false; // corner nodes are not linked directly in a quadratic cell
1196     }
1197   }
1198
1199   switch ( type ) {
1200   case SMDSEntity_Tetra:
1201     return true;
1202   case SMDSEntity_Hexa:
1203     switch ( maxInd - minInd ) {
1204     case 1: return minInd != 3;
1205     case 3: return minInd == 0 || minInd == 4;
1206     case 4: return true;
1207     default:;
1208     }
1209     break;
1210   case SMDSEntity_Pyramid:
1211     if ( maxInd == 4 )
1212       return true;
1213     switch ( maxInd - minInd ) {
1214     case 1:
1215     case 3: return true;
1216     default:;
1217     }
1218     break;
1219   case SMDSEntity_Penta:
1220     switch ( maxInd - minInd ) {
1221     case 1: return minInd != 2;
1222     case 2: return minInd == 0 || minInd == 3;
1223     case 3: return true;
1224     default:;
1225     }
1226     break;
1227   case SMDSEntity_Quad_Tetra:
1228     {
1229       switch ( minInd ) {
1230       case 0: if( maxInd==4 ||  maxInd==6 ||  maxInd==7 ) return true;
1231       case 1: if( maxInd==4 ||  maxInd==5 ||  maxInd==8 ) return true;
1232       case 2: if( maxInd==5 ||  maxInd==6 ||  maxInd==9 ) return true;
1233       case 3: if( maxInd==7 ||  maxInd==8 ||  maxInd==9 ) return true;
1234       default:;
1235       }
1236       break;
1237     }
1238   case SMDSEntity_Quad_Hexa:
1239     {
1240       switch ( minInd ) {
1241       case 0: if( maxInd==8 ||  maxInd==11 ||  maxInd==16 ) return true;
1242       case 1: if( maxInd==8 ||  maxInd==9 ||  maxInd==17 ) return true;
1243       case 2: if( maxInd==9 ||  maxInd==10 ||  maxInd==18 ) return true;
1244       case 3: if( maxInd==10 ||  maxInd==11 ||  maxInd==19 ) return true;
1245       case 4: if( maxInd==12 ||  maxInd==15 ||  maxInd==16 ) return true;
1246       case 5: if( maxInd==12 ||  maxInd==13 ||  maxInd==17 ) return true;
1247       case 6: if( maxInd==13 ||  maxInd==14 ||  maxInd==18 ) return true;
1248       case 7: if( maxInd==14 ||  maxInd==15 ||  maxInd==19 ) return true;
1249       default:;
1250       }
1251       break;
1252     }
1253   case SMDSEntity_Quad_Pyramid:
1254     {
1255       switch ( minInd ) {
1256       case 0: if( maxInd==5 ||  maxInd==8 ||  maxInd==9 ) return true;
1257       case 1: if( maxInd==5 ||  maxInd==6 ||  maxInd==10 ) return true;
1258       case 2: if( maxInd==6 ||  maxInd==7 ||  maxInd==11 ) return true;
1259       case 3: if( maxInd==7 ||  maxInd==8 ||  maxInd==12 ) return true;
1260       case 4: if( maxInd==9 ||  maxInd==10 ||  maxInd==11 ||  maxInd==12 ) return true;
1261       default:;
1262       }
1263       break;
1264     }
1265   case SMDSEntity_Quad_Penta:
1266     {
1267       switch ( minInd ) {
1268       case 0: if( maxInd==6 ||  maxInd==8 ||  maxInd==12 ) return true;
1269       case 1: if( maxInd==6 ||  maxInd==7 ||  maxInd==13 ) return true;
1270       case 2: if( maxInd==7 ||  maxInd==8 ||  maxInd==14 ) return true;
1271       case 3: if( maxInd==9 ||  maxInd==11 ||  maxInd==12 ) return true;
1272       case 4: if( maxInd==9 ||  maxInd==10 ||  maxInd==13 ) return true;
1273       case 5: if( maxInd==10 ||  maxInd==11 ||  maxInd==14 ) return true;
1274       default:;
1275       }
1276       break;
1277     }
1278   default:;
1279   }
1280   return false;
1281 }
1282
1283 //=======================================================================
1284 //function : GetNodeIndex
1285 //purpose  : Return an index of theNode
1286 //=======================================================================
1287
1288 int SMDS_VolumeTool::GetNodeIndex(const SMDS_MeshNode* theNode) const
1289 {
1290   if ( myVolume ) {
1291     for ( int i = 0; i < myVolumeNbNodes; i++ ) {
1292       if ( myVolumeNodes[ i ] == theNode )
1293         return i;
1294     }
1295   }
1296   return -1;
1297 }
1298
1299 //================================================================================
1300 /*!
1301  * \brief Fill vector with boundary faces existing in the mesh
1302   * \param faces - vector of found nodes
1303   * \retval int - nb of found faces
1304  */
1305 //================================================================================
1306
1307 int SMDS_VolumeTool::GetAllExistingFaces(vector<const SMDS_MeshElement*> & faces)
1308 {
1309   faces.clear();
1310   faces.reserve( NbFaces() );
1311   for ( int iF = 0; iF < NbFaces(); ++iF ) {
1312     const SMDS_MeshFace* face = 0;
1313     const SMDS_MeshNode** nodes = GetFaceNodes( iF );
1314     switch ( NbFaceNodes( iF )) {
1315     case 3:
1316       face = SMDS_Mesh::FindFace( nodes[0], nodes[1], nodes[2] ); break;
1317     case 4:
1318       face = SMDS_Mesh::FindFace( nodes[0], nodes[1], nodes[2], nodes[3] ); break;
1319     case 6:
1320       face = SMDS_Mesh::FindFace( nodes[0], nodes[1], nodes[2],
1321                                   nodes[3], nodes[4], nodes[5]); break;
1322     case 8:
1323       face = SMDS_Mesh::FindFace( nodes[0], nodes[1], nodes[2], nodes[3],
1324                                   nodes[4], nodes[5], nodes[6], nodes[7]); break;
1325     }
1326     if ( face )
1327       faces.push_back( face );
1328   }
1329   return faces.size();
1330 }
1331
1332
1333 //================================================================================
1334 /*!
1335  * \brief Fill vector with boundary edges existing in the mesh
1336   * \param edges - vector of found edges
1337   * \retval int - nb of found faces
1338  */
1339 //================================================================================
1340
1341 int SMDS_VolumeTool::GetAllExistingEdges(vector<const SMDS_MeshElement*> & edges) const
1342 {
1343   edges.clear();
1344   edges.reserve( myVolumeNbNodes * 2 );
1345   for ( int i = 0; i < myVolumeNbNodes; ++i ) {
1346     for ( int j = i + 1; j < myVolumeNbNodes; ++j ) {
1347       if ( IsLinked( i, j )) {
1348         const SMDS_MeshElement* edge =
1349           SMDS_Mesh::FindEdge( myVolumeNodes[i], myVolumeNodes[j] );
1350         if ( edge )
1351           edges.push_back( edge );
1352       }
1353     }
1354   }
1355   return edges.size();
1356 }
1357
1358 //================================================================================
1359 /*!
1360  * \brief check that only one volume is build on the face nodes
1361  *
1362  * If a face is shared by one of <ignoreVolumes>, it is considered free
1363  */
1364 //================================================================================
1365
1366 bool SMDS_VolumeTool::IsFreeFace( int faceIndex, const SMDS_MeshElement** otherVol/*=0*/ )
1367 {
1368   const bool isFree = true;
1369
1370   if (!setFace( faceIndex ))
1371     return !isFree;
1372
1373   const SMDS_MeshNode** nodes = GetFaceNodes( faceIndex );
1374   const int nbFaceNodes = myFaceNbNodes;
1375
1376   // evaluate nb of face nodes shared by other volume
1377   int maxNbShared = -1;
1378   typedef map< const SMDS_MeshElement*, int > TElemIntMap;
1379   TElemIntMap volNbShared;
1380   TElemIntMap::iterator vNbIt;
1381   for ( int iNode = 0; iNode < nbFaceNodes; iNode++ ) {
1382     const SMDS_MeshNode* n = nodes[ iNode ];
1383     SMDS_ElemIteratorPtr eIt = n->GetInverseElementIterator( SMDSAbs_Volume );
1384     while ( eIt->more() ) {
1385       const SMDS_MeshElement* elem = eIt->next();
1386       if ( elem != myVolume ) {
1387         vNbIt = volNbShared.insert( make_pair( elem, 0 )).first;
1388         (*vNbIt).second++;
1389         if ( vNbIt->second > maxNbShared )
1390           maxNbShared = vNbIt->second;
1391       }
1392     }
1393   }
1394   if ( maxNbShared < 3 )
1395     return isFree; // is free
1396
1397   // find volumes laying on the opposite side of the face
1398   // and sharing all nodes
1399   XYZ intNormal; // internal normal
1400   GetFaceNormal( faceIndex, intNormal.x, intNormal.y, intNormal.z );
1401   if ( IsFaceExternal( faceIndex ))
1402     intNormal = XYZ( -intNormal.x, -intNormal.y, -intNormal.z );
1403   XYZ p0 ( nodes[0] ), baryCenter;
1404   for ( vNbIt = volNbShared.begin(); vNbIt != volNbShared.end();  ) {
1405     const int& nbShared = (*vNbIt).second;
1406     if ( nbShared >= 3 ) {
1407       SMDS_VolumeTool volume( (*vNbIt).first );
1408       volume.GetBaryCenter( baryCenter.x, baryCenter.y, baryCenter.z );
1409       XYZ intNormal2( baryCenter - p0 );
1410       if ( intNormal.Dot( intNormal2 ) < 0 ) {
1411         // opposite side
1412         if ( nbShared >= nbFaceNodes )
1413         {
1414           // a volume shares the whole facet
1415           if ( otherVol ) *otherVol = vNbIt->first;
1416           return !isFree;
1417         }
1418         ++vNbIt;
1419         continue;
1420       }
1421     }
1422     // remove a volume from volNbShared map
1423     volNbShared.erase( vNbIt++ );
1424   }
1425
1426   // here volNbShared contains only volumes laying on the opposite side of
1427   // the face and sharing 3 or more but not all face nodes with myVolume
1428   if ( volNbShared.size() < 2 ) {
1429     return isFree; // is free
1430   }
1431
1432   // check if the whole area of a face is shared
1433   for ( int iNode = 0; iNode < nbFaceNodes; iNode++ )
1434   {
1435     const SMDS_MeshNode* n = nodes[ iNode ];
1436     // check if n is shared by one of volumes of volNbShared
1437     bool isShared = false;
1438     SMDS_ElemIteratorPtr eIt = n->GetInverseElementIterator( SMDSAbs_Volume );
1439     while ( eIt->more() && !isShared )
1440       isShared = volNbShared.count( eIt->next() );
1441     if ( !isShared )
1442       return isFree;
1443   }
1444   if ( otherVol ) *otherVol = volNbShared.begin()->first;
1445   return !isFree;
1446
1447 //   if ( !myVolume->IsPoly() )
1448 //   {
1449 //     bool isShared[] = { false, false, false, false }; // 4 triangle parts of a quadrangle
1450 //     for ( vNbIt = volNbShared.begin(); vNbIt != volNbShared.end(); vNbIt++ ) {
1451 //       SMDS_VolumeTool volume( (*vNbIt).first );
1452 //       bool prevLinkShared = false;
1453 //       int nbSharedLinks = 0;
1454 //       for ( int iNode = 0; iNode < nbFaceNodes; iNode++ ) {
1455 //         bool linkShared = volume.IsLinked( nodes[ iNode ], nodes[ iNode + 1] );
1456 //         if ( linkShared )
1457 //           nbSharedLinks++;
1458 //         if ( linkShared && prevLinkShared &&
1459 //              volume.IsLinked( nodes[ iNode - 1 ], nodes[ iNode + 1] ))
1460 //           isShared[ iNode ] = true;
1461 //         prevLinkShared = linkShared;
1462 //       }
1463 //       if ( nbSharedLinks == nbFaceNodes )
1464 //         return !free; // is not free
1465 //       if ( nbFaceNodes == 4 ) {
1466 //         // check traingle parts 1 & 3
1467 //         if ( isShared[1] && isShared[3] )
1468 //           return !free; // is not free
1469 //         // check triangle parts 0 & 2;
1470 //         // 0 part could not be checked in the loop; check it here
1471 //         if ( isShared[2] && prevLinkShared &&
1472 //              volume.IsLinked( nodes[ 0 ], nodes[ 1 ] ) &&
1473 //              volume.IsLinked( nodes[ 1 ], nodes[ 3 ] ) )
1474 //           return !free; // is not free
1475 //       }
1476 //     }
1477 //   }
1478 //  return free;
1479 }
1480
1481 //=======================================================================
1482 //function : GetFaceIndex
1483 //purpose  : Return index of a face formed by theFaceNodes
1484 //=======================================================================
1485
1486 int SMDS_VolumeTool::GetFaceIndex( const set<const SMDS_MeshNode*>& theFaceNodes )
1487 {
1488   for ( int iFace = 0; iFace < myNbFaces; iFace++ ) {
1489     const SMDS_MeshNode** nodes = GetFaceNodes( iFace );
1490     int nbFaceNodes = NbFaceNodes( iFace );
1491     set<const SMDS_MeshNode*> nodeSet;
1492     for ( int iNode = 0; iNode < nbFaceNodes; iNode++ )
1493       nodeSet.insert( nodes[ iNode ] );
1494     if ( theFaceNodes == nodeSet )
1495       return iFace;
1496   }
1497   return -1;
1498 }
1499
1500 //=======================================================================
1501 //function : GetFaceIndex
1502 //purpose  : Return index of a face formed by theFaceNodes
1503 //=======================================================================
1504
1505 /*int SMDS_VolumeTool::GetFaceIndex( const set<int>& theFaceNodesIndices )
1506 {
1507   for ( int iFace = 0; iFace < myNbFaces; iFace++ ) {
1508     const int* nodes = GetFaceNodesIndices( iFace );
1509     int nbFaceNodes = NbFaceNodes( iFace );
1510     set<int> nodeSet;
1511     for ( int iNode = 0; iNode < nbFaceNodes; iNode++ )
1512       nodeSet.insert( nodes[ iNode ] );
1513     if ( theFaceNodesIndices == nodeSet )
1514       return iFace;
1515   }
1516   return -1;
1517 }*/
1518
1519 //=======================================================================
1520 //function : setFace
1521 //purpose  :
1522 //=======================================================================
1523
1524 bool SMDS_VolumeTool::setFace( int faceIndex )
1525 {
1526   if ( !myVolume )
1527     return false;
1528
1529   if ( myCurFace == faceIndex )
1530     return true;
1531
1532   myCurFace = -1;
1533
1534   if ( faceIndex < 0 || faceIndex >= NbFaces() )
1535     return false;
1536
1537   if (myFaceNodes != NULL) {
1538     delete [] myFaceNodes;
1539     myFaceNodes = NULL;
1540   }
1541
1542   if (myVolume->IsPoly()) {
1543     if (!myPolyedre) {
1544       MESSAGE("Warning: bad volumic element");
1545       return false;
1546     }
1547
1548     // set face nodes
1549     int iNode;
1550     myFaceNbNodes = myPolyedre->NbFaceNodes(faceIndex + 1);
1551     myFaceNodes = new const SMDS_MeshNode* [myFaceNbNodes + 1];
1552     for ( iNode = 0; iNode < myFaceNbNodes; iNode++ )
1553       myFaceNodes[ iNode ] = myPolyedre->GetFaceNode(faceIndex + 1, iNode + 1);
1554     myFaceNodes[ myFaceNbNodes ] = myFaceNodes[ 0 ]; // last = first
1555
1556     // check orientation
1557     if (myExternalFaces)
1558     {
1559       myCurFace = faceIndex; // avoid infinite recursion in IsFaceExternal()
1560       myExternalFaces = false; // force normal computation by IsFaceExternal()
1561       if ( !IsFaceExternal( faceIndex ))
1562         for ( int i = 0, j = myFaceNbNodes; i < j; ++i, --j )
1563           std::swap( myFaceNodes[i], myFaceNodes[j] );
1564       myExternalFaces = true;
1565     }
1566   }
1567   else {
1568     // choose face node indices
1569     switch ( myVolumeNbNodes ) {
1570     case 4:
1571       myFaceNbNodes = Tetra_nbN[ faceIndex ];
1572       if ( myExternalFaces )
1573         myFaceNodeIndices = myVolForward ? Tetra_F[ faceIndex ] : Tetra_RE[ faceIndex ];
1574       else
1575         myFaceNodeIndices = myVolForward ? Tetra_F[ faceIndex ] : Tetra_R[ faceIndex ];
1576       break;
1577     case 5:
1578       myFaceNbNodes = Pyramid_nbN[ faceIndex ];
1579       if ( myExternalFaces )
1580         myFaceNodeIndices = myVolForward ? Pyramid_F[ faceIndex ] : Pyramid_RE[ faceIndex ];
1581       else
1582         myFaceNodeIndices = myVolForward ? Pyramid_F[ faceIndex ] : Pyramid_R[ faceIndex ];
1583       break;
1584     case 6:
1585       myFaceNbNodes = Penta_nbN[ faceIndex ];
1586       if ( myExternalFaces )
1587         myFaceNodeIndices = myVolForward ? Penta_FE[ faceIndex ] : Penta_RE[ faceIndex ];
1588       else
1589         myFaceNodeIndices = myVolForward ? Penta_F[ faceIndex ] : Penta_R[ faceIndex ];
1590       break;
1591     case 8:
1592       myFaceNbNodes = Hexa_nbN[ faceIndex ];
1593       if ( myExternalFaces )
1594         myFaceNodeIndices = myVolForward ? Hexa_FE[ faceIndex ] : Hexa_RE[ faceIndex ];
1595       else
1596         myFaceNodeIndices = Hexa_F[ faceIndex ];
1597       break;
1598     case 10:
1599       myFaceNbNodes = QuadTetra_nbN[ faceIndex ];
1600       if ( myExternalFaces )
1601         myFaceNodeIndices = myVolForward ? QuadTetra_F[ faceIndex ] : QuadTetra_RE[ faceIndex ];
1602       else
1603         myFaceNodeIndices = myVolForward ? QuadTetra_F[ faceIndex ] : QuadTetra_R[ faceIndex ];
1604       break;
1605     case 13:
1606       myFaceNbNodes = QuadPyram_nbN[ faceIndex ];
1607       if ( myExternalFaces )
1608         myFaceNodeIndices = myVolForward ? QuadPyram_F[ faceIndex ] : QuadPyram_RE[ faceIndex ];
1609       else
1610         myFaceNodeIndices = myVolForward ? QuadPyram_F[ faceIndex ] : QuadPyram_R[ faceIndex ];
1611       break;
1612     case 15:
1613       myFaceNbNodes = QuadPenta_nbN[ faceIndex ];
1614       if ( myExternalFaces )
1615         myFaceNodeIndices = myVolForward ? QuadPenta_FE[ faceIndex ] : QuadPenta_RE[ faceIndex ];
1616       else
1617         myFaceNodeIndices = myVolForward ? QuadPenta_F[ faceIndex ] : QuadPenta_R[ faceIndex ];
1618       break;
1619     case 20:
1620       myFaceNbNodes = QuadHexa_nbN[ faceIndex ];
1621       if ( myExternalFaces )
1622         myFaceNodeIndices = myVolForward ? QuadHexa_FE[ faceIndex ] : QuadHexa_RE[ faceIndex ];
1623       else
1624         myFaceNodeIndices = QuadHexa_F[ faceIndex ];
1625       break;
1626     default:
1627       return false;
1628     }
1629
1630     // set face nodes
1631     myFaceNodes = new const SMDS_MeshNode* [myFaceNbNodes + 1];
1632     for ( int iNode = 0; iNode < myFaceNbNodes; iNode++ )
1633       myFaceNodes[ iNode ] = myVolumeNodes[ myFaceNodeIndices[ iNode ]];
1634     myFaceNodes[ myFaceNbNodes ] = myFaceNodes[ 0 ];
1635   }
1636
1637   myCurFace = faceIndex;
1638
1639   return true;
1640 }
1641
1642 //=======================================================================
1643 //function : GetType
1644 //purpose  : return VolumeType by nb of nodes in a volume
1645 //=======================================================================
1646
1647 SMDS_VolumeTool::VolumeType SMDS_VolumeTool::GetType(int nbNodes)
1648 {
1649   switch ( nbNodes ) {
1650   case 4: return TETRA;
1651   case 5: return PYRAM;
1652   case 6: return PENTA;
1653   case 8: return HEXA;
1654   case 10: return QUAD_TETRA;
1655   case 13: return QUAD_PYRAM;
1656   case 15: return QUAD_PENTA;
1657   case 20: return QUAD_HEXA;
1658   default:return UNKNOWN;
1659   }
1660 }
1661
1662 //=======================================================================
1663 //function : NbFaces
1664 //purpose  : return nb of faces by volume type
1665 //=======================================================================
1666
1667 int SMDS_VolumeTool::NbFaces( VolumeType type )
1668 {
1669   switch ( type ) {
1670   case TETRA     :
1671   case QUAD_TETRA: return 4;
1672   case PYRAM     :
1673   case QUAD_PYRAM: return 5;
1674   case PENTA     :
1675   case QUAD_PENTA: return 5;
1676   case HEXA      :
1677   case QUAD_HEXA : return 6;
1678   default:    return 0;
1679   }
1680 }
1681
1682 //================================================================================
1683 /*!
1684  * \brief Useful to know nb of corner nodes of a quadratic volume
1685   * \param type - volume type
1686   * \retval int - nb of corner nodes
1687  */
1688 //================================================================================
1689
1690 int SMDS_VolumeTool::NbCornerNodes(VolumeType type)
1691 {
1692   switch ( type ) {
1693   case TETRA     :
1694   case QUAD_TETRA: return 4;
1695   case PYRAM     :
1696   case QUAD_PYRAM: return 5;
1697   case PENTA     :
1698   case QUAD_PENTA: return 6;
1699   case HEXA      :
1700   case QUAD_HEXA : return 8;
1701   default:    return 0;
1702   }
1703   return 0;
1704 }
1705   //
1706
1707 //=======================================================================
1708 //function : GetFaceNodesIndices
1709 //purpose  : Return the array of face nodes indices
1710 //           To comfort link iteration, the array
1711 //           length == NbFaceNodes( faceIndex ) + 1 and
1712 //           the last node index == the first one.
1713 //=======================================================================
1714
1715 const int* SMDS_VolumeTool::GetFaceNodesIndices(VolumeType type,
1716                                                 int        faceIndex,
1717                                                 bool       external)
1718 {
1719   switch ( type ) {
1720   case TETRA: return Tetra_F[ faceIndex ];
1721   case PYRAM: return Pyramid_F[ faceIndex ];
1722   case PENTA: return external ? Penta_FE[ faceIndex ] : Penta_F[ faceIndex ];
1723   case HEXA:  return external ? Hexa_FE[ faceIndex ] : Hexa_F[ faceIndex ];
1724   case QUAD_TETRA: return QuadTetra_F[ faceIndex ];
1725   case QUAD_PYRAM: return QuadPyram_F[ faceIndex ];
1726   case QUAD_PENTA: return external ? QuadPenta_FE[ faceIndex ] : QuadPenta_F[ faceIndex ];
1727   case QUAD_HEXA:  return external ? QuadHexa_FE[ faceIndex ] : QuadHexa_F[ faceIndex ];
1728   default:;
1729   }
1730   return 0;
1731 }
1732
1733 //=======================================================================
1734 //function : NbFaceNodes
1735 //purpose  : Return number of nodes in the array of face nodes
1736 //=======================================================================
1737
1738 int SMDS_VolumeTool::NbFaceNodes(VolumeType type,
1739                                  int        faceIndex )
1740 {
1741   switch ( type ) {
1742   case TETRA: return Tetra_nbN[ faceIndex ];
1743   case PYRAM: return Pyramid_nbN[ faceIndex ];
1744   case PENTA: return Penta_nbN[ faceIndex ];
1745   case HEXA:  return Hexa_nbN[ faceIndex ];
1746   case QUAD_TETRA: return QuadTetra_nbN[ faceIndex ];
1747   case QUAD_PYRAM: return QuadPyram_nbN[ faceIndex ];
1748   case QUAD_PENTA: return QuadPenta_nbN[ faceIndex ];
1749   case QUAD_HEXA:  return QuadHexa_nbN[ faceIndex ];
1750   default:;
1751   }
1752   return 0;
1753 }
1754
1755 //=======================================================================
1756 //function : Get
1757 //purpose  : return element
1758 //=======================================================================
1759
1760 const SMDS_MeshVolume* SMDS_VolumeTool::Element() const
1761 {
1762   return static_cast<const SMDS_MeshVolume*>( myVolume );
1763 }
1764
1765 //=======================================================================
1766 //function : ID
1767 //purpose  : return element ID
1768 //=======================================================================
1769
1770 int SMDS_VolumeTool::ID() const
1771 {
1772   return myVolume ? myVolume->GetID() : 0;
1773 }