Salome HOME
Merge from V6_main 01/04/2013
[modules/smesh.git] / src / SMDS / SMDS_VtkVolume.cxx
1 // Copyright (C) 2010-2013  CEA/DEN, EDF R&D, OPEN CASCADE
2 //
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License.
7 //
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
11 // Lesser General Public License for more details.
12 //
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
16 //
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
18 //
19
20 #include "SMDS_VtkVolume.hxx"
21 #include "SMDS_MeshNode.hxx"
22 #include "SMDS_Mesh.hxx"
23 #include "SMDS_VtkCellIterator.hxx"
24
25 #include "utilities.h"
26
27 #include <vector>
28
29 SMDS_VtkVolume::SMDS_VtkVolume()
30 {
31 }
32
33 SMDS_VtkVolume::SMDS_VtkVolume(const std::vector<vtkIdType>& nodeIds, SMDS_Mesh* mesh)
34 {
35   init(nodeIds, mesh);
36 }
37 /*!
38  * typed used are vtk types (@see vtkCellType.h)
39  * see GetEntityType() for conversion in SMDS type (@see SMDSAbs_ElementType.hxx)
40  */
41 void SMDS_VtkVolume::init(const std::vector<vtkIdType>& nodeIds, SMDS_Mesh* mesh)
42 {
43   SMDS_MeshVolume::init();
44   vtkUnstructuredGrid* grid = mesh->getGrid();
45   myMeshId = mesh->getMeshId();
46   vtkIdType aType = VTK_TETRA;
47   switch (nodeIds.size()) // cases are in order of usage frequency
48   {
49     case 4:
50       aType = VTK_TETRA;
51       break;
52     case 8:
53       aType = VTK_HEXAHEDRON;
54       break;
55     case 5:
56       aType = VTK_PYRAMID;
57       break;
58     case 6:
59       aType = VTK_WEDGE;
60       break;
61     case 10:
62       aType = VTK_QUADRATIC_TETRA;
63       break;
64     case 20:
65       aType = VTK_QUADRATIC_HEXAHEDRON;
66       break;
67     case 13:
68       aType = VTK_QUADRATIC_PYRAMID;
69       break;
70     case 15:
71       aType = VTK_QUADRATIC_WEDGE;
72       break;
73     case 12:
74       aType = VTK_HEXAGONAL_PRISM;
75       break;
76     case 27:
77       aType = VTK_TRIQUADRATIC_HEXAHEDRON;
78       break;
79     default:
80       aType = VTK_HEXAHEDRON;
81       break;
82   }
83   myVtkID = grid->InsertNextLinkedCell(aType, nodeIds.size(), (vtkIdType *) &nodeIds[0]);
84   mesh->setMyModified();
85   //MESSAGE("SMDS_VtkVolume::init myVtkID " << myVtkID);
86 }
87
88 //#ifdef VTK_HAVE_POLYHEDRON
89 void SMDS_VtkVolume::initPoly(const std::vector<vtkIdType>& nodeIds,
90                               const std::vector<int>&       nbNodesPerFace,
91                               SMDS_Mesh*                    mesh)
92 {
93   SMDS_MeshVolume::init();
94   //MESSAGE("SMDS_VtkVolume::initPoly");
95   SMDS_UnstructuredGrid* grid = mesh->getGrid();
96   //double center[3];
97   //this->gravityCenter(grid, &nodeIds[0], nodeIds.size(), &center[0]);
98   vector<vtkIdType> ptIds;
99   vtkIdType nbFaces = nbNodesPerFace.size();
100   int k = 0;
101   for (int i = 0; i < nbFaces; i++)
102     {
103       int nf = nbNodesPerFace[i];
104       ptIds.push_back(nf);
105       // EAP: a right approach is:
106       // - either the user should care of order of nodes or
107       // - the user should use a service method arranging nodes if he
108       //   don't want or can't to do it by him-self
109       // The method below works OK only with planar faces and convex polyhedrones
110       //
111       // double a[3];
112       // double b[3];
113       // double c[3];
114       // grid->GetPoints()->GetPoint(nodeIds[k], a);
115       // grid->GetPoints()->GetPoint(nodeIds[k + 1], b);
116       // grid->GetPoints()->GetPoint(nodeIds[k + 2], c);
117       // bool isFaceForward = this->isForward(a, b, c, center);
118       //MESSAGE("isFaceForward " << i << " " << isFaceForward);
119       const vtkIdType *facePts = &nodeIds[k];
120       //if (isFaceForward)
121         for (int n = 0; n < nf; n++)
122           ptIds.push_back(facePts[n]);
123       // else
124       //   for (int n = nf - 1; n >= 0; n--)
125       //     ptIds.push_back(facePts[n]);
126       k += nf;
127     }
128   myVtkID = grid->InsertNextLinkedCell(VTK_POLYHEDRON, nbFaces, &ptIds[0]);
129   mesh->setMyModified();
130 }
131 //#endif
132
133 bool SMDS_VtkVolume::ChangeNodes(const SMDS_MeshNode* nodes[], const int nbNodes)
134 {
135   vtkUnstructuredGrid* grid = SMDS_Mesh::_meshList[myMeshId]->getGrid();
136   vtkIdType npts = 0;
137   vtkIdType* pts = 0;
138   grid->GetCellPoints(myVtkID, npts, pts);
139   if (nbNodes != npts)
140     {
141       MESSAGE("ChangeNodes problem: not the same number of nodes " << npts << " -> " << nbNodes);
142       return false;
143     }
144   for (int i = 0; i < nbNodes; i++)
145     {
146       pts[i] = nodes[i]->getVtkId();
147     }
148   SMDS_Mesh::_meshList[myMeshId]->setMyModified();
149   return true;
150 }
151
152 /*!
153  * Reorder in VTK order a list of nodes given in SMDS order.
154  * To be used before ChangeNodes: lists are given or computed in SMDS order.
155  */
156 bool SMDS_VtkVolume::vtkOrder(const SMDS_MeshNode* nodes[], const int nbNodes)
157 {
158   if (nbNodes != this->NbNodes())
159     {
160       MESSAGE("vtkOrder, wrong number of nodes " << nbNodes << " instead of "<< this->NbNodes());
161       return false;
162     }
163   vtkUnstructuredGrid* grid = SMDS_Mesh::_meshList[myMeshId]->getGrid();
164   vtkIdType aVtkType = grid->GetCellType(this->myVtkID);
165   const std::vector<int>& interlace = SMDS_MeshCell::toVtkOrder( VTKCellType( aVtkType ));
166   if ( !interlace.empty() )
167   {
168     ASSERT( interlace.size() == nbNodes );
169     std::vector<const SMDS_MeshNode*> initNodes( nodes, nodes+nbNodes );
170     for ( size_t i = 0; i < interlace.size(); ++i )
171       nodes[i] = initNodes[ interlace[i] ];
172   }
173   return true;
174 }
175
176 SMDS_VtkVolume::~SMDS_VtkVolume()
177 {
178 }
179
180 void SMDS_VtkVolume::Print(ostream & OS) const
181 {
182   OS << "volume <" << GetID() << "> : ";
183 }
184
185 int SMDS_VtkVolume::NbFaces() const
186 {
187   vtkUnstructuredGrid* grid = SMDS_Mesh::_meshList[myMeshId]->getGrid();
188   vtkIdType aVtkType = grid->GetCellType(this->myVtkID);
189   int nbFaces = 4;
190   switch (aVtkType)
191   {
192     case VTK_TETRA:
193     case VTK_QUADRATIC_TETRA:
194       nbFaces = 4;
195       break;
196     case VTK_PYRAMID:
197     case VTK_WEDGE:
198     case VTK_QUADRATIC_PYRAMID:
199     case VTK_QUADRATIC_WEDGE:
200       nbFaces = 5;
201       break;
202     case VTK_HEXAHEDRON:
203     case VTK_QUADRATIC_HEXAHEDRON:
204     case VTK_TRIQUADRATIC_HEXAHEDRON:
205       nbFaces = 6;
206       break;
207     case VTK_POLYHEDRON:
208       {
209         vtkIdType nFaces = 0;
210         vtkIdType* ptIds = 0;
211         grid->GetFaceStream(this->myVtkID, nFaces, ptIds);
212         nbFaces = nFaces;
213         break;
214       }
215     case VTK_HEXAGONAL_PRISM:
216       nbFaces = 8;
217       break;
218     default:
219       MESSAGE("invalid volume type")
220       ;
221       nbFaces = 0;
222       break;
223   }
224   return nbFaces;
225 }
226
227 int SMDS_VtkVolume::NbNodes() const
228 {
229   vtkUnstructuredGrid* grid = SMDS_Mesh::_meshList[myMeshId]->getGrid();
230   vtkIdType aVtkType = grid->GetCellType(this->myVtkID);
231   int nbPoints = 0;
232   if (aVtkType != VTK_POLYHEDRON)
233     {
234       nbPoints = grid->GetCell(myVtkID)->GetNumberOfPoints();
235     }
236   else
237     {
238       vtkIdType nFaces = 0;
239       vtkIdType* ptIds = 0;
240       grid->GetFaceStream(this->myVtkID, nFaces, ptIds);
241       int id = 0;
242       for (int i = 0; i < nFaces; i++)
243         {
244           int nodesInFace = ptIds[id];
245           nbPoints += nodesInFace;
246           id += (nodesInFace + 1);
247         }
248     }
249   return nbPoints;
250 }
251
252 int SMDS_VtkVolume::NbEdges() const
253 {
254   vtkUnstructuredGrid* grid = SMDS_Mesh::_meshList[myMeshId]->getGrid();
255   vtkIdType aVtkType = grid->GetCellType(this->myVtkID);
256   int nbEdges = 6;
257   switch (aVtkType)
258   {
259     case VTK_TETRA:
260     case VTK_QUADRATIC_TETRA:
261       nbEdges = 6;
262       break;
263     case VTK_PYRAMID:
264     case VTK_QUADRATIC_PYRAMID:
265       nbEdges = 8;
266       break;
267     case VTK_WEDGE:
268     case VTK_QUADRATIC_WEDGE:
269       nbEdges = 9;
270       break;
271     case VTK_HEXAHEDRON:
272     case VTK_QUADRATIC_HEXAHEDRON:
273     case VTK_TRIQUADRATIC_HEXAHEDRON:
274       nbEdges = 12;
275       break;
276     case VTK_POLYHEDRON:
277       {
278         vtkIdType nFaces = 0;
279         vtkIdType* ptIds = 0;
280         grid->GetFaceStream(this->myVtkID, nFaces, ptIds);
281         nbEdges = 0;
282         int id = 0;
283         for (int i = 0; i < nFaces; i++)
284           {
285             int edgesInFace = ptIds[id];
286             id += (edgesInFace + 1);
287             nbEdges += edgesInFace;
288           }
289         nbEdges = nbEdges / 2;
290         break;
291       }
292     case VTK_HEXAGONAL_PRISM:
293       nbEdges = 18;
294       break;
295     default:
296       MESSAGE("invalid volume type")
297       ;
298       nbEdges = 0;
299       break;
300   }
301   return nbEdges;
302 }
303
304 /*! polyhedron only,
305  *  1 <= face_ind <= NbFaces()
306  */
307 int SMDS_VtkVolume::NbFaceNodes(const int face_ind) const
308 {
309   vtkUnstructuredGrid* grid = SMDS_Mesh::_meshList[myMeshId]->getGrid();
310   vtkIdType aVtkType = grid->GetCellType(this->myVtkID);
311   int nbNodes = 0;
312   if (aVtkType == VTK_POLYHEDRON)
313     {
314       vtkIdType nFaces = 0;
315       vtkIdType* ptIds = 0;
316       grid->GetFaceStream(this->myVtkID, nFaces, ptIds);
317       int id = 0;
318       for (int i = 0; i < nFaces; i++)
319         {
320           int nodesInFace = ptIds[id];
321           id += (nodesInFace + 1);
322           if (i == face_ind - 1)
323             {
324               nbNodes = nodesInFace;
325               break;
326             }
327         }
328     }
329   return nbNodes;
330 }
331
332 /*! polyhedron only,
333  *  1 <= face_ind <= NbFaces()
334  *  1 <= node_ind <= NbFaceNodes()
335  */
336 const SMDS_MeshNode* SMDS_VtkVolume::GetFaceNode(const int face_ind, const int node_ind) const
337 {
338   SMDS_Mesh *mesh = SMDS_Mesh::_meshList[myMeshId];
339   vtkUnstructuredGrid* grid = mesh->getGrid();
340   vtkIdType aVtkType = grid->GetCellType(this->myVtkID);
341   const SMDS_MeshNode* node = 0;
342   if (aVtkType == VTK_POLYHEDRON)
343     {
344       vtkIdType nFaces = 0;
345       vtkIdType* ptIds = 0;
346       grid->GetFaceStream(this->myVtkID, nFaces, ptIds);
347       int id = 0;
348       for (int i = 0; i < nFaces; i++)
349         {
350           int nodesInFace = ptIds[id]; // nodeIds in ptIds[id+1 .. id+nodesInFace]
351           if (i == face_ind - 1) // first face is number 1
352             {
353               if ((node_ind > 0) && (node_ind <= nodesInFace))
354                 node = mesh->FindNodeVtk(ptIds[id + node_ind]); // ptIds[id+1] : first node
355               break;
356             }
357           id += (nodesInFace + 1);
358         }
359     }
360   return node;
361 }
362
363 /*! polyhedron only,
364  *  return number of nodes for each face
365  */
366 std::vector<int> SMDS_VtkVolume::GetQuantities() const
367 {
368   vector<int> quantities;
369   quantities.clear();
370   SMDS_Mesh *mesh = SMDS_Mesh::_meshList[myMeshId];
371   vtkUnstructuredGrid* grid = mesh->getGrid();
372   vtkIdType aVtkType = grid->GetCellType(this->myVtkID);
373   if (aVtkType == VTK_POLYHEDRON)
374     {
375       vtkIdType nFaces = 0;
376       vtkIdType* ptIds = 0;
377       grid->GetFaceStream(this->myVtkID, nFaces, ptIds);
378       int id = 0;
379       for (int i = 0; i < nFaces; i++)
380         {
381           int nodesInFace = ptIds[id]; // nodeIds in ptIds[id+1 .. id+nodesInFace]
382           quantities.push_back(nodesInFace);
383           id += (nodesInFace + 1);
384         }
385     }
386   return quantities;
387 }
388
389 SMDS_ElemIteratorPtr SMDS_VtkVolume::elementsIterator(SMDSAbs_ElementType type) const
390 {
391   switch (type)
392   {
393     case SMDSAbs_Node:
394       {
395         SMDSAbs_EntityType aType = this->GetEntityType();
396         if (aType == SMDSEntity_Polyhedra)
397           return SMDS_ElemIteratorPtr(new SMDS_VtkCellIteratorPolyH(SMDS_Mesh::_meshList[myMeshId], myVtkID, aType));
398         else
399           return SMDS_ElemIteratorPtr(new SMDS_VtkCellIterator(SMDS_Mesh::_meshList[myMeshId], myVtkID, aType));
400       }
401     default:
402       MESSAGE("ERROR : Iterator not implemented");
403       return SMDS_ElemIteratorPtr((SMDS_ElemIterator*) NULL);
404   }
405 }
406
407 SMDS_ElemIteratorPtr SMDS_VtkVolume::nodesIteratorToUNV() const
408 {
409   return SMDS_ElemIteratorPtr(new SMDS_VtkCellIteratorToUNV(SMDS_Mesh::_meshList[myMeshId], myVtkID, GetEntityType()));
410 }
411
412 SMDS_ElemIteratorPtr SMDS_VtkVolume::interlacedNodesElemIterator() const
413 {
414   return SMDS_ElemIteratorPtr(new SMDS_VtkCellIteratorToUNV(SMDS_Mesh::_meshList[myMeshId], myVtkID, GetEntityType()));
415 }
416
417 SMDSAbs_ElementType SMDS_VtkVolume::GetType() const
418 {
419   return SMDSAbs_Volume;
420 }
421
422 /*!
423  * \brief Return node by its index
424  * \param ind - node index
425  * \retval const SMDS_MeshNode* - the node
426  */
427 const SMDS_MeshNode* SMDS_VtkVolume::GetNode(const int ind) const
428 {
429   vtkUnstructuredGrid* grid = SMDS_Mesh::_meshList[myMeshId]->getGrid();
430   vtkIdType aVtkType = grid->GetCellType(this->myVtkID);
431   vtkIdType npts, *pts;
432   grid->GetCellPoints( this->myVtkID, npts, pts );
433   const std::vector<int>& interlace = SMDS_MeshCell::fromVtkOrder( VTKCellType( aVtkType ));
434   return SMDS_Mesh::_meshList[myMeshId]->FindNodeVtk( pts[ interlace.empty() ? ind : interlace[ind]] );
435 }
436
437 bool SMDS_VtkVolume::IsQuadratic() const
438 {
439   vtkUnstructuredGrid* grid = SMDS_Mesh::_meshList[myMeshId]->getGrid();
440   vtkIdType aVtkType = grid->GetCellType(this->myVtkID);
441   // TODO quadratic polyhedrons ?
442   switch (aVtkType)
443   {
444     case VTK_QUADRATIC_TETRA:
445     case VTK_QUADRATIC_PYRAMID:
446     case VTK_QUADRATIC_WEDGE:
447     case VTK_QUADRATIC_HEXAHEDRON:
448     case VTK_TRIQUADRATIC_HEXAHEDRON:
449       return true;
450       break;
451     default:
452       return false;
453   }
454 }
455
456 bool SMDS_VtkVolume::IsPoly() const
457 {
458   vtkUnstructuredGrid* grid = SMDS_Mesh::_meshList[myMeshId]->getGrid();
459   vtkIdType aVtkType = grid->GetCellType(this->myVtkID);
460   return (aVtkType == VTK_POLYHEDRON);
461 }
462
463 bool SMDS_VtkVolume::IsMediumNode(const SMDS_MeshNode* node) const
464 {
465   vtkUnstructuredGrid* grid = SMDS_Mesh::_meshList[myMeshId]->getGrid();
466   vtkIdType aVtkType = grid->GetCellType(this->myVtkID);
467   int rankFirstMedium = 0;
468   switch (aVtkType)
469   {
470     case VTK_QUADRATIC_TETRA:
471       rankFirstMedium = 4; // medium nodes are of rank 4 to 9
472       break;
473     case VTK_QUADRATIC_PYRAMID:
474       rankFirstMedium = 5; // medium nodes are of rank 5 to 12
475       break;
476     case VTK_QUADRATIC_WEDGE:
477       rankFirstMedium = 6; // medium nodes are of rank 6 to 14
478       break;
479     case VTK_QUADRATIC_HEXAHEDRON:
480     case VTK_TRIQUADRATIC_HEXAHEDRON:
481       rankFirstMedium = 8; // medium nodes are of rank 8 to 19
482       break;
483     default:
484       return false;
485   }
486   vtkIdType npts = 0;
487   vtkIdType* pts = 0;
488   grid->GetCellPoints(myVtkID, npts, pts);
489   vtkIdType nodeId = node->getVtkId();
490   for (int rank = 0; rank < npts; rank++)
491     {
492       if (pts[rank] == nodeId)
493         {
494           if (rank < rankFirstMedium)
495             return false;
496           else
497             return true;
498         }
499     }
500   //throw SALOME_Exception(LOCALIZED("node does not belong to this element"));
501   MESSAGE("======================================================");
502   MESSAGE("= IsMediumNode: node does not belong to this element =");
503   MESSAGE("======================================================");
504   return false;
505 }
506
507 int SMDS_VtkVolume::NbCornerNodes() const
508 {
509   vtkUnstructuredGrid* grid = SMDS_Mesh::_meshList[myMeshId]->getGrid();
510   int            nbN = grid->GetCell(myVtkID)->GetNumberOfPoints();
511   vtkIdType aVtkType = grid->GetCellType(myVtkID);
512   switch (aVtkType)
513   {
514   case VTK_QUADRATIC_TETRA:         return 4;
515   case VTK_QUADRATIC_PYRAMID:       return 5;
516   case VTK_QUADRATIC_WEDGE:         return 6;
517   case VTK_QUADRATIC_HEXAHEDRON:
518   case VTK_TRIQUADRATIC_HEXAHEDRON: return 8;
519   default:;
520   }
521   return nbN;
522 }
523
524 SMDSAbs_EntityType SMDS_VtkVolume::GetEntityType() const
525 {
526   vtkUnstructuredGrid* grid = SMDS_Mesh::_meshList[myMeshId]->getGrid();
527   vtkIdType aVtkType = grid->GetCellType(this->myVtkID);
528
529   SMDSAbs_EntityType aType = SMDSEntity_Tetra;
530   switch (aVtkType)
531   {
532     case VTK_TETRA:
533       aType = SMDSEntity_Tetra;
534       break;
535     case VTK_PYRAMID:
536       aType = SMDSEntity_Pyramid;
537       break;
538     case VTK_WEDGE:
539       aType = SMDSEntity_Penta;
540       break;
541     case VTK_HEXAHEDRON:
542       aType = SMDSEntity_Hexa;
543       break;
544     case VTK_QUADRATIC_TETRA:
545       aType = SMDSEntity_Quad_Tetra;
546       break;
547     case VTK_QUADRATIC_PYRAMID:
548       aType = SMDSEntity_Quad_Pyramid;
549       break;
550     case VTK_QUADRATIC_WEDGE:
551       aType = SMDSEntity_Quad_Penta;
552       break;
553     case VTK_QUADRATIC_HEXAHEDRON:
554       aType = SMDSEntity_Quad_Hexa;
555       break;
556     case VTK_TRIQUADRATIC_HEXAHEDRON:
557       aType = SMDSEntity_TriQuad_Hexa;
558       break;
559     case VTK_HEXAGONAL_PRISM:
560       aType = SMDSEntity_Hexagonal_Prism;
561       break;
562 //#ifdef VTK_HAVE_POLYHEDRON
563     case VTK_POLYHEDRON:
564       aType = SMDSEntity_Polyhedra;
565       break;
566 //#endif
567     default:
568       aType = SMDSEntity_Polyhedra;
569       break;
570   }
571   return aType;
572 }
573
574 SMDSAbs_GeometryType SMDS_VtkVolume::GetGeomType() const
575 {
576   vtkUnstructuredGrid* grid = SMDS_Mesh::_meshList[myMeshId]->getGrid();
577   vtkIdType aVtkType = grid->GetCellType(this->myVtkID);
578
579   SMDSAbs_GeometryType aType = SMDSGeom_NONE;
580   switch (aVtkType)
581   {
582     case VTK_TETRA:
583     case VTK_QUADRATIC_TETRA:
584       aType = SMDSGeom_TETRA;
585       break;
586     case VTK_PYRAMID:
587     case VTK_QUADRATIC_PYRAMID:
588       aType = SMDSGeom_PYRAMID;
589       break;
590     case VTK_WEDGE:
591     case VTK_QUADRATIC_WEDGE:
592       aType = SMDSGeom_PENTA;
593       break;
594     case VTK_HEXAHEDRON:
595     case VTK_QUADRATIC_HEXAHEDRON:
596     case VTK_TRIQUADRATIC_HEXAHEDRON:
597       aType = SMDSGeom_HEXA;
598       break;
599     case VTK_HEXAGONAL_PRISM:
600       aType = SMDSGeom_HEXAGONAL_PRISM;
601       break;
602 //#ifdef VTK_HAVE_POLYHEDRON
603     case VTK_POLYHEDRON:
604       aType = SMDSGeom_POLYHEDRA;
605       break;
606 //#endif
607     default:
608       aType = SMDSGeom_POLYHEDRA;
609       break;
610   }
611   return aType;
612 }
613
614 vtkIdType SMDS_VtkVolume::GetVtkType() const
615 {
616   vtkUnstructuredGrid* grid = SMDS_Mesh::_meshList[myMeshId]->getGrid();
617   vtkIdType aType = grid->GetCellType(myVtkID);
618   return aType;
619 }
620
621 void SMDS_VtkVolume::gravityCenter(SMDS_UnstructuredGrid* grid,
622                                    const vtkIdType *      nodeIds,
623                                    int                    nbNodes,
624                                    double*                result)
625 {
626   for (int j = 0; j < 3; j++)
627     result[j] = 0;
628   if (nbNodes <= 0)
629     return;
630   for (int i = 0; i < nbNodes; i++)
631     {
632       double *coords = grid->GetPoint(nodeIds[i]);
633       for (int j = 0; j < 3; j++)
634         result[j] += coords[j];
635     }
636   for (int j = 0; j < 3; j++)
637     result[j] = result[j] / nbNodes;
638   //MESSAGE("center " << result[0] << " " << result[1] << " "  << result[2]);
639   return;
640 }
641
642 bool SMDS_VtkVolume::isForward(double* a, double* b, double* c, double* d)
643 {
644   double u[3], v[3], w[3];
645   for (int j = 0; j < 3; j++)
646     {
647       //MESSAGE("a,b,c,d " << a[j] << " " << b[j] << " " << c[j] << " " << d[j]);
648       u[j] = b[j] - a[j];
649       v[j] = c[j] - a[j];
650       w[j] = d[j] - a[j];
651       //MESSAGE("u,v,w " << u[j] << " " << v[j] << " " << w[j]);
652     }
653   double prodmixte = (u[1]*v[2] - u[2]*v[1]) * w[0]
654                    + (u[2]*v[0] - u[0]*v[2]) * w[1]
655                    + (u[0]*v[1] - u[1]*v[0]) * w[2];
656   return (prodmixte < 0);
657 }
658
659 /*! For polyhedron only
660  *  @return actual number of nodes (not the sum of nodes of all faces)
661  */
662 int SMDS_VtkVolume::NbUniqueNodes() const
663 {
664   vtkUnstructuredGrid* grid = SMDS_Mesh::_meshList[myMeshId]->getGrid();
665   return grid->GetCell(myVtkID)->GetNumberOfPoints();
666 }
667
668 /*! For polyhedron use only
669  *  @return iterator on actual nodes (not through the faces)
670  */
671 SMDS_ElemIteratorPtr SMDS_VtkVolume::uniqueNodesIterator() const
672 {
673   MESSAGE("uniqueNodesIterator");
674   return SMDS_ElemIteratorPtr(new SMDS_VtkCellIterator(SMDS_Mesh::_meshList[myMeshId], myVtkID, GetEntityType()));
675 }