Salome HOME
PR: correct some warnings
[modules/smesh.git] / src / SMDS / SMDS_Mesh.cxx
1 //  Copyright (C) 2007-2010  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 //  SMESH SMDS : implementation of Salome mesh data structure
24 //
25 #ifdef _MSC_VER
26 #pragma warning(disable:4786)
27 #endif
28
29 #include "utilities.h"
30 #include "SMDS_Mesh.hxx"
31 #include "SMDS_VolumeOfNodes.hxx"
32 #include "SMDS_VolumeOfFaces.hxx"
33 #include "SMDS_FaceOfNodes.hxx"
34 #include "SMDS_FaceOfEdges.hxx"
35 #include "SMDS_PolyhedralVolumeOfNodes.hxx"
36 #include "SMDS_PolygonalFaceOfNodes.hxx"
37 #include "SMDS_QuadraticEdge.hxx"
38 #include "SMDS_QuadraticFaceOfNodes.hxx"
39 #include "SMDS_QuadraticVolumeOfNodes.hxx"
40 #include "SMDS_SpacePosition.hxx"
41 #include "SMDS_UnstructuredGrid.hxx"
42
43 #include <vtkUnstructuredGrid.h>
44 #include <vtkUnstructuredGridWriter.h>
45 #include <vtkUnsignedCharArray.h>
46 #include <vtkCell.h>
47 #include <vtkCellLinks.h>
48 #include <vtkIdList.h>
49
50 #include <algorithm>
51 #include <map>
52 #include <iostream>
53 #include <fstream>
54 using namespace std;
55
56 #ifndef WIN32
57 #include <sys/sysinfo.h>
58 #endif
59
60 // number of added entities to check memory after
61 #define CHECKMEMORY_INTERVAL 1000
62
63 vector<SMDS_Mesh*> SMDS_Mesh::_meshList = vector<SMDS_Mesh*>();
64 int SMDS_Mesh::chunkSize = 1024;
65
66
67 //================================================================================
68 /*!
69  * \brief Raise an exception if free memory (ram+swap) too low
70  * \param doNotRaise - if true, suppres exception, just return free memory size
71  * \retval int - amount of available memory in MB or negative number in failure case
72  */
73 //================================================================================
74
75 int SMDS_Mesh::CheckMemory(const bool doNotRaise) throw (std::bad_alloc)
76 {
77 #ifndef WIN32
78   struct sysinfo si;
79   int err = sysinfo( &si );
80   if ( err )
81     return -1;
82
83   static int limit = -1;
84   if ( limit < 0 ) {
85     int status = system("SMDS_MemoryLimit"); // it returns lower limit of free RAM
86     if (status >= 0 ) {
87       limit = WEXITSTATUS(status);
88     }
89     if ( limit < 20 )
90       limit = 20;
91     else
92       limit = int( limit * 1.5 );
93 #ifdef _DEBUG_
94     MESSAGE ( "SMDS_Mesh::CheckMemory() memory limit = " << limit << " MB" );
95 #endif
96   }
97
98   const unsigned long Mbyte = 1024 * 1024;
99   // compute separately to avoid overflow
100   int freeMb =
101     ( si.freeram  * si.mem_unit ) / Mbyte +
102     ( si.freeswap * si.mem_unit ) / Mbyte;
103
104   if ( freeMb > limit )
105     return freeMb - limit;
106
107   if ( doNotRaise )
108     return 0;
109 #ifdef _DEBUG_
110   MESSAGE ("SMDS_Mesh::CheckMemory() throws as free memory too low: " << freeMb <<" MB" );
111 #endif
112   throw std::bad_alloc();
113 #else
114   return -1;
115 #endif
116 }
117
118 ///////////////////////////////////////////////////////////////////////////////
119 /// Create a new mesh object
120 ///////////////////////////////////////////////////////////////////////////////
121 SMDS_Mesh::SMDS_Mesh()
122         :myParent(NULL),
123          myNodeIDFactory(new SMDS_MeshNodeIDFactory()),
124          myElementIDFactory(new SMDS_MeshElementIDFactory()),
125          myHasConstructionEdges(false), myHasConstructionFaces(false),
126          myHasInverseElements(true),
127          myNodeMin(0), myNodeMax(0),
128          myNodePool(0), myEdgePool(0), myFacePool(0), myVolumePool(0),
129          myModified(false), myModifTime(0), myCompactTime(0),
130          xmin(0), xmax(0), ymin(0), ymax(0), zmin(0), zmax(0)
131 {
132   myMeshId = _meshList.size();         // --- index of the mesh to push back in the vector
133   MESSAGE("myMeshId=" << myMeshId);
134   MESSAGE("sizeof(SMDS_MeshElement) " << sizeof(SMDS_MeshElement) );
135   MESSAGE("sizeof(SMDS_MeshNode) " << sizeof(SMDS_MeshNode) );
136   MESSAGE("sizeof(SMDS_MeshCell) " << sizeof(SMDS_MeshCell) );
137   MESSAGE("sizeof(SMDS_VtkVolume) " << sizeof(SMDS_VtkVolume) );
138   MESSAGE("sizeof(SMDS_Position) " << sizeof(SMDS_Position) );
139   MESSAGE("sizeof(SMDS_SpacePosition) " << sizeof(SMDS_SpacePosition) );
140   myNodeIDFactory->SetMesh(this);
141   myElementIDFactory->SetMesh(this);
142   _meshList.push_back(this);
143   myNodePool = new ObjectPool<SMDS_MeshNode>(SMDS_Mesh::chunkSize);
144   myEdgePool = new ObjectPool<SMDS_VtkEdge>(SMDS_Mesh::chunkSize);
145   myFacePool = new ObjectPool<SMDS_VtkFace>(SMDS_Mesh::chunkSize);
146   myVolumePool = new ObjectPool<SMDS_VtkVolume>(SMDS_Mesh::chunkSize);
147
148   myNodes.clear();
149   myCells.clear();
150   //myCellIdSmdsToVtk.clear();
151   myCellIdVtkToSmds.clear();
152   myGrid = SMDS_UnstructuredGrid::New();
153   myGrid->setSMDS_mesh(this);
154   myGrid->Initialize();
155   myGrid->Allocate();
156   vtkPoints* points = vtkPoints::New();
157   points->SetNumberOfPoints(SMDS_Mesh::chunkSize);
158   myGrid->SetPoints( points );
159   points->Delete();
160   myGrid->BuildLinks();
161   this->Modified();
162 }
163
164 ///////////////////////////////////////////////////////////////////////////////
165 /// Create a new child mesh
166 /// Note that the tree structure of SMDS_Mesh seems to be unused in this version
167 /// (2003-09-08) of SMESH
168 ///////////////////////////////////////////////////////////////////////////////
169 SMDS_Mesh::SMDS_Mesh(SMDS_Mesh * parent)
170         :myParent(parent), myNodeIDFactory(parent->myNodeIDFactory),
171          myElementIDFactory(parent->myElementIDFactory),
172          myHasConstructionEdges(false), myHasConstructionFaces(false),
173          myHasInverseElements(true),
174          myNodePool(parent->myNodePool),
175          myEdgePool(parent->myEdgePool),
176          myFacePool(parent->myFacePool),
177          myVolumePool(parent->myVolumePool)
178 {
179 }
180
181 ///////////////////////////////////////////////////////////////////////////////
182 ///Create a submesh and add it to the current mesh
183 ///////////////////////////////////////////////////////////////////////////////
184
185 SMDS_Mesh *SMDS_Mesh::AddSubMesh()
186 {
187         SMDS_Mesh *submesh = new SMDS_Mesh(this);
188         myChildren.insert(myChildren.end(), submesh);
189         return submesh;
190 }
191
192 ///////////////////////////////////////////////////////////////////////////////
193 ///create a MeshNode and add it to the current Mesh
194 ///An ID is automatically assigned to the node.
195 ///@return : The created node
196 ///////////////////////////////////////////////////////////////////////////////
197
198 SMDS_MeshNode * SMDS_Mesh::AddNode(double x, double y, double z)
199 {
200   return SMDS_Mesh::AddNodeWithID(x,y,z,myNodeIDFactory->GetFreeID());
201 }
202
203 ///////////////////////////////////////////////////////////////////////////////
204 ///create a MeshNode and add it to the current Mesh
205 ///@param ID : The ID of the MeshNode to create
206 ///@return : The created node or NULL if a node with this ID already exists
207 ///////////////////////////////////////////////////////////////////////////////
208 SMDS_MeshNode * SMDS_Mesh::AddNodeWithID(double x, double y, double z, int ID)
209 {
210   // find the MeshNode corresponding to ID
211   const SMDS_MeshElement *node = myNodeIDFactory->MeshElement(ID);
212   if(!node){
213     // TODO ID < 1
214     if (ID <= 0)
215       {
216         MESSAGE("=============>  Bad Node Id: " << ID);
217         ID = myNodeIDFactory->GetFreeID();
218       }
219     myNodeIDFactory->adjustMaxId(ID);
220     SMDS_MeshNode * node = myNodePool->getNew();
221     node->init(ID, myMeshId, 0, x, y, z);
222     if (ID >= myNodes.size())
223     {
224         myNodes.resize(ID+SMDS_Mesh::chunkSize, 0);
225         MESSAGE(" ------------------ myNodes resize " << ID << " --> " << ID+SMDS_Mesh::chunkSize);
226     }
227     myNodes[ID] = node;
228     myNodeIDFactory->BindID(ID,node);
229     myInfo.myNbNodes++;
230     myModified = true;
231     this->adjustBoundingBox(x, y, z);
232     return node;
233   }else
234     return NULL;
235 }
236
237 ///////////////////////////////////////////////////////////////////////////////
238 /// create a Mesh0DElement and add it to the current Mesh
239 /// @return : The created Mesh0DElement
240 ///////////////////////////////////////////////////////////////////////////////
241 SMDS_Mesh0DElement* SMDS_Mesh::Add0DElementWithID(int idnode, int ID)
242 {
243   SMDS_MeshNode * node = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode);
244   if (!node) return NULL;
245   return SMDS_Mesh::Add0DElementWithID(node, ID);
246 }
247
248 ///////////////////////////////////////////////////////////////////////////////
249 /// create a Mesh0DElement and add it to the current Mesh
250 /// @return : The created Mesh0DElement
251 ///////////////////////////////////////////////////////////////////////////////
252 SMDS_Mesh0DElement* SMDS_Mesh::Add0DElement(const SMDS_MeshNode * node)
253 {
254   return SMDS_Mesh::Add0DElementWithID(node, myElementIDFactory->GetFreeID());
255 }
256
257 ///////////////////////////////////////////////////////////////////////////////
258 /// Create a new Mesh0DElement and at it to the mesh
259 /// @param idnode ID of the node
260 /// @param ID ID of the 0D element to create
261 /// @return The created 0D element or NULL if an element with this
262 ///         ID already exists or if input node is not found.
263 ///////////////////////////////////////////////////////////////////////////////
264 SMDS_Mesh0DElement* SMDS_Mesh::Add0DElementWithID(const SMDS_MeshNode * n, int ID)
265 {
266   if (!n) return 0;
267
268   //if (my0DElements.Extent() % CHECKMEMORY_INTERVAL == 0) CheckMemory();
269   //MESSAGE("Add0DElementWithID" << ID)
270   SMDS_Mesh0DElement * el0d = new SMDS_Mesh0DElement(n);
271   if (myElementIDFactory->BindID(ID, el0d)) {
272     //SMDS_MeshNode *node = const_cast<SMDS_MeshNode*>(n);
273     //node->AddInverseElement(el0d);// --- fait avec BindID
274     adjustmyCellsCapacity(ID);
275     myCells[ID] = el0d;
276     myInfo.myNb0DElements++;
277     return el0d;
278   }
279
280   delete el0d;
281   return NULL;
282 }
283
284 ///////////////////////////////////////////////////////////////////////////////
285 /// create a MeshEdge and add it to the current Mesh
286 /// @return : The created MeshEdge
287 ///////////////////////////////////////////////////////////////////////////////
288
289 SMDS_MeshEdge* SMDS_Mesh::AddEdgeWithID(int idnode1, int idnode2, int ID)
290 {
291   SMDS_MeshNode * node1 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode1);
292   SMDS_MeshNode * node2 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode2);
293   if(!node1 || !node2) return NULL;
294   return SMDS_Mesh::AddEdgeWithID(node1, node2, ID);
295 }
296
297 ///////////////////////////////////////////////////////////////////////////////
298 /// create a MeshEdge and add it to the current Mesh
299 /// @return : The created MeshEdge
300 ///////////////////////////////////////////////////////////////////////////////
301
302 SMDS_MeshEdge* SMDS_Mesh::AddEdge(const SMDS_MeshNode * node1,
303                                   const SMDS_MeshNode * node2)
304 {
305   return SMDS_Mesh::AddEdgeWithID(node1, node2, myElementIDFactory->GetFreeID());
306 }
307
308 ///////////////////////////////////////////////////////////////////////////////
309 /// Create a new edge and at it to the mesh
310 /// @param idnode1 ID of the first node
311 /// @param idnode2 ID of the second node
312 /// @param ID ID of the edge to create
313 /// @return The created edge or NULL if an element with this ID already exists or
314 /// if input nodes are not found.
315 ///////////////////////////////////////////////////////////////////////////////
316
317 SMDS_MeshEdge* SMDS_Mesh::AddEdgeWithID(const SMDS_MeshNode * n1,
318                                         const SMDS_MeshNode * n2,
319                                         int ID)
320 {
321   if ( !n1 || !n2 ) return 0;
322   SMDS_MeshEdge * edge = 0;
323
324   // --- retreive nodes ID
325   vector<vtkIdType> nodeIds;
326   nodeIds.clear();
327   nodeIds.push_back(n1->getVtkId());
328   nodeIds.push_back(n2->getVtkId());
329
330   SMDS_VtkEdge *edgevtk = myEdgePool->getNew();
331   edgevtk->init(nodeIds, this);
332   if (!this->registerElement(ID,edgevtk))
333     {
334       this->myGrid->GetCellTypesArray()->SetValue(edgevtk->getVtkId(), VTK_EMPTY_CELL);
335       myEdgePool->destroy(edgevtk);
336       return 0;
337     }
338   edge = edgevtk;
339   adjustmyCellsCapacity(ID);
340   myCells[ID] = edge;
341   myInfo.myNbEdges++;
342
343 //  if (edge && !registerElement(ID, edge))
344 //  {
345 //    RemoveElement(edge, false);
346 //    edge = NULL;
347 //  }
348   return edge;
349 }
350
351 ///////////////////////////////////////////////////////////////////////////////
352 /// Add a triangle defined by its nodes. An ID is automatically affected to the
353 /// Created face
354 ///////////////////////////////////////////////////////////////////////////////
355
356 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshNode * n1,
357                                   const SMDS_MeshNode * n2,
358                                   const SMDS_MeshNode * n3)
359 {
360   return SMDS_Mesh::AddFaceWithID(n1,n2,n3, myElementIDFactory->GetFreeID());
361 }
362
363 ///////////////////////////////////////////////////////////////////////////////
364 /// Add a triangle defined by its nodes IDs
365 ///////////////////////////////////////////////////////////////////////////////
366
367 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(int idnode1, int idnode2, int idnode3, int ID)
368 {
369   SMDS_MeshNode * node1 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode1);
370   SMDS_MeshNode * node2 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode2);
371   SMDS_MeshNode * node3 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode3);
372   if(!node1 || !node2 || !node3) return NULL;
373   return SMDS_Mesh::AddFaceWithID(node1, node2, node3, ID);
374 }
375
376 ///////////////////////////////////////////////////////////////////////////////
377 /// Add a triangle defined by its nodes
378 ///////////////////////////////////////////////////////////////////////////////
379
380 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshNode * n1,
381                                         const SMDS_MeshNode * n2,
382                                         const SMDS_MeshNode * n3,
383                                         int ID)
384 {
385   //MESSAGE("AddFaceWithID " << ID)
386   SMDS_MeshFace * face=createTriangle(n1, n2, n3, ID);
387
388 //  if (face && !registerElement(ID, face)) {
389 //    RemoveElement(face, false);
390 //    face = NULL;
391 //  }
392   return face;
393 }
394
395 ///////////////////////////////////////////////////////////////////////////////
396 /// Add a quadrangle defined by its nodes. An ID is automatically affected to the
397 /// created face
398 ///////////////////////////////////////////////////////////////////////////////
399
400 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshNode * n1,
401                                   const SMDS_MeshNode * n2,
402                                   const SMDS_MeshNode * n3,
403                                   const SMDS_MeshNode * n4)
404 {
405   return SMDS_Mesh::AddFaceWithID(n1,n2,n3, n4, myElementIDFactory->GetFreeID());
406 }
407
408 ///////////////////////////////////////////////////////////////////////////////
409 /// Add a quadrangle defined by its nodes IDs
410 ///////////////////////////////////////////////////////////////////////////////
411
412 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(int idnode1,
413                                         int idnode2,
414                                         int idnode3,
415                                         int idnode4,
416                                         int ID)
417 {
418   SMDS_MeshNode *node1, *node2, *node3, *node4;
419   node1 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode1);
420   node2 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode2);
421   node3 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode3);
422   node4 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode4);
423   if(!node1 || !node2 || !node3 || !node4) return NULL;
424   return SMDS_Mesh::AddFaceWithID(node1, node2, node3, node4, ID);
425 }
426
427 ///////////////////////////////////////////////////////////////////////////////
428 /// Add a quadrangle defined by its nodes
429 ///////////////////////////////////////////////////////////////////////////////
430
431 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshNode * n1,
432                                         const SMDS_MeshNode * n2,
433                                         const SMDS_MeshNode * n3,
434                                         const SMDS_MeshNode * n4,
435                                         int ID)
436 {
437   //MESSAGE("AddFaceWithID " << ID);
438   SMDS_MeshFace * face=createQuadrangle(n1, n2, n3, n4, ID);
439
440 //  if (face && !registerElement(ID, face)) {
441 //    RemoveElement(face, false);
442 //    face = NULL;
443 //  }
444   return face;
445 }
446
447 ///////////////////////////////////////////////////////////////////////////////
448 /// Add a triangle defined by its edges. An ID is automatically assigned to the
449 /// Created face
450 ///////////////////////////////////////////////////////////////////////////////
451
452 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshEdge * e1,
453                                   const SMDS_MeshEdge * e2,
454                                   const SMDS_MeshEdge * e3)
455 {
456   if (!hasConstructionEdges())
457     return NULL;
458      //MESSAGE("AddFaceWithID");
459  return AddFaceWithID(e1,e2,e3, myElementIDFactory->GetFreeID());
460 }
461
462 ///////////////////////////////////////////////////////////////////////////////
463 /// Add a triangle defined by its edges
464 ///////////////////////////////////////////////////////////////////////////////
465
466 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshEdge * e1,
467                                         const SMDS_MeshEdge * e2,
468                                         const SMDS_MeshEdge * e3,
469                                         int ID)
470 {
471   if (!hasConstructionEdges())
472     return NULL;
473   if ( !e1 || !e2 || !e3 ) return 0;
474
475   //if ( myFaces.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
476   MESSAGE("AddFaceWithID" << ID);
477
478   SMDS_MeshFace * face = new SMDS_FaceOfEdges(e1,e2,e3);
479   adjustmyCellsCapacity(ID);
480   myCells[ID] = face;
481   myInfo.myNbTriangles++;
482
483   if (!registerElement(ID, face)) {
484     registerElement(myElementIDFactory->GetFreeID(), face);
485     //RemoveElement(face, false);
486     //face = NULL;
487   }
488   return face;
489 }
490
491 ///////////////////////////////////////////////////////////////////////////////
492 /// Add a quadrangle defined by its edges. An ID is automatically assigned to the
493 /// Created face
494 ///////////////////////////////////////////////////////////////////////////////
495
496 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshEdge * e1,
497                                   const SMDS_MeshEdge * e2,
498                                   const SMDS_MeshEdge * e3,
499                                   const SMDS_MeshEdge * e4)
500 {
501   if (!hasConstructionEdges())
502     return NULL;
503      //MESSAGE("AddFaceWithID" );
504  return AddFaceWithID(e1,e2,e3,e4, myElementIDFactory->GetFreeID());
505 }
506
507 ///////////////////////////////////////////////////////////////////////////////
508 /// Add a quadrangle defined by its edges
509 ///////////////////////////////////////////////////////////////////////////////
510
511 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshEdge * e1,
512                                         const SMDS_MeshEdge * e2,
513                                         const SMDS_MeshEdge * e3,
514                                         const SMDS_MeshEdge * e4,
515                                         int ID)
516 {
517   if (!hasConstructionEdges())
518     return NULL;
519   MESSAGE("AddFaceWithID" << ID);
520   if ( !e1 || !e2 || !e3 || !e4 ) return 0;
521   //if ( myFaces.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
522   SMDS_MeshFace * face = new SMDS_FaceOfEdges(e1,e2,e3,e4);
523   adjustmyCellsCapacity(ID);
524   myCells[ID] = face;
525   myInfo.myNbQuadrangles++;
526
527   if (!registerElement(ID, face))
528   {
529     registerElement(myElementIDFactory->GetFreeID(), face);
530     //RemoveElement(face, false);
531     //face = NULL;
532   }
533   return face;
534 }
535
536 ///////////////////////////////////////////////////////////////////////////////
537 ///Create a new tetrahedron and add it to the mesh.
538 ///@return The created tetrahedron
539 ///////////////////////////////////////////////////////////////////////////////
540
541 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
542                                       const SMDS_MeshNode * n2,
543                                       const SMDS_MeshNode * n3,
544                                       const SMDS_MeshNode * n4)
545 {
546   int ID = myElementIDFactory->GetFreeID();
547     //MESSAGE("AddVolumeWithID " << ID);
548   SMDS_MeshVolume * v = SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, ID);
549   if(v==NULL) myElementIDFactory->ReleaseID(ID);
550   return v;
551 }
552
553 ///////////////////////////////////////////////////////////////////////////////
554 ///Create a new tetrahedron and add it to the mesh.
555 ///@param ID The ID of the new volume
556 ///@return The created tetrahedron or NULL if an element with this ID already exists
557 ///or if input nodes are not found.
558 ///////////////////////////////////////////////////////////////////////////////
559
560 SMDS_MeshVolume * SMDS_Mesh::AddVolumeWithID(int idnode1,
561                                              int idnode2,
562                                              int idnode3,
563                                              int idnode4,
564                                              int ID)
565 {
566     //MESSAGE("AddVolumeWithID" << ID);
567   SMDS_MeshNode *node1, *node2, *node3, *node4;
568   node1 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode1);
569   node2 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode2);
570   node3 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode3);
571   node4 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode4);
572   if(!node1 || !node2 || !node3 || !node4) return NULL;
573   return SMDS_Mesh::AddVolumeWithID(node1, node2, node3, node4, ID);
574 }
575
576 ///////////////////////////////////////////////////////////////////////////////
577 ///Create a new tetrahedron and add it to the mesh.
578 ///@param ID The ID of the new volume
579 ///@return The created tetrahedron
580 ///////////////////////////////////////////////////////////////////////////////
581
582 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
583                                             const SMDS_MeshNode * n2,
584                                             const SMDS_MeshNode * n3,
585                                             const SMDS_MeshNode * n4,
586                                             int ID)
587 {
588     //MESSAGE("AddVolumeWithID " << ID);
589   SMDS_MeshVolume* volume = 0;
590   if ( !n1 || !n2 || !n3 || !n4) return volume;
591   //if ( myVolumes.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
592   if(hasConstructionFaces()) {
593     SMDS_MeshFace * f1=FindFaceOrCreate(n1,n2,n3);
594     SMDS_MeshFace * f2=FindFaceOrCreate(n1,n2,n4);
595     SMDS_MeshFace * f3=FindFaceOrCreate(n1,n3,n4);
596     SMDS_MeshFace * f4=FindFaceOrCreate(n2,n3,n4);
597     volume=new SMDS_VolumeOfFaces(f1,f2,f3,f4);
598     adjustmyCellsCapacity(ID);
599     myCells[ID] = volume;
600     myInfo.myNbTetras++;
601   }
602   else if(hasConstructionEdges()) {
603     MESSAGE("Error : Not implemented");
604     return NULL;
605   }
606   else {
607     // --- retrieve nodes ID
608     vector<vtkIdType> nodeIds;
609     nodeIds.clear();
610     nodeIds.push_back(n1->getVtkId());
611     nodeIds.push_back(n3->getVtkId()); // order SMDS-->VTK
612     nodeIds.push_back(n2->getVtkId());
613     nodeIds.push_back(n4->getVtkId());
614
615     SMDS_VtkVolume *volvtk = myVolumePool->getNew();
616     volvtk->init(nodeIds, this);
617     if (!this->registerElement(ID,volvtk))
618       {
619         this->myGrid->GetCellTypesArray()->SetValue(volvtk->getVtkId(), VTK_EMPTY_CELL);
620         myVolumePool->destroy(volvtk);
621         return 0;
622       }
623     volume = volvtk;
624     adjustmyCellsCapacity(ID);
625     myCells[ID] = volume;
626     myInfo.myNbTetras++;
627   }
628
629 //  if (!registerElement(ID, volume)) {
630 //    RemoveElement(volume, false);
631 //    volume = NULL;
632 //  }
633   return volume;
634 }
635
636 ///////////////////////////////////////////////////////////////////////////////
637 ///Create a new pyramid and add it to the mesh.
638 ///Nodes 1,2,3 and 4 define the base of the pyramid
639 ///@return The created pyramid
640 ///////////////////////////////////////////////////////////////////////////////
641
642 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
643                                       const SMDS_MeshNode * n2,
644                                       const SMDS_MeshNode * n3,
645                                       const SMDS_MeshNode * n4,
646                                       const SMDS_MeshNode * n5)
647 {
648   int ID = myElementIDFactory->GetFreeID();
649     //MESSAGE("AddVolumeWithID " << ID);
650   SMDS_MeshVolume * v = SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, ID);
651   if(v==NULL) myElementIDFactory->ReleaseID(ID);
652   return v;
653 }
654
655 ///////////////////////////////////////////////////////////////////////////////
656 ///Create a new pyramid and add it to the mesh.
657 ///Nodes 1,2,3 and 4 define the base of the pyramid
658 ///@param ID The ID of the new volume
659 ///@return The created pyramid or NULL if an element with this ID already exists
660 ///or if input nodes are not found.
661 ///////////////////////////////////////////////////////////////////////////////
662
663 SMDS_MeshVolume * SMDS_Mesh::AddVolumeWithID(int idnode1,
664                                              int idnode2,
665                                              int idnode3,
666                                              int idnode4,
667                                              int idnode5,
668                                              int ID)
669 {
670     //MESSAGE("AddVolumeWithID " << ID);
671   SMDS_MeshNode *node1, *node2, *node3, *node4, *node5;
672   node1 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode1);
673   node2 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode2);
674   node3 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode3);
675   node4 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode4);
676   node5 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode5);
677   if(!node1 || !node2 || !node3 || !node4 || !node5) return NULL;
678   return SMDS_Mesh::AddVolumeWithID(node1, node2, node3, node4, node5, ID);
679 }
680
681 ///////////////////////////////////////////////////////////////////////////////
682 ///Create a new pyramid and add it to the mesh.
683 ///Nodes 1,2,3 and 4 define the base of the pyramid
684 ///@param ID The ID of the new volume
685 ///@return The created pyramid
686 ///////////////////////////////////////////////////////////////////////////////
687
688 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
689                                             const SMDS_MeshNode * n2,
690                                             const SMDS_MeshNode * n3,
691                                             const SMDS_MeshNode * n4,
692                                             const SMDS_MeshNode * n5,
693                                             int ID)
694 {
695     //MESSAGE("AddVolumeWithID " << ID);
696   SMDS_MeshVolume* volume = 0;
697   if ( !n1 || !n2 || !n3 || !n4 || !n5) return volume;
698   //if ( myVolumes.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
699   if(hasConstructionFaces()) {
700     SMDS_MeshFace * f1=FindFaceOrCreate(n1,n2,n3,n4);
701     SMDS_MeshFace * f2=FindFaceOrCreate(n1,n2,n5);
702     SMDS_MeshFace * f3=FindFaceOrCreate(n2,n3,n5);
703     SMDS_MeshFace * f4=FindFaceOrCreate(n3,n4,n5);
704     volume=new SMDS_VolumeOfFaces(f1,f2,f3,f4);
705     adjustmyCellsCapacity(ID);
706     myCells[ID] = volume;
707     myInfo.myNbPyramids++;
708   }
709   else if(hasConstructionEdges()) {
710     MESSAGE("Error : Not implemented");
711     return NULL;
712   }
713   else {
714     // --- retrieve nodes ID
715     vector<vtkIdType> nodeIds;
716     nodeIds.clear();
717     nodeIds.push_back(n1->getVtkId());
718     nodeIds.push_back(n4->getVtkId());
719     nodeIds.push_back(n3->getVtkId());
720     nodeIds.push_back(n2->getVtkId());
721     nodeIds.push_back(n5->getVtkId());
722
723     SMDS_VtkVolume *volvtk = myVolumePool->getNew();
724     volvtk->init(nodeIds, this);
725     if (!this->registerElement(ID,volvtk))
726       {
727         this->myGrid->GetCellTypesArray()->SetValue(volvtk->getVtkId(), VTK_EMPTY_CELL);
728         myVolumePool->destroy(volvtk);
729         return 0;
730       }
731     volume = volvtk;
732     adjustmyCellsCapacity(ID);
733     myCells[ID] = volume;
734     myInfo.myNbPyramids++;
735   }
736
737 //  if (!registerElement(ID, volume)) {
738 //    RemoveElement(volume, false);
739 //    volume = NULL;
740 //  }
741   return volume;
742 }
743
744 ///////////////////////////////////////////////////////////////////////////////
745 ///Create a new prism and add it to the mesh.
746 ///Nodes 1,2,3 is a triangle and 1,2,5,4 a quadrangle.
747 ///@return The created prism
748 ///////////////////////////////////////////////////////////////////////////////
749
750 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
751                                       const SMDS_MeshNode * n2,
752                                       const SMDS_MeshNode * n3,
753                                       const SMDS_MeshNode * n4,
754                                       const SMDS_MeshNode * n5,
755                                       const SMDS_MeshNode * n6)
756 {
757   int ID = myElementIDFactory->GetFreeID();
758     //MESSAGE("AddVolumeWithID " << ID);
759   SMDS_MeshVolume * v = SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n6, ID);
760   if(v==NULL) myElementIDFactory->ReleaseID(ID);
761   return v;
762 }
763
764 ///////////////////////////////////////////////////////////////////////////////
765 ///Create a new prism and add it to the mesh.
766 ///Nodes 1,2,3 is a triangle and 1,2,5,4 a quadrangle.
767 ///@param ID The ID of the new volume
768 ///@return The created prism or NULL if an element with this ID already exists
769 ///or if input nodes are not found.
770 ///////////////////////////////////////////////////////////////////////////////
771
772 SMDS_MeshVolume * SMDS_Mesh::AddVolumeWithID(int idnode1,
773                                              int idnode2,
774                                              int idnode3,
775                                              int idnode4,
776                                              int idnode5,
777                                              int idnode6,
778                                              int ID)
779 {
780     //MESSAGE("AddVolumeWithID " << ID);
781   SMDS_MeshNode *node1, *node2, *node3, *node4, *node5, *node6;
782   node1 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode1);
783   node2 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode2);
784   node3 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode3);
785   node4 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode4);
786   node5 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode5);
787   node6 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode6);
788   if(!node1 || !node2 || !node3 || !node4 || !node5 || !node6) return NULL;
789   return SMDS_Mesh::AddVolumeWithID(node1, node2, node3, node4, node5, node6, ID);
790 }
791
792 ///////////////////////////////////////////////////////////////////////////////
793 ///Create a new prism and add it to the mesh.
794 ///Nodes 1,2,3 is a triangle and 1,2,5,4 a quadrangle.
795 ///@param ID The ID of the new volume
796 ///@return The created prism
797 ///////////////////////////////////////////////////////////////////////////////
798
799 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
800                                             const SMDS_MeshNode * n2,
801                                             const SMDS_MeshNode * n3,
802                                             const SMDS_MeshNode * n4,
803                                             const SMDS_MeshNode * n5,
804                                             const SMDS_MeshNode * n6,
805                                             int ID)
806 {
807     //MESSAGE("AddVolumeWithID " << ID);
808   SMDS_MeshVolume* volume = 0;
809   if ( !n1 || !n2 || !n3 || !n4 || !n5 || !n6) return volume;
810   //if ( myVolumes.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
811   if(hasConstructionFaces()) {
812     SMDS_MeshFace * f1=FindFaceOrCreate(n1,n2,n3);
813     SMDS_MeshFace * f2=FindFaceOrCreate(n4,n5,n6);
814     SMDS_MeshFace * f3=FindFaceOrCreate(n1,n4,n5,n2);
815     SMDS_MeshFace * f4=FindFaceOrCreate(n2,n5,n6,n3);
816     SMDS_MeshFace * f5=FindFaceOrCreate(n3,n6,n4,n1);
817     volume=new SMDS_VolumeOfFaces(f1,f2,f3,f4,f5);
818     adjustmyCellsCapacity(ID);
819     myCells[ID] = volume;
820     myInfo.myNbPrisms++;
821   }
822   else if(hasConstructionEdges()) {
823     MESSAGE("Error : Not implemented");
824     return NULL;
825   }
826   else {
827     // --- retrieve nodes ID
828     vector<vtkIdType> nodeIds;
829     nodeIds.clear();
830     nodeIds.push_back(n1->getVtkId());
831     nodeIds.push_back(n2->getVtkId());
832     nodeIds.push_back(n3->getVtkId());
833     nodeIds.push_back(n4->getVtkId());
834     nodeIds.push_back(n5->getVtkId());
835     nodeIds.push_back(n6->getVtkId());
836
837     SMDS_VtkVolume *volvtk = myVolumePool->getNew();
838     volvtk->init(nodeIds, this);
839     if (!this->registerElement(ID,volvtk))
840       {
841         this->myGrid->GetCellTypesArray()->SetValue(volvtk->getVtkId(), VTK_EMPTY_CELL);
842         myVolumePool->destroy(volvtk);
843         return 0;
844       }
845     volume = volvtk;
846     adjustmyCellsCapacity(ID);
847     myCells[ID] = volume;
848     myInfo.myNbPrisms++;
849   }
850
851 //  if (!registerElement(ID, volume)) {
852 //    RemoveElement(volume, false);
853 //    volume = NULL;
854 //  }
855   return volume;
856 }
857
858 ///////////////////////////////////////////////////////////////////////////////
859 ///Create a new hexahedron and add it to the mesh.
860 ///Nodes 1,2,3,4 and 5,6,7,8 are quadrangle and 5,1 and 7,3 are an edges.
861 ///@return The created hexahedron
862 ///////////////////////////////////////////////////////////////////////////////
863
864 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
865                                       const SMDS_MeshNode * n2,
866                                       const SMDS_MeshNode * n3,
867                                       const SMDS_MeshNode * n4,
868                                       const SMDS_MeshNode * n5,
869                                       const SMDS_MeshNode * n6,
870                                       const SMDS_MeshNode * n7,
871                                       const SMDS_MeshNode * n8)
872 {
873   int ID = myElementIDFactory->GetFreeID();
874      //MESSAGE("AddVolumeWithID " << ID);
875  SMDS_MeshVolume * v = SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n6, n7, n8, ID);
876   if(v==NULL) myElementIDFactory->ReleaseID(ID);
877   return v;
878 }
879
880 ///////////////////////////////////////////////////////////////////////////////
881 ///Create a new hexahedron and add it to the mesh.
882 ///Nodes 1,2,3,4 and 5,6,7,8 are quadrangle and 5,1 and 7,3 are an edges.
883 ///@param ID The ID of the new volume
884 ///@return The created hexahedron or NULL if an element with this ID already
885 ///exists or if input nodes are not found.
886 ///////////////////////////////////////////////////////////////////////////////
887
888 SMDS_MeshVolume * SMDS_Mesh::AddVolumeWithID(int idnode1,
889                                              int idnode2,
890                                              int idnode3,
891                                              int idnode4,
892                                              int idnode5,
893                                              int idnode6,
894                                              int idnode7,
895                                              int idnode8,
896                                              int ID)
897 {
898     //MESSAGE("AddVolumeWithID " << ID);
899   SMDS_MeshNode *node1, *node2, *node3, *node4, *node5, *node6, *node7, *node8;
900   node1 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode1);
901   node2 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode2);
902   node3 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode3);
903   node4 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode4);
904   node5 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode5);
905   node6 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode6);
906   node7 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode7);
907   node8 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode8);
908   if(!node1 || !node2 || !node3 || !node4 || !node5 || !node6 || !node7 || !node8)
909     return NULL;
910   return SMDS_Mesh::AddVolumeWithID(node1, node2, node3, node4, node5, node6,
911                                     node7, node8, ID);
912 }
913
914 ///////////////////////////////////////////////////////////////////////////////
915 ///Create a new hexahedron and add it to the mesh.
916 ///Nodes 1,2,3,4 and 5,6,7,8 are quadrangle and 5,1 and 7,3 are an edges.
917 ///@param ID The ID of the new volume
918 ///@return The created prism or NULL if an element with this ID already exists
919 ///or if input nodes are not found.
920 ///////////////////////////////////////////////////////////////////////////////
921
922 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
923                                             const SMDS_MeshNode * n2,
924                                             const SMDS_MeshNode * n3,
925                                             const SMDS_MeshNode * n4,
926                                             const SMDS_MeshNode * n5,
927                                             const SMDS_MeshNode * n6,
928                                             const SMDS_MeshNode * n7,
929                                             const SMDS_MeshNode * n8,
930                                             int ID)
931 {
932     //MESSAGE("AddVolumeWithID " << ID);
933   SMDS_MeshVolume* volume = 0;
934   if ( !n1 || !n2 || !n3 || !n4 || !n5 || !n6 || !n7 || !n8) return volume;
935   //if ( myVolumes.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
936   if(hasConstructionFaces()) {
937     SMDS_MeshFace * f1=FindFaceOrCreate(n1,n2,n3,n4);
938     SMDS_MeshFace * f2=FindFaceOrCreate(n5,n6,n7,n8);
939     SMDS_MeshFace * f3=FindFaceOrCreate(n1,n4,n8,n5);
940     SMDS_MeshFace * f4=FindFaceOrCreate(n1,n2,n6,n5);
941     SMDS_MeshFace * f5=FindFaceOrCreate(n2,n3,n7,n6);
942     SMDS_MeshFace * f6=FindFaceOrCreate(n3,n4,n8,n7);
943     volume=new SMDS_VolumeOfFaces(f1,f2,f3,f4,f5,f6);
944     adjustmyCellsCapacity(ID);
945     myCells[ID] = volume;
946     myInfo.myNbHexas++;
947   }
948   else if(hasConstructionEdges()) {
949     MESSAGE("Error : Not implemented");
950     return NULL;
951   }
952   else {
953     // --- retrieve nodes ID
954     vector<vtkIdType> nodeIds;
955     nodeIds.clear();
956     nodeIds.push_back(n1->getVtkId());
957     nodeIds.push_back(n4->getVtkId());
958     nodeIds.push_back(n3->getVtkId());
959     nodeIds.push_back(n2->getVtkId());
960     nodeIds.push_back(n5->getVtkId());
961     nodeIds.push_back(n8->getVtkId());
962     nodeIds.push_back(n7->getVtkId());
963     nodeIds.push_back(n6->getVtkId());
964
965     SMDS_VtkVolume *volvtk = myVolumePool->getNew();
966     volvtk->init(nodeIds, this);
967     if (!this->registerElement(ID,volvtk))
968       {
969         this->myGrid->GetCellTypesArray()->SetValue(volvtk->getVtkId(), VTK_EMPTY_CELL);
970         myVolumePool->destroy(volvtk);
971         return 0;
972       }
973     volume = volvtk;
974     adjustmyCellsCapacity(ID);
975     myCells[ID] = volume;
976     myInfo.myNbHexas++;
977   }
978  
979 //  if (!registerElement(ID, volume)) {
980 //    RemoveElement(volume, false);
981 //    volume = NULL;
982 //  }
983   return volume;
984 }
985
986 ///////////////////////////////////////////////////////////////////////////////
987 ///Create a new tetrahedron defined by its faces and add it to the mesh.
988 ///@return The created tetrahedron
989 ///////////////////////////////////////////////////////////////////////////////
990
991 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshFace * f1,
992                                       const SMDS_MeshFace * f2,
993                                       const SMDS_MeshFace * f3,
994                                       const SMDS_MeshFace * f4)
995 {
996     //MESSAGE("AddVolumeWithID");
997   if (!hasConstructionFaces())
998     return NULL;
999   return AddVolumeWithID(f1,f2,f3,f4, myElementIDFactory->GetFreeID());
1000 }
1001
1002 ///////////////////////////////////////////////////////////////////////////////
1003 ///Create a new tetrahedron defined by its faces and add it to the mesh.
1004 ///@param ID The ID of the new volume
1005 ///@return The created tetrahedron
1006 ///////////////////////////////////////////////////////////////////////////////
1007
1008 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshFace * f1,
1009                                             const SMDS_MeshFace * f2,
1010                                             const SMDS_MeshFace * f3,
1011                                             const SMDS_MeshFace * f4,
1012                                             int ID)
1013 {
1014   MESSAGE("AddVolumeWithID" << ID);
1015   if (!hasConstructionFaces())
1016     return NULL;
1017   if ( !f1 || !f2 || !f3 || !f4) return 0;
1018   //if ( myVolumes.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
1019   SMDS_MeshVolume * volume = new SMDS_VolumeOfFaces(f1,f2,f3,f4);
1020   adjustmyCellsCapacity(ID);
1021   myCells[ID] = volume;
1022   myInfo.myNbTetras++;
1023
1024   if (!registerElement(ID, volume)) {
1025     registerElement(myElementIDFactory->GetFreeID(), volume);
1026     //RemoveElement(volume, false);
1027     //volume = NULL;
1028   }
1029   return volume;
1030 }
1031
1032 ///////////////////////////////////////////////////////////////////////////////
1033 ///Create a new pyramid defined by its faces and add it to the mesh.
1034 ///@return The created pyramid
1035 ///////////////////////////////////////////////////////////////////////////////
1036
1037 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshFace * f1,
1038                                       const SMDS_MeshFace * f2,
1039                                       const SMDS_MeshFace * f3,
1040                                       const SMDS_MeshFace * f4,
1041                                       const SMDS_MeshFace * f5)
1042 {
1043      //MESSAGE("AddVolumeWithID");
1044  if (!hasConstructionFaces())
1045     return NULL;
1046   return AddVolumeWithID(f1,f2,f3,f4,f5, myElementIDFactory->GetFreeID());
1047 }
1048
1049 ///////////////////////////////////////////////////////////////////////////////
1050 ///Create a new pyramid defined by its faces and add it to the mesh.
1051 ///@param ID The ID of the new volume
1052 ///@return The created pyramid
1053 ///////////////////////////////////////////////////////////////////////////////
1054
1055 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshFace * f1,
1056                                             const SMDS_MeshFace * f2,
1057                                             const SMDS_MeshFace * f3,
1058                                             const SMDS_MeshFace * f4,
1059                                             const SMDS_MeshFace * f5,
1060                                             int ID)
1061 {
1062   MESSAGE("AddVolumeWithID" << ID);
1063   if (!hasConstructionFaces())
1064     return NULL;
1065   if ( !f1 || !f2 || !f3 || !f4 || !f5) return 0;
1066   //if ( myVolumes.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
1067   SMDS_MeshVolume * volume = new SMDS_VolumeOfFaces(f1,f2,f3,f4,f5);
1068   adjustmyCellsCapacity(ID);
1069   myCells[ID] = volume;
1070   myInfo.myNbPyramids++;
1071
1072   if (!registerElement(ID, volume)) {
1073     registerElement(myElementIDFactory->GetFreeID(), volume);
1074     //RemoveElement(volume, false);
1075     //volume = NULL;
1076   }
1077   return volume;
1078 }
1079
1080 ///////////////////////////////////////////////////////////////////////////////
1081 ///Create a new prism defined by its faces and add it to the mesh.
1082 ///@return The created prism
1083 ///////////////////////////////////////////////////////////////////////////////
1084
1085 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshFace * f1,
1086                                       const SMDS_MeshFace * f2,
1087                                       const SMDS_MeshFace * f3,
1088                                       const SMDS_MeshFace * f4,
1089                                       const SMDS_MeshFace * f5,
1090                                       const SMDS_MeshFace * f6)
1091 {
1092      //MESSAGE("AddVolumeWithID" );
1093  if (!hasConstructionFaces())
1094     return NULL;
1095   return AddVolumeWithID(f1,f2,f3,f4,f5,f6, myElementIDFactory->GetFreeID());
1096 }
1097
1098 ///////////////////////////////////////////////////////////////////////////////
1099 ///Create a new prism defined by its faces and add it to the mesh.
1100 ///@param ID The ID of the new volume
1101 ///@return The created prism
1102 ///////////////////////////////////////////////////////////////////////////////
1103
1104 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshFace * f1,
1105                                             const SMDS_MeshFace * f2,
1106                                             const SMDS_MeshFace * f3,
1107                                             const SMDS_MeshFace * f4,
1108                                             const SMDS_MeshFace * f5,
1109                                             const SMDS_MeshFace * f6,
1110                                             int ID)
1111 {
1112   MESSAGE("AddVolumeWithID" << ID);
1113   if (!hasConstructionFaces())
1114     return NULL;
1115   if ( !f1 || !f2 || !f3 || !f4 || !f5 || !f6) return 0;
1116   //if ( myVolumes.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
1117   SMDS_MeshVolume * volume = new SMDS_VolumeOfFaces(f1,f2,f3,f4,f5,f6);
1118   adjustmyCellsCapacity(ID);
1119   myCells[ID] = volume;
1120   myInfo.myNbPrisms++;
1121
1122   if (!registerElement(ID, volume)) {
1123     registerElement(myElementIDFactory->GetFreeID(), volume);
1124     //RemoveElement(volume, false);
1125     //volume = NULL;
1126   }
1127   return volume;
1128 }
1129
1130 ///////////////////////////////////////////////////////////////////////////////
1131 /// Add a polygon defined by its nodes IDs
1132 ///////////////////////////////////////////////////////////////////////////////
1133
1134 SMDS_MeshFace* SMDS_Mesh::AddPolygonalFaceWithID (vector<int> nodes_ids,
1135                                                   const int        ID)
1136 {
1137   int nbNodes = nodes_ids.size();
1138   vector<const SMDS_MeshNode*> nodes (nbNodes);
1139   for (int i = 0; i < nbNodes; i++) {
1140     nodes[i] = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(nodes_ids[i]);
1141     if (!nodes[i]) return NULL;
1142   }
1143   return SMDS_Mesh::AddPolygonalFaceWithID(nodes, ID);
1144 }
1145
1146 ///////////////////////////////////////////////////////////////////////////////
1147 /// Add a polygon defined by its nodes
1148 ///////////////////////////////////////////////////////////////////////////////
1149
1150 SMDS_MeshFace* SMDS_Mesh::AddPolygonalFaceWithID
1151                           (vector<const SMDS_MeshNode*> nodes,
1152                            const int                         ID)
1153 {
1154   SMDS_MeshFace * face;
1155
1156   //if ( myFaces.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
1157   if (hasConstructionEdges())
1158     {
1159       MESSAGE("Error : Not implemented");
1160       return NULL;
1161     }
1162   else
1163     {
1164 //#ifdef VTK_HAVE_POLYHEDRON
1165     MESSAGE("AddPolygonalFaceWithID vtk " << ID);
1166     vector<vtkIdType> nodeIds;
1167     nodeIds.clear();
1168     vector<const SMDS_MeshNode*>::iterator it = nodes.begin();
1169     for ( ; it != nodes.end(); ++it)
1170       nodeIds.push_back((*it)->getVtkId());
1171
1172     SMDS_VtkFace *facevtk = myFacePool->getNew();
1173     facevtk->initPoly(nodeIds, this);
1174     if (!this->registerElement(ID,facevtk))
1175       {
1176         this->myGrid->GetCellTypesArray()->SetValue(facevtk->getVtkId(), VTK_EMPTY_CELL);
1177         myFacePool->destroy(facevtk);
1178         return 0;
1179       }
1180     face = facevtk;
1181 //#else
1182 //    MESSAGE("AddPolygonalFaceWithID smds " << ID);
1183 //     for ( int i = 0; i < nodes.size(); ++i )
1184 //      if ( !nodes[ i ] ) return 0;
1185 //      face = new SMDS_PolygonalFaceOfNodes(nodes);
1186 //#endif
1187       adjustmyCellsCapacity(ID);
1188       myCells[ID] = face;
1189       myInfo.myNbPolygons++;
1190     }
1191
1192 //#ifndef VTK_HAVE_POLYHEDRON
1193 //  if (!registerElement(ID, face))
1194 //    {
1195 //      registerElement(myElementIDFactory->GetFreeID(), face);
1196 //      //RemoveElement(face, false);
1197 //      //face = NULL;
1198 //    }
1199 //#endif
1200  return face;
1201 }
1202
1203 ///////////////////////////////////////////////////////////////////////////////
1204 /// Add a polygon defined by its nodes.
1205 /// An ID is automatically affected to the created face.
1206 ///////////////////////////////////////////////////////////////////////////////
1207
1208 SMDS_MeshFace* SMDS_Mesh::AddPolygonalFace (vector<const SMDS_MeshNode*> nodes)
1209 {
1210   return SMDS_Mesh::AddPolygonalFaceWithID(nodes, myElementIDFactory->GetFreeID());
1211 }
1212
1213 ///////////////////////////////////////////////////////////////////////////////
1214 /// Create a new polyhedral volume and add it to the mesh.
1215 /// @param ID The ID of the new volume
1216 /// @return The created volume or NULL if an element with this ID already exists
1217 /// or if input nodes are not found.
1218 ///////////////////////////////////////////////////////////////////////////////
1219
1220 SMDS_MeshVolume * SMDS_Mesh::AddPolyhedralVolumeWithID
1221                              (vector<int> nodes_ids,
1222                               vector<int> quantities,
1223                               const int        ID)
1224 {
1225   int nbNodes = nodes_ids.size();
1226   vector<const SMDS_MeshNode*> nodes (nbNodes);
1227   for (int i = 0; i < nbNodes; i++) {
1228     nodes[i] = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(nodes_ids[i]);
1229     if (!nodes[i]) return NULL;
1230   }
1231   return SMDS_Mesh::AddPolyhedralVolumeWithID(nodes, quantities, ID);
1232 }
1233
1234 ///////////////////////////////////////////////////////////////////////////////
1235 /// Create a new polyhedral volume and add it to the mesh.
1236 /// @param ID The ID of the new volume
1237 /// @return The created  volume
1238 ///////////////////////////////////////////////////////////////////////////////
1239
1240 SMDS_MeshVolume* SMDS_Mesh::AddPolyhedralVolumeWithID
1241                             (vector<const SMDS_MeshNode*> nodes,
1242                              vector<int>                  quantities,
1243                              const int                    ID)
1244 {
1245   SMDS_MeshVolume* volume;
1246   //if ( myVolumes.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
1247   if (hasConstructionFaces())
1248     {
1249       MESSAGE("Error : Not implemented");
1250       return NULL;
1251     }
1252   else if (hasConstructionEdges())
1253     {
1254       MESSAGE("Error : Not implemented");
1255       return NULL;
1256     }
1257   else
1258     {
1259 //#ifdef VTK_HAVE_POLYHEDRON
1260       MESSAGE("AddPolyhedralVolumeWithID vtk " << ID);
1261       vector<vtkIdType> nodeIds;
1262       nodeIds.clear();
1263       vector<const SMDS_MeshNode*>::iterator it = nodes.begin();
1264       for (; it != nodes.end(); ++it)
1265         nodeIds.push_back((*it)->getVtkId());
1266
1267       SMDS_VtkVolume *volvtk = myVolumePool->getNew();
1268       volvtk->initPoly(nodeIds, quantities, this);
1269       if (!this->registerElement(ID, volvtk))
1270         {
1271           this->myGrid->GetCellTypesArray()->SetValue(volvtk->getVtkId(), VTK_EMPTY_CELL);
1272           myVolumePool->destroy(volvtk);
1273           return 0;
1274         }
1275       volume = volvtk;
1276 //#else
1277 //      MESSAGE("AddPolyhedralVolumeWithID smds " << ID);
1278 //      for ( int i = 0; i < nodes.size(); ++i )
1279 //      if ( !nodes[ i ] ) return 0;
1280 //      volume = new SMDS_PolyhedralVolumeOfNodes(nodes, quantities);
1281 //#endif
1282       adjustmyCellsCapacity(ID);
1283       myCells[ID] = volume;
1284       myInfo.myNbPolyhedrons++;
1285     }
1286
1287 //#ifndef VTK_HAVE_POLYHEDRON
1288 //  if (!registerElement(ID, volume))
1289 //    {
1290 //      registerElement(myElementIDFactory->GetFreeID(), volume);
1291 //      //RemoveElement(volume, false);
1292 //      //volume = NULL;
1293 //    }
1294 //#endif
1295   return volume;
1296 }
1297
1298 ///////////////////////////////////////////////////////////////////////////////
1299 /// Create a new polyhedral volume and add it to the mesh.
1300 /// @return The created  volume
1301 ///////////////////////////////////////////////////////////////////////////////
1302
1303 SMDS_MeshVolume* SMDS_Mesh::AddPolyhedralVolume
1304                             (vector<const SMDS_MeshNode*> nodes,
1305                              vector<int>                  quantities)
1306 {
1307   int ID = myElementIDFactory->GetFreeID();
1308   SMDS_MeshVolume * v = SMDS_Mesh::AddPolyhedralVolumeWithID(nodes, quantities, ID);
1309   if (v == NULL) myElementIDFactory->ReleaseID(ID);
1310   return v;
1311 }
1312
1313 SMDS_MeshVolume* SMDS_Mesh::AddVolumeFromVtkIds(const std::vector<int>& vtkNodeIds)
1314 {
1315   int ID = myElementIDFactory->GetFreeID();
1316   SMDS_MeshVolume * v = SMDS_Mesh::AddVolumeFromVtkIdsWithID(vtkNodeIds, ID);
1317   if (v == NULL) myElementIDFactory->ReleaseID(ID);
1318   return v;
1319 }
1320
1321 SMDS_MeshVolume* SMDS_Mesh::AddVolumeFromVtkIdsWithID(const std::vector<int>& vtkNodeIds, const int ID)
1322 {
1323   SMDS_VtkVolume *volvtk = myVolumePool->getNew();
1324   volvtk->init(vtkNodeIds, this);
1325   if (!this->registerElement(ID,volvtk))
1326     {
1327       this->myGrid->GetCellTypesArray()->SetValue(volvtk->getVtkId(), VTK_EMPTY_CELL);
1328       myVolumePool->destroy(volvtk);
1329       return 0;
1330     }
1331   adjustmyCellsCapacity(ID);
1332   myCells[ID] = volvtk;
1333   vtkIdType aVtkType = volvtk->GetVtkType();
1334   switch (aVtkType)
1335   {
1336     case VTK_TETRA:
1337       myInfo.myNbTetras++;
1338       break;
1339     case VTK_PYRAMID:
1340       myInfo.myNbPyramids++;
1341       break;
1342     case VTK_WEDGE:
1343       myInfo.myNbPrisms++;
1344       break;
1345     case VTK_HEXAHEDRON:
1346       myInfo.myNbHexas++;
1347       break;
1348     case VTK_QUADRATIC_TETRA:
1349       myInfo.myNbQuadTetras++;
1350       break;
1351     case VTK_QUADRATIC_PYRAMID:
1352       myInfo.myNbQuadPyramids++;
1353       break;
1354     case VTK_QUADRATIC_WEDGE:
1355       myInfo.myNbQuadPrisms++;
1356       break;
1357     case VTK_QUADRATIC_HEXAHEDRON:
1358       myInfo.myNbQuadHexas++;
1359       break;
1360 //#ifdef VTK_HAVE_POLYHEDRON
1361     case VTK_POLYHEDRON:
1362       myInfo.myNbPolyhedrons++;
1363       break;
1364 //#endif
1365     default:
1366       myInfo.myNbPolyhedrons++;
1367       break;
1368   }
1369   return volvtk;
1370 }
1371
1372 ///////////////////////////////////////////////////////////////////////////////
1373 /// Registers element with the given ID, maintains inverse connections
1374 ///////////////////////////////////////////////////////////////////////////////
1375 bool SMDS_Mesh::registerElement(int ID, SMDS_MeshElement* element)
1376 {
1377   //MESSAGE("registerElement " << ID);
1378   if ((ID >=0) && (ID < myCells.size()) && myCells[ID]) // --- already bound
1379   {
1380     MESSAGE(" ------------------ already bound "<< ID << " " << myCells[ID]->getVtkId());
1381     return false;
1382   }
1383
1384   element->myID = ID;
1385   element->myMeshId = myMeshId;
1386
1387   SMDS_MeshCell *cell = dynamic_cast<SMDS_MeshCell*>(element);
1388   MYASSERT(cell);
1389   int vtkId = cell->getVtkId();  
1390   if (vtkId == -1)
1391     vtkId = myElementIDFactory->SetInVtkGrid(element);
1392
1393   if (vtkId >= myCellIdVtkToSmds.size()) // --- resize local vector
1394   {
1395     MESSAGE(" --------------------- resize myCellIdVtkToSmds " << vtkId << " --> " << vtkId + SMDS_Mesh::chunkSize);
1396     myCellIdVtkToSmds.resize(vtkId + SMDS_Mesh::chunkSize, -1);
1397   }
1398   myCellIdVtkToSmds[vtkId] = ID;
1399
1400   myElementIDFactory->updateMinMax(ID);
1401   return true;
1402 }
1403
1404 ///////////////////////////////////////////////////////////////////////////////
1405 /// Return the node whose SMDS ID is 'ID'.
1406 ///////////////////////////////////////////////////////////////////////////////
1407 const SMDS_MeshNode * SMDS_Mesh::FindNode(int ID) const
1408 {
1409   if (ID < 1 || ID >= myNodes.size())
1410   {
1411     MESSAGE("------------------------------------------------------------------------- ");
1412     MESSAGE("----------------------------------- bad ID " << ID << " " << myNodes.size());
1413     MESSAGE("------------------------------------------------------------------------- ");
1414     return 0;
1415   }
1416   return (const SMDS_MeshNode *)myNodes[ID];
1417 }
1418
1419 ///////////////////////////////////////////////////////////////////////////////
1420 /// Return the node whose VTK ID is 'vtkId'.
1421 ///////////////////////////////////////////////////////////////////////////////
1422 const SMDS_MeshNode * SMDS_Mesh::FindNodeVtk(int vtkId) const
1423 {
1424   // TODO if needed use mesh->nodeIdFromVtkToSmds
1425   if (vtkId < 0 || vtkId >= (myNodes.size() -1))
1426   {
1427     MESSAGE("------------------------------------------------------------------------- ");
1428     MESSAGE("---------------------------- bad VTK ID " << vtkId << " " << myNodes.size());
1429     MESSAGE("------------------------------------------------------------------------- ");
1430     return 0;
1431   }
1432   return (const SMDS_MeshNode *)myNodes[vtkId+1];
1433 }
1434
1435 ///////////////////////////////////////////////////////////////////////////////
1436 ///Create a triangle and add it to the current mesh. This method do not bind an
1437 ///ID to the create triangle.
1438 ///////////////////////////////////////////////////////////////////////////////
1439 SMDS_MeshFace * SMDS_Mesh::createTriangle(const SMDS_MeshNode * node1,
1440                                           const SMDS_MeshNode * node2,
1441                                           const SMDS_MeshNode * node3,
1442                                           int ID)
1443 {
1444   if ( !node1 || !node2 || !node3) return 0;
1445 //  if ( myFaces.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
1446   if(hasConstructionEdges())
1447   {
1448     SMDS_MeshEdge *edge1, *edge2, *edge3;
1449     edge1=FindEdgeOrCreate(node1,node2);
1450     edge2=FindEdgeOrCreate(node2,node3);
1451     edge3=FindEdgeOrCreate(node3,node1);
1452
1453     //int ID = myElementIDFactory->GetFreeID(); // -PR- voir si on range cet element
1454     SMDS_MeshFace * face = new SMDS_FaceOfEdges(edge1,edge2,edge3);
1455     adjustmyCellsCapacity(ID);
1456     myCells[ID] = face;
1457     myInfo.myNbTriangles++;
1458     return face;
1459   }
1460   else
1461   {
1462     // --- retrieve nodes ID
1463     vector<vtkIdType> nodeIds;
1464     nodeIds.clear();
1465     nodeIds.push_back(node1->getVtkId());
1466     nodeIds.push_back(node2->getVtkId());
1467     nodeIds.push_back(node3->getVtkId());
1468
1469     SMDS_MeshFace * face = 0;
1470     SMDS_VtkFace *facevtk = myFacePool->getNew();
1471     facevtk->init(nodeIds, this); // put in vtkUnstructuredGrid
1472     if (!this->registerElement(ID,facevtk))
1473       {
1474         this->myGrid->GetCellTypesArray()->SetValue(facevtk->getVtkId(), VTK_EMPTY_CELL);
1475         myFacePool->destroy(facevtk);
1476         return 0;
1477       }
1478     face = facevtk;
1479     adjustmyCellsCapacity(ID);
1480     myCells[ID] = face;
1481     //MESSAGE("createTriangle " << ID << " " << face);
1482     myInfo.myNbTriangles++;
1483     return face;
1484   }
1485 }
1486
1487 ///////////////////////////////////////////////////////////////////////////////
1488 ///Create a quadrangle and add it to the current mesh. This methode do not bind
1489 ///a ID to the create triangle.
1490 ///////////////////////////////////////////////////////////////////////////////
1491 SMDS_MeshFace * SMDS_Mesh::createQuadrangle(const SMDS_MeshNode * node1,
1492                                             const SMDS_MeshNode * node2,
1493                                             const SMDS_MeshNode * node3,
1494                                             const SMDS_MeshNode * node4,
1495                                             int ID)
1496 {
1497   if ( !node1 || !node2 || !node3 || !node4 ) return 0;
1498 //  if ( myFaces.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
1499   if(hasConstructionEdges())
1500   {
1501       //MESSAGE("createQuadrangle hasConstructionEdges "<< ID);
1502     SMDS_MeshEdge *edge1, *edge2, *edge3, *edge4;
1503     edge1=FindEdgeOrCreate(node1,node2);
1504     edge2=FindEdgeOrCreate(node2,node3);
1505     edge3=FindEdgeOrCreate(node3,node4);
1506     edge4=FindEdgeOrCreate(node4,node1);
1507
1508     SMDS_MeshFace * face = new SMDS_FaceOfEdges(edge1,edge2,edge3,edge4);
1509     adjustmyCellsCapacity(ID);
1510     myCells[ID] = face;
1511     myInfo.myNbQuadrangles++;
1512     return face;
1513   }
1514   else
1515   {
1516     // --- retrieve nodes ID
1517     vector<vtkIdType> nodeIds;
1518     nodeIds.clear();
1519     nodeIds.push_back(node1->getVtkId());
1520     nodeIds.push_back(node2->getVtkId());
1521     nodeIds.push_back(node3->getVtkId());
1522     nodeIds.push_back(node4->getVtkId());
1523
1524     SMDS_MeshFace * face = 0;
1525     SMDS_VtkFace *facevtk = myFacePool->getNew();
1526     facevtk->init(nodeIds, this);
1527     if (!this->registerElement(ID,facevtk))
1528       {
1529         this->myGrid->GetCellTypesArray()->SetValue(facevtk->getVtkId(), VTK_EMPTY_CELL);
1530         myFacePool->destroy(facevtk);
1531         return 0;
1532       }
1533     face = facevtk;
1534     adjustmyCellsCapacity(ID);
1535     myCells[ID] = face;
1536     myInfo.myNbQuadrangles++;
1537     return face;
1538   }
1539 }
1540
1541 ///////////////////////////////////////////////////////////////////////////////
1542 /// Remove a node and all the elements which own this node
1543 ///////////////////////////////////////////////////////////////////////////////
1544
1545 void SMDS_Mesh::RemoveNode(const SMDS_MeshNode * node)
1546 {
1547     MESSAGE("RemoveNode");
1548         RemoveElement(node, true);
1549 }
1550
1551 ///////////////////////////////////////////////////////////////////////////////
1552 /// Remove an edge and all the elements which own this edge
1553 ///////////////////////////////////////////////////////////////////////////////
1554
1555 void SMDS_Mesh::Remove0DElement(const SMDS_Mesh0DElement * elem0d)
1556 {
1557     MESSAGE("Remove0DElement");
1558   RemoveElement(elem0d,true);
1559 }
1560
1561 ///////////////////////////////////////////////////////////////////////////////
1562 /// Remove an edge and all the elements which own this edge
1563 ///////////////////////////////////////////////////////////////////////////////
1564
1565 void SMDS_Mesh::RemoveEdge(const SMDS_MeshEdge * edge)
1566 {
1567     MESSAGE("RemoveEdge");
1568         RemoveElement(edge,true);
1569 }
1570
1571 ///////////////////////////////////////////////////////////////////////////////
1572 /// Remove an face and all the elements which own this face
1573 ///////////////////////////////////////////////////////////////////////////////
1574
1575 void SMDS_Mesh::RemoveFace(const SMDS_MeshFace * face)
1576 {
1577     MESSAGE("RemoveFace");
1578         RemoveElement(face, true);
1579 }
1580
1581 ///////////////////////////////////////////////////////////////////////////////
1582 /// Remove a volume
1583 ///////////////////////////////////////////////////////////////////////////////
1584
1585 void SMDS_Mesh::RemoveVolume(const SMDS_MeshVolume * volume)
1586 {
1587     MESSAGE("RemoveVolume");
1588         RemoveElement(volume, true);
1589 }
1590
1591 //=======================================================================
1592 //function : RemoveFromParent
1593 //purpose  :
1594 //=======================================================================
1595
1596 bool SMDS_Mesh::RemoveFromParent()
1597 {
1598         if (myParent==NULL) return false;
1599         else return (myParent->RemoveSubMesh(this));
1600 }
1601
1602 //=======================================================================
1603 //function : RemoveSubMesh
1604 //purpose  :
1605 //=======================================================================
1606
1607 bool SMDS_Mesh::RemoveSubMesh(const SMDS_Mesh * aMesh)
1608 {
1609         bool found = false;
1610
1611         list<SMDS_Mesh *>::iterator itmsh=myChildren.begin();
1612         for (; itmsh!=myChildren.end() && !found; itmsh++)
1613         {
1614                 SMDS_Mesh * submesh = *itmsh;
1615                 if (submesh == aMesh)
1616                 {
1617                         found = true;
1618                         myChildren.erase(itmsh);
1619                 }
1620         }
1621
1622         return found;
1623 }
1624
1625 //=======================================================================
1626 //function : ChangeElementNodes
1627 //purpose  :
1628 //=======================================================================
1629
1630 bool SMDS_Mesh::ChangeElementNodes(const SMDS_MeshElement * element,
1631                                    const SMDS_MeshNode    * nodes[],
1632                                    const int                nbnodes)
1633 {
1634   MESSAGE("SMDS_Mesh::ChangeElementNodes");
1635   // keep current nodes of elem
1636   set<const SMDS_MeshElement*> oldNodes;
1637   SMDS_ElemIteratorPtr itn = element->nodesIterator();
1638   while(itn->more())
1639     oldNodes.insert(itn->next());
1640
1641   // change nodes
1642   bool Ok = false;
1643   SMDS_MeshCell* cell = dynamic_cast<SMDS_MeshCell*>((SMDS_MeshElement*) element);
1644   if (cell)
1645     {
1646       Ok = cell->vtkOrder(nodes, nbnodes);
1647       Ok = cell->ChangeNodes(nodes, nbnodes);
1648     }
1649
1650   if ( Ok ) { // update InverseElements
1651
1652     set<const SMDS_MeshElement*>::iterator it;
1653
1654     // AddInverseElement to new nodes
1655     for ( int i = 0; i < nbnodes; i++ ) {
1656       it = oldNodes.find( nodes[i] );
1657       if ( it == oldNodes.end() )
1658         // new node
1659         const_cast<SMDS_MeshNode*>( nodes[i] )->AddInverseElement( cell );
1660       else
1661         // remove from oldNodes a node that remains in elem
1662         oldNodes.erase( it );
1663     }
1664     // RemoveInverseElement from the nodes removed from elem
1665     for ( it = oldNodes.begin(); it != oldNodes.end(); it++ )
1666     {
1667       SMDS_MeshNode * n = static_cast<SMDS_MeshNode *>
1668         (const_cast<SMDS_MeshElement *>( *it ));
1669       n->RemoveInverseElement( cell );
1670     }
1671   }
1672
1673   return Ok;
1674 }
1675
1676 //=======================================================================
1677 //function : ChangePolyhedronNodes
1678 //purpose  : to change nodes of polyhedral volume
1679 //=======================================================================
1680 bool SMDS_Mesh::ChangePolyhedronNodes (const SMDS_MeshElement *            elem,
1681                                        const vector<const SMDS_MeshNode*>& nodes,
1682                                        const vector<int>                 & quantities)
1683 {
1684   if (elem->GetType() != SMDSAbs_Volume) {
1685     MESSAGE("WRONG ELEM TYPE");
1686     return false;
1687   }
1688
1689   const SMDS_VtkVolume* vol = dynamic_cast<const SMDS_VtkVolume*>(elem);
1690   if (!vol) {
1691     return false;
1692   }
1693
1694   // keep current nodes of elem
1695   set<const SMDS_MeshElement*> oldNodes;
1696   SMDS_ElemIteratorPtr itn = elem->nodesIterator();
1697   while (itn->more()) {
1698     oldNodes.insert(itn->next());
1699   }
1700
1701   // change nodes
1702   // TODO remove this function
1703   //bool Ok = const_cast<SMDS_VtkVolume*>(vol)->ChangeNodes(nodes, quantities);
1704   bool Ok = false;
1705   if (!Ok) {
1706     return false;
1707   }
1708
1709   // update InverseElements
1710
1711   // AddInverseElement to new nodes
1712   int nbnodes = nodes.size();
1713   set<const SMDS_MeshElement*>::iterator it;
1714   for (int i = 0; i < nbnodes; i++) {
1715     it = oldNodes.find(nodes[i]);
1716     if (it == oldNodes.end()) {
1717       // new node
1718       const_cast<SMDS_MeshNode*>(nodes[i])->AddInverseElement(elem);
1719     } else {
1720       // remove from oldNodes a node that remains in elem
1721       oldNodes.erase(it);
1722     }
1723   }
1724
1725   // RemoveInverseElement from the nodes removed from elem
1726   for (it = oldNodes.begin(); it != oldNodes.end(); it++) {
1727     SMDS_MeshNode * n = static_cast<SMDS_MeshNode *>
1728       (const_cast<SMDS_MeshElement *>( *it ));
1729     n->RemoveInverseElement(elem);
1730   }
1731
1732   return Ok;
1733 }
1734
1735
1736 //=======================================================================
1737 //function : Find0DElement
1738 //purpose  :
1739 //=======================================================================
1740 const SMDS_Mesh0DElement* SMDS_Mesh::Find0DElement(int idnode) const
1741 {
1742   const SMDS_MeshNode * node = FindNode(idnode);
1743   if(node == NULL) return NULL;
1744   return Find0DElement(node);
1745 }
1746
1747 const SMDS_Mesh0DElement* SMDS_Mesh::Find0DElement(const SMDS_MeshNode * node)
1748 {
1749   if (!node) return 0;
1750   const SMDS_Mesh0DElement* toReturn = NULL;
1751   SMDS_ElemIteratorPtr it1 = node->GetInverseElementIterator(SMDSAbs_0DElement);
1752   while (it1->more() && (toReturn == NULL)) {
1753     const SMDS_MeshElement* e = it1->next();
1754     if (e->NbNodes() == 1) {
1755       toReturn = static_cast<const SMDS_Mesh0DElement*>(e);
1756     }
1757   }
1758   return toReturn;
1759 }
1760
1761 //=======================================================================
1762 //function : Find0DElementOrCreate
1763 //purpose  :
1764 //=======================================================================
1765 //SMDS_Mesh0DElement* SMDS_Mesh::Find0DElementOrCreate(const SMDS_MeshNode * node)
1766 //{
1767 //  if (!node) return 0;
1768 //  SMDS_Mesh0DElement * toReturn = NULL;
1769 //  toReturn = const_cast<SMDS_Mesh0DElement*>(Find0DElement(node));
1770 //  if (toReturn == NULL) {
1771 //    //if (my0DElements.Extent() % CHECKMEMORY_INTERVAL == 0) CheckMemory();
1772 //    toReturn = new SMDS_Mesh0DElement(node);
1773 //    my0DElements.Add(toReturn);
1774 //    myInfo.myNb0DElements++;
1775 //  }
1776 //  return toReturn;
1777 //}
1778
1779
1780 //=======================================================================
1781 //function : FindEdge
1782 //purpose  :
1783 //=======================================================================
1784
1785 const SMDS_MeshEdge* SMDS_Mesh::FindEdge(int idnode1, int idnode2) const
1786 {
1787   const SMDS_MeshNode * node1=FindNode(idnode1);
1788   const SMDS_MeshNode * node2=FindNode(idnode2);
1789   if((node1==NULL)||(node2==NULL)) return NULL;
1790   return FindEdge(node1,node2);
1791 }
1792
1793 //#include "Profiler.h"
1794 const SMDS_MeshEdge* SMDS_Mesh::FindEdge(const SMDS_MeshNode * node1,
1795                                          const SMDS_MeshNode * node2)
1796 {
1797   if ( !node1 ) return 0;
1798   const SMDS_MeshEdge * toReturn=NULL;
1799   //PROFILER_Init();
1800   //PROFILER_Set();
1801   SMDS_ElemIteratorPtr it1=node1->GetInverseElementIterator(SMDSAbs_Edge);
1802   //PROFILER_Get(0);
1803   //PROFILER_Set();
1804   while(it1->more()) {
1805     const SMDS_MeshElement * e = it1->next();
1806     if ( e->NbNodes() == 2 && e->GetNodeIndex( node2 ) >= 0 ) {
1807       toReturn = static_cast<const SMDS_MeshEdge*>( e );
1808       break;
1809     }
1810   }
1811   //PROFILER_Get(1);
1812   return toReturn;
1813 }
1814
1815
1816 //=======================================================================
1817 //function : FindEdgeOrCreate
1818 //purpose  :
1819 //=======================================================================
1820
1821 SMDS_MeshEdge* SMDS_Mesh::FindEdgeOrCreate(const SMDS_MeshNode * node1,
1822                                            const SMDS_MeshNode * node2)
1823 {
1824   if ( !node1 || !node2) return 0;
1825   SMDS_MeshEdge * toReturn=NULL;
1826   toReturn=const_cast<SMDS_MeshEdge*>(FindEdge(node1,node2));
1827   if(toReturn==NULL) {
1828     //if ( myEdges.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
1829     int ID = myElementIDFactory->GetFreeID(); // -PR- voir si on range cet element
1830     adjustmyCellsCapacity(ID);
1831     vector<vtkIdType> nodeIds;
1832     nodeIds.clear();
1833     nodeIds.push_back(node1->getVtkId());
1834     nodeIds.push_back(node2->getVtkId());
1835
1836     SMDS_VtkEdge *edgevtk = myEdgePool->getNew();
1837     edgevtk->init(nodeIds, this);
1838     if (!this->registerElement(ID,edgevtk))
1839       {
1840         this->myGrid->GetCellTypesArray()->SetValue(edgevtk->getVtkId(), VTK_EMPTY_CELL);
1841         myEdgePool->destroy(edgevtk);
1842         return 0;
1843       }
1844     toReturn = edgevtk;
1845     myCells[ID] = toReturn;
1846     myInfo.myNbEdges++;
1847   }
1848   return toReturn;
1849 }
1850
1851
1852 //=======================================================================
1853 //function : FindEdge
1854 //purpose  :
1855 //=======================================================================
1856
1857 const SMDS_MeshEdge* SMDS_Mesh::FindEdge(int idnode1, int idnode2,
1858                                          int idnode3) const
1859 {
1860   const SMDS_MeshNode * node1=FindNode(idnode1);
1861   const SMDS_MeshNode * node2=FindNode(idnode2);
1862   const SMDS_MeshNode * node3=FindNode(idnode3);
1863   return FindEdge(node1,node2,node3);
1864 }
1865
1866 const SMDS_MeshEdge* SMDS_Mesh::FindEdge(const SMDS_MeshNode * node1,
1867                                          const SMDS_MeshNode * node2,
1868                                          const SMDS_MeshNode * node3)
1869 {
1870   if ( !node1 ) return 0;
1871   SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Edge);
1872   while(it1->more()) {
1873     const SMDS_MeshElement * e = it1->next();
1874     if ( e->NbNodes() == 3 ) {
1875       SMDS_ElemIteratorPtr it2 = e->nodesIterator();
1876       while(it2->more()) {
1877         const SMDS_MeshElement* n = it2->next();
1878         if( n!=node1 &&
1879             n!=node2 &&
1880             n!=node3 )
1881         {
1882           e = 0;
1883           break;
1884         }
1885       }
1886       if ( e )
1887         return static_cast<const SMDS_MeshEdge *> (e);
1888     }
1889   }
1890   return 0;
1891 }
1892
1893
1894 //=======================================================================
1895 //function : FindFace
1896 //purpose  :
1897 //=======================================================================
1898
1899 const SMDS_MeshFace* SMDS_Mesh::FindFace(int idnode1, int idnode2,
1900         int idnode3) const
1901 {
1902   const SMDS_MeshNode * node1=FindNode(idnode1);
1903   const SMDS_MeshNode * node2=FindNode(idnode2);
1904   const SMDS_MeshNode * node3=FindNode(idnode3);
1905   return FindFace(node1, node2, node3);
1906 }
1907
1908 const SMDS_MeshFace* SMDS_Mesh::FindFace(const SMDS_MeshNode *node1,
1909                                          const SMDS_MeshNode *node2,
1910                                          const SMDS_MeshNode *node3)
1911 {
1912   if ( !node1 ) return 0;
1913   SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Face);
1914   while(it1->more()) {
1915     const SMDS_MeshElement * e = it1->next();
1916     if ( e->NbNodes() == 3 ) {
1917       SMDS_ElemIteratorPtr it2 = e->nodesIterator();
1918       while(it2->more()) {
1919         const SMDS_MeshElement* n = it2->next();
1920         if( n!=node1 &&
1921             n!=node2 &&
1922             n!=node3 )
1923         {
1924           e = 0;
1925           break;
1926         }
1927       }
1928       if ( e )
1929         return static_cast<const SMDS_MeshFace *> (e);
1930     }
1931   }
1932   return 0;
1933 }
1934
1935 SMDS_MeshFace* SMDS_Mesh::FindFaceOrCreate(const SMDS_MeshNode *node1,
1936                                            const SMDS_MeshNode *node2,
1937                                            const SMDS_MeshNode *node3)
1938 {
1939   SMDS_MeshFace * toReturn=NULL;
1940   toReturn = const_cast<SMDS_MeshFace*>(FindFace(node1,node2,node3));
1941   if(toReturn==NULL) {
1942     int ID = myElementIDFactory->GetFreeID();
1943     toReturn = createTriangle(node1,node2,node3, ID);
1944   }
1945   return toReturn;
1946 }
1947
1948
1949 //=======================================================================
1950 //function : FindFace
1951 //purpose  :
1952 //=======================================================================
1953
1954 const SMDS_MeshFace* SMDS_Mesh::FindFace(int idnode1, int idnode2,
1955                                          int idnode3, int idnode4) const
1956 {
1957   const SMDS_MeshNode * node1=FindNode(idnode1);
1958   const SMDS_MeshNode * node2=FindNode(idnode2);
1959   const SMDS_MeshNode * node3=FindNode(idnode3);
1960   const SMDS_MeshNode * node4=FindNode(idnode4);
1961   return FindFace(node1, node2, node3, node4);
1962 }
1963
1964 const SMDS_MeshFace* SMDS_Mesh::FindFace(const SMDS_MeshNode *node1,
1965                                          const SMDS_MeshNode *node2,
1966                                          const SMDS_MeshNode *node3,
1967                                          const SMDS_MeshNode *node4)
1968 {
1969   if ( !node1 ) return 0;
1970   SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Face);
1971   while(it1->more()) {
1972     const SMDS_MeshElement * e = it1->next();
1973     if ( e->NbNodes() == 4 ) {
1974       SMDS_ElemIteratorPtr it2 = e->nodesIterator();
1975       while(it2->more()) {
1976         const SMDS_MeshElement* n = it2->next();
1977         if( n!=node1 &&
1978             n!=node2 &&
1979             n!=node3 &&
1980             n!=node4 )
1981         {
1982           e = 0;
1983           break;
1984         }
1985       }
1986       if ( e )
1987         return static_cast<const SMDS_MeshFace *> (e);
1988     }
1989   }
1990   return 0;
1991 }
1992
1993 SMDS_MeshFace* SMDS_Mesh::FindFaceOrCreate(const SMDS_MeshNode *node1,
1994                                            const SMDS_MeshNode *node2,
1995                                            const SMDS_MeshNode *node3,
1996                                            const SMDS_MeshNode *node4)
1997 {
1998   SMDS_MeshFace * toReturn=NULL;
1999   toReturn=const_cast<SMDS_MeshFace*>(FindFace(node1,node2,node3,node4));
2000   if(toReturn==NULL) {
2001     int ID = myElementIDFactory->GetFreeID();
2002     toReturn=createQuadrangle(node1,node2,node3,node4,ID);
2003   }
2004   return toReturn;
2005 }
2006
2007
2008 //=======================================================================
2009 //function : FindFace
2010 //purpose  :quadratic triangle
2011 //=======================================================================
2012
2013 const SMDS_MeshFace* SMDS_Mesh::FindFace(int idnode1, int idnode2,
2014                                          int idnode3, int idnode4,
2015                                          int idnode5, int idnode6) const
2016 {
2017   const SMDS_MeshNode * node1 = FindNode(idnode1);
2018   const SMDS_MeshNode * node2 = FindNode(idnode2);
2019   const SMDS_MeshNode * node3 = FindNode(idnode3);
2020   const SMDS_MeshNode * node4 = FindNode(idnode4);
2021   const SMDS_MeshNode * node5 = FindNode(idnode5);
2022   const SMDS_MeshNode * node6 = FindNode(idnode6);
2023   return FindFace(node1, node2, node3, node4, node5, node6);
2024 }
2025
2026 const SMDS_MeshFace* SMDS_Mesh::FindFace(const SMDS_MeshNode *node1,
2027                                          const SMDS_MeshNode *node2,
2028                                          const SMDS_MeshNode *node3,
2029                                          const SMDS_MeshNode *node4,
2030                                          const SMDS_MeshNode *node5,
2031                                          const SMDS_MeshNode *node6)
2032 {
2033   if ( !node1 ) return 0;
2034   SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Face);
2035   while(it1->more()) {
2036     const SMDS_MeshElement * e = it1->next();
2037     if ( e->NbNodes() == 6 ) {
2038       SMDS_ElemIteratorPtr it2 = e->nodesIterator();
2039       while(it2->more()) {
2040         const SMDS_MeshElement* n = it2->next();
2041         if( n!=node1 &&
2042             n!=node2 &&
2043             n!=node3 &&
2044             n!=node4 &&
2045             n!=node5 &&
2046             n!=node6 )
2047         {
2048           e = 0;
2049           break;
2050         }
2051       }
2052       if ( e )
2053         return static_cast<const SMDS_MeshFace *> (e);
2054     }
2055   }
2056   return 0;
2057 }
2058
2059
2060 //=======================================================================
2061 //function : FindFace
2062 //purpose  : quadratic quadrangle
2063 //=======================================================================
2064
2065 const SMDS_MeshFace* SMDS_Mesh::FindFace(int idnode1, int idnode2,
2066                                          int idnode3, int idnode4,
2067                                          int idnode5, int idnode6,
2068                                          int idnode7, int idnode8) const
2069 {
2070   const SMDS_MeshNode * node1 = FindNode(idnode1);
2071   const SMDS_MeshNode * node2 = FindNode(idnode2);
2072   const SMDS_MeshNode * node3 = FindNode(idnode3);
2073   const SMDS_MeshNode * node4 = FindNode(idnode4);
2074   const SMDS_MeshNode * node5 = FindNode(idnode5);
2075   const SMDS_MeshNode * node6 = FindNode(idnode6);
2076   const SMDS_MeshNode * node7 = FindNode(idnode7);
2077   const SMDS_MeshNode * node8 = FindNode(idnode8);
2078   return FindFace(node1, node2, node3, node4, node5, node6, node7, node8);
2079 }
2080
2081 const SMDS_MeshFace* SMDS_Mesh::FindFace(const SMDS_MeshNode *node1,
2082                                          const SMDS_MeshNode *node2,
2083                                          const SMDS_MeshNode *node3,
2084                                          const SMDS_MeshNode *node4,
2085                                          const SMDS_MeshNode *node5,
2086                                          const SMDS_MeshNode *node6,
2087                                          const SMDS_MeshNode *node7,
2088                                          const SMDS_MeshNode *node8)
2089 {
2090   if ( !node1 ) return 0;
2091   SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Face);
2092   while(it1->more()) {
2093     const SMDS_MeshElement * e = it1->next();
2094     if ( e->NbNodes() == 8 ) {
2095       SMDS_ElemIteratorPtr it2 = e->nodesIterator();
2096       while(it2->more()) {
2097         const SMDS_MeshElement* n = it2->next();
2098         if( n!=node1 &&
2099             n!=node2 &&
2100             n!=node3 &&
2101             n!=node4 &&
2102             n!=node5 &&
2103             n!=node6 &&
2104             n!=node7 &&
2105             n!=node8 )
2106         {
2107           e = 0;
2108           break;
2109         }
2110       }
2111       if ( e )
2112         return static_cast<const SMDS_MeshFace *> (e);
2113     }
2114   }
2115   return 0;
2116 }
2117
2118
2119 //=======================================================================
2120 //function : FindElement
2121 //purpose  :
2122 //=======================================================================
2123
2124 const SMDS_MeshElement* SMDS_Mesh::FindElement(int IDelem) const
2125 {
2126   if ((IDelem <= 0) || IDelem >= myCells.size())
2127   {
2128     MESSAGE("--------------------------------------------------------------------------------- ");
2129     MESSAGE("----------------------------------- bad IDelem " << IDelem << " " << myCells.size());
2130     MESSAGE("--------------------------------------------------------------------------------- ");
2131     // TODO raise an exception
2132     //assert(0);
2133     return 0;
2134   }
2135   return myCells[IDelem];
2136 }
2137
2138 //=======================================================================
2139 //function : FindFace
2140 //purpose  : find polygon
2141 //=======================================================================
2142
2143 const SMDS_MeshFace* SMDS_Mesh::FindFace (const vector<int>& nodes_ids) const
2144 {
2145   int nbnodes = nodes_ids.size();
2146   vector<const SMDS_MeshNode *> poly_nodes (nbnodes);
2147   for (int inode = 0; inode < nbnodes; inode++) {
2148     const SMDS_MeshNode * node = FindNode(nodes_ids[inode]);
2149     if (node == NULL) return NULL;
2150     poly_nodes[inode] = node;
2151   }
2152   return FindFace(poly_nodes);
2153 }
2154
2155 const SMDS_MeshFace* SMDS_Mesh::FindFace (const vector<const SMDS_MeshNode *>& nodes)
2156 {
2157   return (const SMDS_MeshFace*) FindElement( nodes, SMDSAbs_Face );
2158 }
2159
2160
2161 //================================================================================
2162 /*!
2163  * \brief Return element based on all given nodes
2164  *  \param nodes - node of element
2165  *  \param type - type of element
2166  *  \param noMedium - true if medium nodes of quadratic element are not included in <nodes>
2167  *  \retval const SMDS_MeshElement* - found element or NULL
2168  */
2169 //================================================================================
2170
2171 const SMDS_MeshElement* SMDS_Mesh::FindElement (const vector<const SMDS_MeshNode *>& nodes,
2172                                                 const SMDSAbs_ElementType            type,
2173                                                 const bool                           noMedium)
2174 {
2175   if ( nodes.size() > 0 && nodes[0] )
2176   {
2177     SMDS_ElemIteratorPtr itF = nodes[0]->GetInverseElementIterator(type);
2178     while (itF->more())
2179     {
2180       const SMDS_MeshElement* e = itF->next();
2181       int nbNodesToCheck = noMedium ? e->NbCornerNodes() : e->NbNodes();
2182       if ( nbNodesToCheck == nodes.size() )
2183       {
2184         for ( int i = 1; e && i < nodes.size(); ++ i )
2185         {
2186           int nodeIndex = e->GetNodeIndex( nodes[ i ]);
2187           if ( nodeIndex < 0 || nodeIndex >= nbNodesToCheck )
2188             e = 0;
2189         }
2190         if ( e )
2191           return static_cast<const SMDS_MeshFace *> (e);
2192       }
2193     }
2194   }
2195   return NULL;
2196 }
2197
2198 //=======================================================================
2199 //function : DumpNodes
2200 //purpose  :
2201 //=======================================================================
2202
2203 void SMDS_Mesh::DumpNodes() const
2204 {
2205         MESSAGE("dump nodes of mesh : ");
2206         SMDS_NodeIteratorPtr itnode=nodesIterator();
2207         while(itnode->more()) ; //MESSAGE(itnode->next());
2208 }
2209
2210 //=======================================================================
2211 //function : Dump0DElements
2212 //purpose  :
2213 //=======================================================================
2214 void SMDS_Mesh::Dump0DElements() const
2215 {
2216   MESSAGE("dump 0D elements of mesh : ");
2217   SMDS_0DElementIteratorPtr it0d = elements0dIterator();
2218   while(it0d->more()) ; //MESSAGE(it0d->next());
2219 }
2220
2221 //=======================================================================
2222 //function : DumpEdges
2223 //purpose  :
2224 //=======================================================================
2225
2226 void SMDS_Mesh::DumpEdges() const
2227 {
2228         MESSAGE("dump edges of mesh : ");
2229         SMDS_EdgeIteratorPtr itedge=edgesIterator();
2230         while(itedge->more()) ; //MESSAGE(itedge->next());
2231 }
2232
2233 //=======================================================================
2234 //function : DumpFaces
2235 //purpose  :
2236 //=======================================================================
2237
2238 void SMDS_Mesh::DumpFaces() const
2239 {
2240         MESSAGE("dump faces of mesh : ");
2241         SMDS_FaceIteratorPtr itface=facesIterator();
2242         while(itface->more()) ; //MESSAGE(itface->next());
2243 }
2244
2245 //=======================================================================
2246 //function : DumpVolumes
2247 //purpose  :
2248 //=======================================================================
2249
2250 void SMDS_Mesh::DumpVolumes() const
2251 {
2252         MESSAGE("dump volumes of mesh : ");
2253         SMDS_VolumeIteratorPtr itvol=volumesIterator();
2254         while(itvol->more()) ; //MESSAGE(itvol->next());
2255 }
2256
2257 //=======================================================================
2258 //function : DebugStats
2259 //purpose  :
2260 //=======================================================================
2261
2262 void SMDS_Mesh::DebugStats() const
2263 {
2264   MESSAGE("Debug stats of mesh : ");
2265
2266   MESSAGE("===== NODES ====="<<NbNodes());
2267   MESSAGE("===== 0DELEMS ====="<<Nb0DElements());
2268   MESSAGE("===== EDGES ====="<<NbEdges());
2269   MESSAGE("===== FACES ====="<<NbFaces());
2270   MESSAGE("===== VOLUMES ====="<<NbVolumes());
2271
2272   MESSAGE("End Debug stats of mesh ");
2273
2274   //#ifdef DEB
2275
2276   SMDS_NodeIteratorPtr itnode=nodesIterator();
2277   int sizeofnodes = 0;
2278   int sizeoffaces = 0;
2279
2280   while(itnode->more())
2281   {
2282     const SMDS_MeshNode *node = itnode->next();
2283
2284     sizeofnodes += sizeof(*node);
2285
2286     SMDS_ElemIteratorPtr it = node->GetInverseElementIterator();
2287     while(it->more())
2288     {
2289       const SMDS_MeshElement *me = it->next();
2290       sizeofnodes += sizeof(me);
2291     }
2292   }
2293
2294   SMDS_FaceIteratorPtr itface=facesIterator();
2295   while(itface->more())
2296   {
2297     const SMDS_MeshElement *face = itface->next();
2298     sizeoffaces += sizeof(*face);
2299   }
2300
2301   MESSAGE("total size of node elements = " << sizeofnodes);;
2302   MESSAGE("total size of face elements = " << sizeoffaces);;
2303
2304   //#endif
2305 }
2306
2307 ///////////////////////////////////////////////////////////////////////////////
2308 /// Return the number of nodes
2309 ///////////////////////////////////////////////////////////////////////////////
2310 int SMDS_Mesh::NbNodes() const
2311 {
2312         //MESSAGE(myGrid->GetNumberOfPoints());
2313         //MESSAGE(myInfo.NbNodes());
2314         //MESSAGE(myNodeMax);
2315     return myInfo.NbNodes();
2316 }
2317
2318 ///////////////////////////////////////////////////////////////////////////////
2319 /// Return the number of 0D elements
2320 ///////////////////////////////////////////////////////////////////////////////
2321 int SMDS_Mesh::Nb0DElements() const
2322 {
2323   return myInfo.Nb0DElements(); // -PR- a verfier
2324 }
2325
2326 ///////////////////////////////////////////////////////////////////////////////
2327 /// Return the number of edges (including construction edges)
2328 ///////////////////////////////////////////////////////////////////////////////
2329 int SMDS_Mesh::NbEdges() const
2330 {
2331         return myInfo.NbEdges(); // -PR- a verfier
2332 }
2333
2334 ///////////////////////////////////////////////////////////////////////////////
2335 /// Return the number of faces (including construction faces)
2336 ///////////////////////////////////////////////////////////////////////////////
2337 int SMDS_Mesh::NbFaces() const
2338 {
2339         return myInfo.NbFaces();  // -PR- a verfier
2340 }
2341
2342 ///////////////////////////////////////////////////////////////////////////////
2343 /// Return the number of volumes
2344 ///////////////////////////////////////////////////////////////////////////////
2345 int SMDS_Mesh::NbVolumes() const
2346 {
2347         return myInfo.NbVolumes(); // -PR- a verfier
2348 }
2349
2350 ///////////////////////////////////////////////////////////////////////////////
2351 /// Return the number of child mesh of this mesh.
2352 /// Note that the tree structure of SMDS_Mesh seems to be unused in this version
2353 /// (2003-09-08) of SMESH
2354 ///////////////////////////////////////////////////////////////////////////////
2355 int SMDS_Mesh::NbSubMesh() const
2356 {
2357         return myChildren.size();
2358 }
2359
2360 ///////////////////////////////////////////////////////////////////////////////
2361 /// Destroy the mesh and all its elements
2362 /// All pointer on elements owned by this mesh become illegals.
2363 ///////////////////////////////////////////////////////////////////////////////
2364 SMDS_Mesh::~SMDS_Mesh()
2365 {
2366   list<SMDS_Mesh*>::iterator itc=myChildren.begin();
2367   while(itc!=myChildren.end())
2368   {
2369     delete *itc;
2370     itc++;
2371   }
2372
2373   if(myParent==NULL)
2374   {
2375     delete myNodeIDFactory;
2376     delete myElementIDFactory;
2377   }
2378   else
2379   {
2380     SMDS_ElemIteratorPtr eIt = elementsIterator();
2381     while ( eIt->more() )
2382       {
2383         const SMDS_MeshElement *elem = eIt->next();
2384         myElementIDFactory->ReleaseID(elem->GetID(), elem->getVtkId());
2385       }
2386     SMDS_NodeIteratorPtr itn = nodesIterator();
2387     while (itn->more())
2388       {
2389         const SMDS_MeshNode *node = itn->next();
2390         myNodeIDFactory->ReleaseID(node->GetID(), node->getVtkId());
2391       }
2392   }
2393
2394 //   SetOfNodes::Iterator itn(myNodes);
2395 //   for (; itn.More(); itn.Next())
2396 //     delete itn.Value();
2397
2398 //   SetOf0DElements::Iterator it0d (my0DElements);
2399 //   for (; it0d.More(); it0d.Next())
2400 //   {
2401 //     SMDS_MeshElement* elem = it0d.Value();
2402 //     delete elem;
2403 //   }
2404
2405 //   SetOfEdges::Iterator ite(myEdges);
2406 //   for (; ite.More(); ite.Next())
2407 //   {
2408 //     SMDS_MeshElement* elem = ite.Value();
2409 //     delete elem;
2410 //   }
2411
2412 //   SetOfFaces::Iterator itf(myFaces);
2413 //   for (; itf.More(); itf.Next())
2414 //   {
2415 //     SMDS_MeshElement* elem = itf.Value();
2416 //     delete elem;
2417 //   }
2418
2419 //   SetOfVolumes::Iterator itv(myVolumes);
2420 //   for (; itv.More(); itv.Next())
2421 //   {
2422 //     SMDS_MeshElement* elem = itv.Value();
2423 //     delete elem;
2424 //   }
2425 }
2426
2427 //================================================================================
2428 /*!
2429  * \brief Clear all data
2430  */
2431 //================================================================================
2432
2433 void SMDS_Mesh::Clear()
2434 {
2435   MESSAGE("SMDS_Mesh::Clear");
2436   if (myParent!=NULL)
2437     {
2438     SMDS_ElemIteratorPtr eIt = elementsIterator();
2439     while ( eIt->more() )
2440       {
2441         const SMDS_MeshElement *elem = eIt->next();
2442         myElementIDFactory->ReleaseID(elem->GetID(), elem->getVtkId());
2443       }
2444     SMDS_NodeIteratorPtr itn = nodesIterator();
2445     while (itn->more())
2446       {
2447         const SMDS_MeshNode *node = itn->next();
2448         myNodeIDFactory->ReleaseID(node->GetID(), node->getVtkId());
2449       }
2450     }
2451   else
2452     {
2453     myNodeIDFactory->Clear();
2454     myElementIDFactory->Clear();
2455     }
2456
2457   SMDS_ElemIteratorPtr itv = elementsIterator();
2458   while (itv->more())
2459     {
2460       SMDS_MeshElement* elem = (SMDS_MeshElement*)(itv->next());
2461       SMDSAbs_ElementType aType = elem->GetType();
2462       switch (aType)
2463       {
2464         case SMDSAbs_0DElement:
2465           delete elem;
2466           break;
2467         case SMDSAbs_Edge:
2468            myEdgePool->destroy(static_cast<SMDS_VtkEdge*>(elem));
2469           break;
2470         case SMDSAbs_Face:
2471           myFacePool->destroy(static_cast<SMDS_VtkFace*>(elem));
2472           break;
2473         case SMDSAbs_Volume:
2474           myVolumePool->destroy(static_cast<SMDS_VtkVolume*>(elem));
2475           break;
2476         default:
2477           break;
2478       }
2479     }
2480   myCells.clear();
2481   myCellIdVtkToSmds.clear();
2482   //myCellIdSmdsToVtk.clear();
2483
2484   SMDS_NodeIteratorPtr itn = nodesIterator();
2485   while (itn->more())
2486     {
2487       SMDS_MeshNode *node = (SMDS_MeshNode*)(itn->next());
2488       myNodePool->destroy(node);
2489     }
2490   myNodes.clear();
2491
2492   list<SMDS_Mesh*>::iterator itc=myChildren.begin();
2493   while(itc!=myChildren.end())
2494     (*itc)->Clear();
2495
2496   myModified = false;
2497   xmin = 0; xmax = 0;
2498   ymin = 0; ymax = 0;
2499   zmin = 0; zmax = 0;
2500
2501   myInfo.Clear();
2502
2503   myGrid->Initialize();
2504   myGrid->Allocate();
2505   vtkPoints* points = vtkPoints::New();
2506   points->SetNumberOfPoints(SMDS_Mesh::chunkSize);
2507   myGrid->SetPoints( points );
2508   points->Delete();
2509   myGrid->BuildLinks();
2510 }
2511
2512 ///////////////////////////////////////////////////////////////////////////////
2513 /// Return true if this mesh create faces with edges.
2514 /// A false returned value mean that faces are created with nodes. A concequence
2515 /// is, iteration on edges (SMDS_Element::edgesIterator) will be unavailable.
2516 ///////////////////////////////////////////////////////////////////////////////
2517 bool SMDS_Mesh::hasConstructionEdges()
2518 {
2519         return myHasConstructionEdges;
2520 }
2521
2522 ///////////////////////////////////////////////////////////////////////////////
2523 /// Return true if this mesh create volumes with faces
2524 /// A false returned value mean that volumes are created with nodes or edges.
2525 /// (see hasConstructionEdges)
2526 /// A concequence is, iteration on faces (SMDS_Element::facesIterator) will be
2527 /// unavailable.
2528 ///////////////////////////////////////////////////////////////////////////////
2529 bool SMDS_Mesh::hasConstructionFaces()
2530 {
2531         return myHasConstructionFaces;
2532 }
2533
2534 ///////////////////////////////////////////////////////////////////////////////
2535 /// Return true if nodes are linked to the finit elements, they are belonging to.
2536 /// Currently, It always return true.
2537 ///////////////////////////////////////////////////////////////////////////////
2538 bool SMDS_Mesh::hasInverseElements()
2539 {
2540         return myHasInverseElements;
2541 }
2542
2543 ///////////////////////////////////////////////////////////////////////////////
2544 /// Make this mesh creating construction edges (see hasConstructionEdges)
2545 /// @param b true to have construction edges, else false.
2546 ///////////////////////////////////////////////////////////////////////////////
2547 void SMDS_Mesh::setConstructionEdges(bool b)
2548 {
2549         myHasConstructionEdges=b;
2550 }
2551
2552 ///////////////////////////////////////////////////////////////////////////////
2553 /// Make this mesh creating construction faces (see hasConstructionFaces)
2554 /// @param b true to have construction faces, else false.
2555 ///////////////////////////////////////////////////////////////////////////////
2556 void SMDS_Mesh::setConstructionFaces(bool b)
2557 {
2558          myHasConstructionFaces=b;
2559 }
2560
2561 ///////////////////////////////////////////////////////////////////////////////
2562 /// Make this mesh creating link from nodes to elements (see hasInverseElements)
2563 /// @param b true to link nodes to elements, else false.
2564 ///////////////////////////////////////////////////////////////////////////////
2565 void SMDS_Mesh::setInverseElements(bool b)
2566 {
2567   if(!b) MESSAGE("Error : inverseElement=false not implemented");
2568   myHasInverseElements=b;
2569 }
2570
2571 namespace {
2572
2573 ///////////////////////////////////////////////////////////////////////////////
2574 ///Iterator on NCollection_Map
2575 ///////////////////////////////////////////////////////////////////////////////
2576 template <class MAP, typename ELEM=const SMDS_MeshElement*, class FATHER=SMDS_ElemIterator>
2577 struct MYNode_Map_Iterator: public FATHER
2578 {
2579   int _ctr;
2580   const MAP& _map;
2581   MYNode_Map_Iterator(const MAP& map): _map(map) // map is a std::vector<ELEM>
2582   {
2583       _ctr = 0;
2584   }
2585
2586   bool more()
2587   {
2588       while (_ctr < _map.size())
2589       {
2590           if (_map[_ctr])
2591               return true;
2592           _ctr++;
2593       }
2594           return false;
2595   }
2596
2597   ELEM next()
2598   {
2599     ELEM current = _map[_ctr];
2600     _ctr++;
2601     return current;
2602   }
2603 };
2604
2605 template <class MAP, typename ELEM=const SMDS_MeshElement*, class FATHER=SMDS_ElemIterator>
2606 struct MYElem_Map_Iterator: public FATHER
2607 {
2608   int _ctr;
2609   int _type;
2610   const MAP& _map;
2611   MYElem_Map_Iterator(const MAP& map, int typ): _map(map) // map is a std::vector<ELEM>
2612   {
2613       _ctr = 0;
2614       _type = typ;
2615       while (_ctr < _map.size()) // go to the first valid element
2616       {
2617           if (_map[_ctr])
2618             if ( (_type == SMDSAbs_All) || (_map[_ctr]->GetType() == _type))
2619               break;
2620           _ctr++;
2621       }
2622   }
2623
2624 bool more()
2625   {
2626       while (_ctr < _map.size())
2627       {
2628           if (_map[_ctr])
2629             if ( (_type == SMDSAbs_All) || (_map[_ctr]->GetType() == _type))
2630               return true;
2631           _ctr++;
2632       }
2633           return false;
2634   }
2635
2636   ELEM next()
2637   {
2638     ELEM current = dynamic_cast<ELEM> (_map[_ctr]);
2639     _ctr++;
2640     return current;
2641   }
2642 };
2643
2644 //================================================================================
2645   /*!
2646    * \brief Iterator on elements in id increasing order
2647    */
2648   //================================================================================
2649
2650   template <typename ELEM=const SMDS_MeshElement*>
2651   class IdSortedIterator : public SMDS_Iterator<ELEM>
2652   {
2653     const SMDS_MeshElementIDFactory& myIDFact;
2654     int                              myID, myMaxID, myNbFound, myTotalNb;
2655     SMDSAbs_ElementType              myType;
2656     ELEM                             myElem;
2657
2658   public:
2659     IdSortedIterator(const SMDS_MeshElementIDFactory& fact,
2660                      const SMDSAbs_ElementType        type, // SMDSAbs_All NOT allowed!!! 
2661                      const int                        totalNb)
2662       :myIDFact( fact ),
2663        myID(1), myMaxID( myIDFact.GetMaxID() ),myNbFound(0), myTotalNb( totalNb ),
2664        myType( type ),
2665        myElem(0)
2666     {
2667       next();
2668     }
2669     bool more()
2670     {
2671       return myElem;
2672     }
2673     ELEM next()
2674     {
2675       ELEM current = myElem;
2676
2677       for ( myElem = 0;  !myElem && myNbFound < myTotalNb && myID <= myMaxID; ++myID )
2678         if ((myElem = (ELEM) myIDFact.MeshElement( myID ))
2679             && myElem->GetType() != myType )
2680           myElem = 0;
2681
2682       myNbFound += bool(myElem);
2683
2684       return current;
2685     }
2686   };
2687 }
2688
2689 ///////////////////////////////////////////////////////////////////////////////
2690 /// Return an iterator on nodes of the current mesh factory
2691 ///////////////////////////////////////////////////////////////////////////////
2692
2693 SMDS_NodeIteratorPtr SMDS_Mesh::nodesIterator(bool idInceasingOrder) const
2694 {
2695   typedef MYNode_Map_Iterator
2696     < SetOfNodes, const SMDS_MeshNode*, SMDS_NodeIterator > TIterator;
2697   return SMDS_NodeIteratorPtr( new TIterator(myNodes)); // naturally always sorted by ID
2698
2699 //  typedef IdSortedIterator< const SMDS_MeshNode* >          TSortedIterator;
2700 //  return ( idInceasingOrder ?
2701 //           SMDS_NodeIteratorPtr( new TSortedIterator( *myNodeIDFactory, SMDSAbs_Node, NbNodes())) :
2702 //           SMDS_NodeIteratorPtr( new TIterator(myNodes)));
2703 }
2704
2705 ///////////////////////////////////////////////////////////////////////////////
2706 ///Return an iterator on 0D elements of the current mesh.
2707 ///////////////////////////////////////////////////////////////////////////////
2708
2709 SMDS_0DElementIteratorPtr SMDS_Mesh::elements0dIterator(bool idInceasingOrder) const
2710 {
2711   typedef MYElem_Map_Iterator
2712     < SetOfCells, const SMDS_Mesh0DElement*, SMDS_0DElementIterator > TIterator;
2713   return SMDS_0DElementIteratorPtr(new TIterator(myCells, SMDSAbs_0DElement)); // naturally always sorted by ID
2714
2715 //  typedef MYNCollection_Map_Iterator
2716 //    < SetOf0DElements, const SMDS_Mesh0DElement*, SMDS_0DElementIterator > TIterator;
2717 //  typedef IdSortedIterator< const SMDS_Mesh0DElement* >                    TSortedIterator;
2718 //  return ( idInceasingOrder ?
2719 //           SMDS_0DElementIteratorPtr( new TSortedIterator( *myElementIDFactory,
2720 //                                                           SMDSAbs_0DElement,
2721 //                                                           Nb0DElements() )) :
2722 //           SMDS_0DElementIteratorPtr( new TIterator(my0DElements)));
2723 }
2724
2725 ///////////////////////////////////////////////////////////////////////////////
2726 ///Return an iterator on edges of the current mesh.
2727 ///////////////////////////////////////////////////////////////////////////////
2728
2729 SMDS_EdgeIteratorPtr SMDS_Mesh::edgesIterator(bool idInceasingOrder) const
2730 {
2731   typedef MYElem_Map_Iterator
2732     < SetOfCells, const SMDS_MeshEdge*, SMDS_EdgeIterator > TIterator;
2733   return SMDS_EdgeIteratorPtr(new TIterator(myCells, SMDSAbs_Edge)); // naturally always sorted by ID
2734
2735 //  typedef MYNCollection_Map_Iterator
2736 //    < SetOfEdges, const SMDS_MeshEdge*, SMDS_EdgeIterator > TIterator;
2737 //  typedef IdSortedIterator< const SMDS_MeshEdge* >          TSortedIterator;
2738 //  return ( idInceasingOrder ?
2739 //           SMDS_EdgeIteratorPtr( new TSortedIterator( *myElementIDFactory,
2740 //                                                      SMDSAbs_Edge,
2741 //                                                      NbEdges() )) :
2742 //           SMDS_EdgeIteratorPtr(new TIterator(myEdges)));
2743 }
2744
2745 ///////////////////////////////////////////////////////////////////////////////
2746 ///Return an iterator on faces of the current mesh.
2747 ///////////////////////////////////////////////////////////////////////////////
2748
2749 SMDS_FaceIteratorPtr SMDS_Mesh::facesIterator(bool idInceasingOrder) const
2750 {
2751   typedef MYElem_Map_Iterator
2752     < SetOfCells, const SMDS_MeshFace*, SMDS_FaceIterator > TIterator;
2753   return SMDS_FaceIteratorPtr(new TIterator(myCells, SMDSAbs_Face)); // naturally always sorted by ID
2754
2755 //  typedef MYNCollection_Map_Iterator
2756 //    < SetOfFaces, const SMDS_MeshFace*, SMDS_FaceIterator > TIterator;
2757 //  typedef IdSortedIterator< const SMDS_MeshFace* >          TSortedIterator;
2758 //  return ( idInceasingOrder ?
2759 //           SMDS_FaceIteratorPtr( new TSortedIterator( *myElementIDFactory,
2760 //                                                      SMDSAbs_Face,
2761 //                                                      NbFaces() )) :
2762 //           SMDS_FaceIteratorPtr(new TIterator(myFaces)));
2763 }
2764
2765 ///////////////////////////////////////////////////////////////////////////////
2766 ///Return an iterator on volumes of the current mesh.
2767 ///////////////////////////////////////////////////////////////////////////////
2768
2769 SMDS_VolumeIteratorPtr SMDS_Mesh::volumesIterator(bool idInceasingOrder) const
2770 {
2771   typedef MYElem_Map_Iterator
2772     < SetOfCells, const SMDS_MeshVolume*, SMDS_VolumeIterator > TIterator;
2773   return SMDS_VolumeIteratorPtr(new TIterator(myCells, SMDSAbs_Volume)); // naturally always sorted by ID
2774
2775   //  typedef MYNCollection_Map_Iterator
2776 //    < SetOfVolumes, const SMDS_MeshVolume*, SMDS_VolumeIterator > TIterator;
2777 //  typedef IdSortedIterator< const SMDS_MeshVolume* >              TSortedIterator;
2778 //  return ( idInceasingOrder ?
2779 //           SMDS_VolumeIteratorPtr( new TSortedIterator( *myElementIDFactory,
2780 //                                                        SMDSAbs_Volume,
2781 //                                                        NbVolumes() )) :
2782 //           SMDS_VolumeIteratorPtr(new TIterator(myVolumes)));
2783 }
2784
2785 ///////////////////////////////////////////////////////////////////////////////
2786 /// Return an iterator on elements of the current mesh factory
2787 ///////////////////////////////////////////////////////////////////////////////
2788 SMDS_ElemIteratorPtr SMDS_Mesh::elementsIterator(SMDSAbs_ElementType type) const
2789 {
2790   switch (type) {
2791   case SMDSAbs_All:
2792     return SMDS_ElemIteratorPtr (new MYElem_Map_Iterator< SetOfCells >(myCells, SMDSAbs_All));
2793     break;
2794   case SMDSAbs_Volume:
2795     return SMDS_ElemIteratorPtr (new MYElem_Map_Iterator< SetOfCells >(myCells, SMDSAbs_Volume));
2796   case SMDSAbs_Face:
2797     return SMDS_ElemIteratorPtr (new MYElem_Map_Iterator< SetOfCells >(myCells, SMDSAbs_Face));
2798   case SMDSAbs_Edge:
2799     return SMDS_ElemIteratorPtr (new MYElem_Map_Iterator< SetOfCells >(myCells, SMDSAbs_Edge));
2800   case SMDSAbs_0DElement:
2801     return SMDS_ElemIteratorPtr (new MYElem_Map_Iterator< SetOfCells >(myCells, SMDSAbs_0DElement));
2802   case SMDSAbs_Node:
2803     return SMDS_ElemIteratorPtr (new MYElem_Map_Iterator< SetOfNodes >(myNodes, SMDSAbs_All));
2804     //return myNodeIDFactory->elementsIterator();
2805   default:;
2806   }
2807   return myElementIDFactory->elementsIterator();
2808 }
2809
2810 ///////////////////////////////////////////////////////////////////////////////
2811 /// Do intersection of sets (more than 2)
2812 ///////////////////////////////////////////////////////////////////////////////
2813 static set<const SMDS_MeshElement*> * intersectionOfSets(
2814         set<const SMDS_MeshElement*> vs[], int numberOfSets)
2815 {
2816         set<const SMDS_MeshElement*>* rsetA=new set<const SMDS_MeshElement*>(vs[0]);
2817         set<const SMDS_MeshElement*>* rsetB;
2818
2819         for(int i=0; i<numberOfSets-1; i++)
2820         {
2821                 rsetB=new set<const SMDS_MeshElement*>();
2822                 set_intersection(
2823                         rsetA->begin(), rsetA->end(),
2824                         vs[i+1].begin(), vs[i+1].end(),
2825                         inserter(*rsetB, rsetB->begin()));
2826                 delete rsetA;
2827                 rsetA=rsetB;
2828         }
2829         return rsetA;
2830 }
2831
2832 ///////////////////////////////////////////////////////////////////////////////
2833 /// Return the list of finite elements owning the given element: elements
2834 /// containing all the nodes of the given element, for instance faces and
2835 /// volumes containing a given edge.
2836 ///////////////////////////////////////////////////////////////////////////////
2837 static set<const SMDS_MeshElement*> * getFinitElements(const SMDS_MeshElement * element)
2838 {
2839         int numberOfSets=element->NbNodes();
2840         set<const SMDS_MeshElement*> *initSet = new set<const SMDS_MeshElement*>[numberOfSets];
2841
2842         SMDS_ElemIteratorPtr itNodes=element->nodesIterator();
2843
2844         int i=0;
2845         while(itNodes->more())
2846         {
2847           const SMDS_MeshElement* node = itNodes->next();
2848           MYASSERT(node);
2849                 const SMDS_MeshNode * n=static_cast<const SMDS_MeshNode*>(node);
2850                 SMDS_ElemIteratorPtr itFe = n->GetInverseElementIterator();
2851
2852                 //initSet[i]=set<const SMDS_MeshElement*>();
2853                 while(itFe->more())
2854                 {
2855                   const SMDS_MeshElement* elem = itFe->next();
2856                   MYASSERT(elem);
2857                   initSet[i].insert(elem);
2858
2859                 }
2860
2861                 i++;
2862         }
2863         set<const SMDS_MeshElement*> *retSet=intersectionOfSets(initSet, numberOfSets);
2864         MESSAGE("nb elems " << i << " intersection " << retSet->size());
2865         delete [] initSet;
2866         return retSet;
2867 }
2868
2869 ///////////////////////////////////////////////////////////////////////////////
2870 /// Return the list of nodes used only by the given elements
2871 ///////////////////////////////////////////////////////////////////////////////
2872 static set<const SMDS_MeshElement*> * getExclusiveNodes(
2873         set<const SMDS_MeshElement*>& elements)
2874 {
2875         set<const SMDS_MeshElement*> * toReturn=new set<const SMDS_MeshElement*>();
2876         set<const SMDS_MeshElement*>::iterator itElements=elements.begin();
2877
2878         while(itElements!=elements.end())
2879         {
2880                 SMDS_ElemIteratorPtr itNodes = (*itElements)->nodesIterator();
2881                 itElements++;
2882
2883                 while(itNodes->more())
2884                 {
2885                         const SMDS_MeshNode * n=static_cast<const SMDS_MeshNode*>(itNodes->next());
2886                         SMDS_ElemIteratorPtr itFe = n->GetInverseElementIterator();
2887                         set<const SMDS_MeshElement*> s;
2888                         while(itFe->more())
2889                           s.insert(itFe->next());
2890                         if(s==elements) toReturn->insert(n);
2891                 }
2892         }
2893         return toReturn;
2894 }
2895
2896 ///////////////////////////////////////////////////////////////////////////////
2897 ///Find the children of an element that are made of given nodes
2898 ///@param setOfChildren The set in which matching children will be inserted
2899 ///@param element The element were to search matching children
2900 ///@param nodes The nodes that the children must have to be selected
2901 ///////////////////////////////////////////////////////////////////////////////
2902 void SMDS_Mesh::addChildrenWithNodes(set<const SMDS_MeshElement*>& setOfChildren,
2903                                      const SMDS_MeshElement *      element,
2904                                      set<const SMDS_MeshElement*>& nodes)
2905 {
2906   switch(element->GetType())
2907     {
2908     case SMDSAbs_Node:
2909       MESSAGE("Internal Error: This should not happen");
2910       break;
2911     case SMDSAbs_0DElement:
2912       {
2913       }
2914       break;
2915     case SMDSAbs_Edge:
2916         {
2917                 SMDS_ElemIteratorPtr itn=element->nodesIterator();
2918                 while(itn->more())
2919                 {
2920                         const SMDS_MeshElement * e=itn->next();
2921                         if(nodes.find(e)!=nodes.end())
2922                         {
2923                           setOfChildren.insert(element);
2924                           break;
2925                         }
2926                 }
2927         } break;
2928     case SMDSAbs_Face:
2929         {
2930                 SMDS_ElemIteratorPtr itn=element->nodesIterator();
2931                 while(itn->more())
2932                 {
2933                         const SMDS_MeshElement * e=itn->next();
2934                         if(nodes.find(e)!=nodes.end())
2935                         {
2936                           setOfChildren.insert(element);
2937                           break;
2938                         }
2939                 }
2940                 if(hasConstructionEdges())
2941                 {
2942                         SMDS_ElemIteratorPtr ite=element->edgesIterator();
2943                         while(ite->more())
2944                                 addChildrenWithNodes(setOfChildren, ite->next(), nodes);
2945                 }
2946         } break;
2947     case SMDSAbs_Volume:
2948         {
2949                 if(hasConstructionFaces())
2950                 {
2951                         SMDS_ElemIteratorPtr ite=element->facesIterator();
2952                         while(ite->more())
2953                                 addChildrenWithNodes(setOfChildren, ite->next(), nodes);
2954                 }
2955                 else if(hasConstructionEdges())
2956                 {
2957                         SMDS_ElemIteratorPtr ite=element->edgesIterator();
2958                         while(ite->more())
2959                                 addChildrenWithNodes(setOfChildren, ite->next(), nodes);
2960                 }
2961         }
2962     }
2963 }
2964
2965 ///////////////////////////////////////////////////////////////////////////////
2966 ///@param elem The element to delete
2967 ///@param removenodes if true remaining nodes will be removed
2968 ///////////////////////////////////////////////////////////////////////////////
2969 void SMDS_Mesh::RemoveElement(const SMDS_MeshElement * elem,
2970                               const bool removenodes)
2971 {
2972   list<const SMDS_MeshElement *> removedElems;
2973   list<const SMDS_MeshElement *> removedNodes;
2974   RemoveElement( elem, removedElems, removedNodes, removenodes );
2975 }
2976
2977 ///////////////////////////////////////////////////////////////////////////////
2978 ///@param elem The element to delete
2979 ///@param removedElems to be filled with all removed elements
2980 ///@param removedNodes to be filled with all removed nodes
2981 ///@param removenodes if true remaining nodes will be removed
2982 ///////////////////////////////////////////////////////////////////////////////
2983 void SMDS_Mesh::RemoveElement(const SMDS_MeshElement *        elem,
2984                               list<const SMDS_MeshElement *>& removedElems,
2985                               list<const SMDS_MeshElement *>& removedNodes,
2986                               bool                            removenodes)
2987 {
2988   //MESSAGE("SMDS_Mesh::RemoveElement " << elem->getVtkId() << " " << removenodes);
2989   // get finite elements built on elem
2990   set<const SMDS_MeshElement*> * s1;
2991   if (    (elem->GetType() == SMDSAbs_0DElement)
2992       || ((elem->GetType() == SMDSAbs_Edge) && !hasConstructionEdges())
2993       || ((elem->GetType() == SMDSAbs_Face) && !hasConstructionFaces())
2994       ||  (elem->GetType() == SMDSAbs_Volume) )
2995     {
2996       s1 = new set<const SMDS_MeshElement*> ();
2997       s1->insert(elem);
2998     }
2999   else
3000     s1 = getFinitElements(elem);
3001
3002   // get exclusive nodes (which would become free afterwards)
3003   set<const SMDS_MeshElement*> * s2;
3004   if (elem->GetType() == SMDSAbs_Node) // a node is removed
3005     {
3006       // do not remove nodes except elem
3007       s2 = new set<const SMDS_MeshElement*> ();
3008       s2->insert(elem);
3009       removenodes = true;
3010     }
3011   else
3012     s2 = getExclusiveNodes(*s1);
3013
3014   // form the set of finite and construction elements to remove
3015   set<const SMDS_MeshElement*> s3;
3016   set<const SMDS_MeshElement*>::iterator it = s1->begin();
3017   while (it != s1->end())
3018     {
3019       addChildrenWithNodes(s3, *it, *s2);
3020       s3.insert(*it);
3021       it++;
3022     }
3023   if (elem->GetType() != SMDSAbs_Node)
3024     s3.insert(elem);
3025
3026   // remove finite and construction elements
3027   it = s3.begin();
3028   while (it != s3.end())
3029     {
3030       // Remove element from <InverseElements> of its nodes
3031       SMDS_ElemIteratorPtr itn = (*it)->nodesIterator();
3032       while (itn->more())
3033         {
3034           SMDS_MeshNode * n = static_cast<SMDS_MeshNode *> (const_cast<SMDS_MeshElement *> (itn->next()));
3035           n->RemoveInverseElement((*it));
3036         }
3037       int IdToRemove = (*it)->GetID();
3038       int vtkid = (*it)->getVtkId();
3039       //MESSAGE("elem Id to remove " << IdToRemove << " vtkid " << vtkid <<
3040       //        " vtktype " << (*it)->GetVtkType() << " type " << (*it)->GetType());
3041       switch ((*it)->GetType())
3042       {
3043         case SMDSAbs_Node:
3044           MYASSERT("Internal Error: This should not happen")
3045           ;
3046           break;
3047         case SMDSAbs_0DElement:
3048           if (IdToRemove >= 0)
3049             {
3050               myCells[IdToRemove] = 0; // -PR- ici ou dans myElementIDFactory->ReleaseID ?
3051               myInfo.remove(*it);
3052             }
3053           removedElems.push_back((*it));
3054           myElementIDFactory->ReleaseID(IdToRemove, vtkid);
3055           delete (*it);
3056           break;
3057         case SMDSAbs_Edge:
3058           if (IdToRemove >= 0)
3059             {
3060               myCells[IdToRemove] = 0;
3061               myInfo.RemoveEdge(*it);
3062             }
3063           removedElems.push_back((*it));
3064           myElementIDFactory->ReleaseID(IdToRemove, vtkid);
3065           if (const SMDS_VtkEdge* vtkElem = dynamic_cast<const SMDS_VtkEdge*>(*it))
3066             myEdgePool->destroy((SMDS_VtkEdge*) vtkElem);
3067           else
3068             delete (*it);
3069           break;
3070         case SMDSAbs_Face:
3071           if (IdToRemove >= 0)
3072             {
3073               myCells[IdToRemove] = 0;
3074               myInfo.RemoveFace(*it);
3075             }
3076           removedElems.push_back((*it));
3077           myElementIDFactory->ReleaseID(IdToRemove, vtkid);
3078           if (const SMDS_VtkFace* vtkElem = dynamic_cast<const SMDS_VtkFace*>(*it))
3079             myFacePool->destroy((SMDS_VtkFace*) vtkElem);
3080           else
3081             delete (*it);
3082           break;
3083         case SMDSAbs_Volume:
3084           if (IdToRemove >= 0)
3085             {
3086               myCells[IdToRemove] = 0;
3087               myInfo.RemoveVolume(*it);
3088             }
3089           removedElems.push_back((*it));
3090           myElementIDFactory->ReleaseID(IdToRemove, vtkid);
3091           if (const SMDS_VtkVolume* vtkElem = dynamic_cast<const SMDS_VtkVolume*>(*it))
3092             myVolumePool->destroy((SMDS_VtkVolume*) vtkElem);
3093           else
3094             delete (*it);
3095           break;
3096       }
3097       if (vtkid >= 0)
3098         {
3099           //MESSAGE("VTK_EMPTY_CELL in " << vtkid);
3100           this->myGrid->GetCellTypesArray()->SetValue(vtkid, VTK_EMPTY_CELL);
3101         }
3102       it++;
3103     }
3104
3105   // remove exclusive (free) nodes
3106   if (removenodes)
3107     {
3108       it = s2->begin();
3109       while (it != s2->end())
3110         {
3111           int IdToRemove = (*it)->GetID();
3112           //MESSAGE( "SMDS: RM node " << IdToRemove);
3113           if (IdToRemove >= 0)
3114             {
3115               myNodes[IdToRemove] = 0;
3116               myInfo.myNbNodes--;
3117             }
3118           myNodeIDFactory->ReleaseID((*it)->GetID(), (*it)->getVtkId());
3119           removedNodes.push_back((*it));
3120           if (const SMDS_MeshNode* vtkElem = dynamic_cast<const SMDS_MeshNode*>(*it))
3121             myNodePool->destroy((SMDS_MeshNode*) vtkElem);
3122           else
3123             delete (*it);
3124           it++;
3125         }
3126     }
3127
3128   delete s2;
3129   delete s1;
3130 }
3131
3132
3133 ///////////////////////////////////////////////////////////////////////////////
3134 ///@param elem The element to delete
3135 ///////////////////////////////////////////////////////////////////////////////
3136 void SMDS_Mesh::RemoveFreeElement(const SMDS_MeshElement * elem)
3137 {
3138   int elemId = elem->GetID();
3139   int vtkId = elem->getVtkId();
3140   //MESSAGE("RemoveFreeElement " << elemId);
3141   SMDSAbs_ElementType aType = elem->GetType();
3142   SMDS_MeshElement* todest = (SMDS_MeshElement*)(elem);
3143   if (aType == SMDSAbs_Node) {
3144     //MESSAGE("Remove free node " << elemId);
3145     // only free node can be removed by this method
3146     const SMDS_MeshNode* n = static_cast<SMDS_MeshNode*>(todest);
3147     SMDS_ElemIteratorPtr itFe = n->GetInverseElementIterator();
3148     if (!itFe->more()) { // free node
3149       myNodes[elemId] = 0;
3150       myInfo.myNbNodes--;
3151       myNodePool->destroy(static_cast<SMDS_MeshNode*>(todest));
3152       myNodeIDFactory->ReleaseID(elemId, vtkId);
3153     }
3154   } else {
3155     if (hasConstructionEdges() || hasConstructionFaces())
3156       // this methods is only for meshes without descendants
3157       return;
3158
3159     //MESSAGE("Remove free element " << elemId);
3160     // Remove element from <InverseElements> of its nodes
3161     SMDS_ElemIteratorPtr itn = elem->nodesIterator();
3162     while (itn->more()) {
3163       SMDS_MeshNode * n = static_cast<SMDS_MeshNode *>
3164         (const_cast<SMDS_MeshElement *>(itn->next()));
3165       n->RemoveInverseElement(elem);
3166     }
3167
3168     // in meshes without descendants elements are always free
3169      switch (aType) {
3170     case SMDSAbs_0DElement:
3171       myCells[elemId] = 0;
3172       myInfo.remove(elem);
3173       delete elem;
3174       break;
3175     case SMDSAbs_Edge:
3176       myCells[elemId] = 0;
3177       myInfo.RemoveEdge(elem);
3178       myEdgePool->destroy(static_cast<SMDS_VtkEdge*>(todest));
3179       break;
3180     case SMDSAbs_Face:
3181       myCells[elemId] = 0;
3182       myInfo.RemoveFace(elem);
3183       myFacePool->destroy(static_cast<SMDS_VtkFace*>(todest));
3184       break;
3185     case SMDSAbs_Volume:
3186       myCells[elemId] = 0;
3187       myInfo.RemoveVolume(elem);
3188       myVolumePool->destroy(static_cast<SMDS_VtkVolume*>(todest));
3189       break;
3190     default:
3191       break;
3192     }
3193     myElementIDFactory->ReleaseID(elemId, vtkId);
3194
3195     this->myGrid->GetCellTypesArray()->SetValue(vtkId, VTK_EMPTY_CELL);
3196     // --- to do: keep vtkid in a list of reusable cells
3197   }
3198 }
3199
3200 /*!
3201  * Checks if the element is present in mesh.
3202  * Useful to determine dead pointers.
3203  */
3204 bool SMDS_Mesh::Contains (const SMDS_MeshElement* elem) const
3205 {
3206   // we should not imply on validity of *elem, so iterate on containers
3207   // of all types in the hope of finding <elem> somewhere there
3208   SMDS_NodeIteratorPtr itn = nodesIterator();
3209   while (itn->more())
3210     if (elem == itn->next())
3211       return true;
3212   SMDS_0DElementIteratorPtr it0d = elements0dIterator();
3213   while (it0d->more())
3214     if (elem == it0d->next())
3215       return true;
3216   SMDS_EdgeIteratorPtr ite = edgesIterator();
3217   while (ite->more())
3218     if (elem == ite->next())
3219       return true;
3220   SMDS_FaceIteratorPtr itf = facesIterator();
3221   while (itf->more())
3222     if (elem == itf->next())
3223       return true;
3224   SMDS_VolumeIteratorPtr itv = volumesIterator();
3225   while (itv->more())
3226     if (elem == itv->next())
3227       return true;
3228   return false;
3229 }
3230
3231 //=======================================================================
3232 //function : MaxNodeID
3233 //purpose  :
3234 //=======================================================================
3235
3236 int SMDS_Mesh::MaxNodeID() const
3237 {
3238   return myNodeMax;
3239 }
3240
3241 //=======================================================================
3242 //function : MinNodeID
3243 //purpose  :
3244 //=======================================================================
3245
3246 int SMDS_Mesh::MinNodeID() const
3247 {
3248   return myNodeMin;
3249 }
3250
3251 //=======================================================================
3252 //function : MaxElementID
3253 //purpose  :
3254 //=======================================================================
3255
3256 int SMDS_Mesh::MaxElementID() const
3257 {
3258   return myElementIDFactory->GetMaxID();
3259 }
3260
3261 //=======================================================================
3262 //function : MinElementID
3263 //purpose  :
3264 //=======================================================================
3265
3266 int SMDS_Mesh::MinElementID() const
3267 {
3268   return myElementIDFactory->GetMinID();
3269 }
3270
3271 //=======================================================================
3272 //function : Renumber
3273 //purpose  : Renumber all nodes or elements.
3274 //=======================================================================
3275
3276 void SMDS_Mesh::Renumber (const bool isNodes, const int  startID, const int  deltaID)
3277 {
3278     MESSAGE("Renumber");
3279   if ( deltaID == 0 )
3280     return;
3281
3282   SMDS_MeshNodeIDFactory * idFactory =
3283     isNodes ? myNodeIDFactory : myElementIDFactory;
3284
3285   // get existing elements in the order of ID increasing
3286   map<int,SMDS_MeshElement*> elemMap;
3287   SMDS_ElemIteratorPtr idElemIt = idFactory->elementsIterator();
3288   while ( idElemIt->more() ) {
3289     SMDS_MeshElement* elem = const_cast<SMDS_MeshElement*>(idElemIt->next());
3290     int id = elem->GetID();
3291     elemMap.insert(map<int,SMDS_MeshElement*>::value_type(id, elem));
3292   }
3293   // release their ids
3294   map<int,SMDS_MeshElement*>::iterator elemIt = elemMap.begin();
3295   idFactory->Clear();
3296 //   for ( ; elemIt != elemMap.end(); elemIt++ )
3297 //   {
3298 //     int id = (*elemIt).first;
3299 //     idFactory->ReleaseID( id );
3300 //   }
3301   // set new IDs
3302   int ID = startID;
3303   elemIt = elemMap.begin();
3304   for ( ; elemIt != elemMap.end(); elemIt++ )
3305   {
3306     idFactory->BindID( ID, (*elemIt).second );
3307     ID += deltaID;
3308   }
3309 }
3310
3311 //=======================================================================
3312 //function : GetElementType
3313 //purpose  : Return type of element or node with id
3314 //=======================================================================
3315
3316 SMDSAbs_ElementType SMDS_Mesh::GetElementType( const int id, const bool iselem ) const
3317 {
3318   SMDS_MeshElement* elem = 0;
3319   if( iselem )
3320     elem = myElementIDFactory->MeshElement( id );
3321   else
3322     elem = myNodeIDFactory->MeshElement( id );
3323
3324   if( !elem )
3325   {
3326     //throw SALOME_Exception(LOCALIZED ("this element isn't exist"));
3327     return SMDSAbs_All;
3328   }
3329   else
3330     return elem->GetType();
3331 }
3332
3333
3334
3335 //********************************************************************
3336 //********************************************************************
3337 //********                                                   *********
3338 //*****       Methods for addition of quadratic elements        ******
3339 //********                                                   *********
3340 //********************************************************************
3341 //********************************************************************
3342
3343 //=======================================================================
3344 //function : AddEdgeWithID
3345 //purpose  :
3346 //=======================================================================
3347 SMDS_MeshEdge* SMDS_Mesh::AddEdgeWithID(int n1, int n2, int n12, int ID)
3348 {
3349   return SMDS_Mesh::AddEdgeWithID
3350     ((SMDS_MeshNode*) myNodeIDFactory->MeshElement(n1),
3351      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n2),
3352      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n12),
3353      ID);
3354 }
3355
3356 //=======================================================================
3357 //function : AddEdge
3358 //purpose  :
3359 //=======================================================================
3360 SMDS_MeshEdge* SMDS_Mesh::AddEdge(const SMDS_MeshNode* n1,
3361                                   const SMDS_MeshNode* n2,
3362                                   const SMDS_MeshNode* n12)
3363 {
3364   return SMDS_Mesh::AddEdgeWithID(n1, n2, n12, myElementIDFactory->GetFreeID());
3365 }
3366
3367 //=======================================================================
3368 //function : AddEdgeWithID
3369 //purpose  :
3370 //=======================================================================
3371 SMDS_MeshEdge* SMDS_Mesh::AddEdgeWithID(const SMDS_MeshNode * n1,
3372                                         const SMDS_MeshNode * n2,
3373                                         const SMDS_MeshNode * n12,
3374                                         int ID)
3375 {
3376   if ( !n1 || !n2 || !n12 ) return 0;
3377
3378   // --- retrieve nodes ID
3379   vector<vtkIdType> nodeIds;
3380   nodeIds.clear();
3381   nodeIds.push_back(n1->getVtkId());
3382   nodeIds.push_back(n2->getVtkId());
3383   nodeIds.push_back(n12->getVtkId());
3384
3385   SMDS_MeshEdge * edge = 0;
3386   SMDS_VtkEdge *edgevtk = myEdgePool->getNew();
3387   edgevtk->init(nodeIds, this);
3388   if (!this->registerElement(ID,edgevtk))
3389     {
3390       this->myGrid->GetCellTypesArray()->SetValue(edgevtk->getVtkId(), VTK_EMPTY_CELL);
3391       myEdgePool->destroy(edgevtk);
3392       return 0;
3393     }
3394   edge = edgevtk;
3395   adjustmyCellsCapacity(ID);
3396   myCells[ID] = edge;
3397   myInfo.myNbQuadEdges++;
3398
3399 //  if (!registerElement(ID, edge)) {
3400 //        RemoveElement(edge, false);
3401 //        edge = NULL;
3402 //  }
3403   return edge;
3404
3405 }
3406
3407
3408 //=======================================================================
3409 //function : AddFace
3410 //purpose  :
3411 //=======================================================================
3412 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshNode * n1,
3413                                   const SMDS_MeshNode * n2,
3414                                   const SMDS_MeshNode * n3,
3415                                   const SMDS_MeshNode * n12,
3416                                   const SMDS_MeshNode * n23,
3417                                   const SMDS_MeshNode * n31)
3418 {
3419   return SMDS_Mesh::AddFaceWithID(n1,n2,n3,n12,n23,n31,
3420                                   myElementIDFactory->GetFreeID());
3421 }
3422
3423 //=======================================================================
3424 //function : AddFaceWithID
3425 //purpose  :
3426 //=======================================================================
3427 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(int n1, int n2, int n3,
3428                                         int n12,int n23,int n31, int ID)
3429 {
3430   return SMDS_Mesh::AddFaceWithID
3431     ((SMDS_MeshNode *)myNodeIDFactory->MeshElement(n1) ,
3432      (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n2) ,
3433      (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n3) ,
3434      (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n12),
3435      (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n23),
3436      (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n31),
3437      ID);
3438 }
3439
3440 //=======================================================================
3441 //function : AddFaceWithID
3442 //purpose  :
3443 //=======================================================================
3444 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshNode * n1,
3445                                         const SMDS_MeshNode * n2,
3446                                         const SMDS_MeshNode * n3,
3447                                         const SMDS_MeshNode * n12,
3448                                         const SMDS_MeshNode * n23,
3449                                         const SMDS_MeshNode * n31,
3450                                         int ID)
3451 {
3452   if ( !n1 || !n2 || !n3 || !n12 || !n23 || !n31) return 0;
3453   if(hasConstructionEdges()) {
3454     // creation quadratic edges - not implemented
3455     return 0;
3456   }
3457   else
3458   {
3459     // --- retrieve nodes ID
3460     vector<vtkIdType> nodeIds;
3461     nodeIds.clear();
3462     nodeIds.push_back(n1->getVtkId());
3463     nodeIds.push_back(n2->getVtkId());
3464     nodeIds.push_back(n3->getVtkId());
3465     nodeIds.push_back(n12->getVtkId());
3466     nodeIds.push_back(n23->getVtkId());
3467     nodeIds.push_back(n31->getVtkId());
3468
3469     SMDS_MeshFace * face = 0;
3470     SMDS_VtkFace *facevtk = myFacePool->getNew();
3471     facevtk->init(nodeIds, this);
3472     if (!this->registerElement(ID,facevtk))
3473       {
3474         this->myGrid->GetCellTypesArray()->SetValue(facevtk->getVtkId(), VTK_EMPTY_CELL);
3475         myFacePool->destroy(facevtk);
3476         return 0;
3477       }
3478     face = facevtk;
3479     adjustmyCellsCapacity(ID);
3480     myCells[ID] = face;
3481     myInfo.myNbQuadTriangles++;
3482
3483 //    if (!registerElement(ID, face)) {
3484 //      RemoveElement(face, false);
3485 //      face = NULL;
3486 //    }
3487     return face;
3488   }
3489 }
3490
3491
3492 //=======================================================================
3493 //function : AddFace
3494 //purpose  :
3495 //=======================================================================
3496 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshNode * n1,
3497                                   const SMDS_MeshNode * n2,
3498                                   const SMDS_MeshNode * n3,
3499                                   const SMDS_MeshNode * n4,
3500                                   const SMDS_MeshNode * n12,
3501                                   const SMDS_MeshNode * n23,
3502                                   const SMDS_MeshNode * n34,
3503                                   const SMDS_MeshNode * n41)
3504 {
3505   return SMDS_Mesh::AddFaceWithID(n1,n2,n3,n4,n12,n23,n34,n41,
3506                                   myElementIDFactory->GetFreeID());
3507 }
3508
3509 //=======================================================================
3510 //function : AddFaceWithID
3511 //purpose  :
3512 //=======================================================================
3513 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(int n1, int n2, int n3, int n4,
3514                                         int n12,int n23,int n34,int n41, int ID)
3515 {
3516   return SMDS_Mesh::AddFaceWithID
3517     ((SMDS_MeshNode *)myNodeIDFactory->MeshElement(n1) ,
3518      (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n2) ,
3519      (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n3) ,
3520      (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n4) ,
3521      (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n12),
3522      (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n23),
3523      (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n34),
3524      (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n41),
3525      ID);
3526 }
3527
3528 //=======================================================================
3529 //function : AddFaceWithID
3530 //purpose  :
3531 //=======================================================================
3532 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshNode * n1,
3533                                         const SMDS_MeshNode * n2,
3534                                         const SMDS_MeshNode * n3,
3535                                         const SMDS_MeshNode * n4,
3536                                         const SMDS_MeshNode * n12,
3537                                         const SMDS_MeshNode * n23,
3538                                         const SMDS_MeshNode * n34,
3539                                         const SMDS_MeshNode * n41,
3540                                         int ID)
3541 {
3542   if ( !n1 || !n2 || !n3 || !n4 || !n12 || !n23 || !n34 || !n41) return 0;
3543   if(hasConstructionEdges()) {
3544     // creation quadratic edges - not implemented
3545         return 0;
3546   }
3547   else
3548   {
3549     // --- retrieve nodes ID
3550     vector<vtkIdType> nodeIds;
3551     nodeIds.clear();
3552     nodeIds.push_back(n1->getVtkId());
3553     nodeIds.push_back(n2->getVtkId());
3554     nodeIds.push_back(n3->getVtkId());
3555     nodeIds.push_back(n4->getVtkId());
3556     nodeIds.push_back(n12->getVtkId());
3557     nodeIds.push_back(n23->getVtkId());
3558     nodeIds.push_back(n34->getVtkId());
3559     nodeIds.push_back(n41->getVtkId());
3560
3561     SMDS_MeshFace * face = 0;
3562     SMDS_VtkFace *facevtk = myFacePool->getNew();
3563     facevtk->init(nodeIds, this);
3564     if (!this->registerElement(ID,facevtk))
3565       {
3566         this->myGrid->GetCellTypesArray()->SetValue(facevtk->getVtkId(), VTK_EMPTY_CELL);
3567         myFacePool->destroy(facevtk);
3568         return 0;
3569       }
3570     face = facevtk;
3571     adjustmyCellsCapacity(ID);
3572     myCells[ID] = face;
3573     myInfo.myNbQuadQuadrangles++;
3574
3575 //    if (!registerElement(ID, face)) {
3576 //      RemoveElement(face, false);
3577 //      face = NULL;
3578 //    }
3579     return face;
3580   }
3581 }
3582
3583
3584 //=======================================================================
3585 //function : AddVolume
3586 //purpose  :
3587 //=======================================================================
3588 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
3589                                       const SMDS_MeshNode * n2,
3590                                       const SMDS_MeshNode * n3,
3591                                       const SMDS_MeshNode * n4,
3592                                       const SMDS_MeshNode * n12,
3593                                       const SMDS_MeshNode * n23,
3594                                       const SMDS_MeshNode * n31,
3595                                       const SMDS_MeshNode * n14,
3596                                       const SMDS_MeshNode * n24,
3597                                       const SMDS_MeshNode * n34)
3598 {
3599   int ID = myElementIDFactory->GetFreeID();
3600   SMDS_MeshVolume * v = SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n12, n23,
3601                                                    n31, n14, n24, n34, ID);
3602   if(v==NULL) myElementIDFactory->ReleaseID(ID);
3603   return v;
3604 }
3605
3606 //=======================================================================
3607 //function : AddVolumeWithID
3608 //purpose  :
3609 //=======================================================================
3610 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3, int n4,
3611                                             int n12,int n23,int n31,
3612                                             int n14,int n24,int n34, int ID)
3613 {
3614   return SMDS_Mesh::AddVolumeWithID
3615     ((SMDS_MeshNode*) myNodeIDFactory->MeshElement(n1) ,
3616      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n2) ,
3617      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n3) ,
3618      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n4) ,
3619      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n12),
3620      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n23),
3621      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n31),
3622      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n14),
3623      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n24),
3624      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n34),
3625      ID);
3626 }
3627
3628 //=======================================================================
3629 //function : AddVolumeWithID
3630 //purpose  : 2d order tetrahedron of 10 nodes
3631 //=======================================================================
3632 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
3633                                             const SMDS_MeshNode * n2,
3634                                             const SMDS_MeshNode * n3,
3635                                             const SMDS_MeshNode * n4,
3636                                             const SMDS_MeshNode * n12,
3637                                             const SMDS_MeshNode * n23,
3638                                             const SMDS_MeshNode * n31,
3639                                             const SMDS_MeshNode * n14,
3640                                             const SMDS_MeshNode * n24,
3641                                             const SMDS_MeshNode * n34,
3642                                             int ID)
3643 {
3644   if ( !n1 || !n2 || !n3 || !n4 || !n12 || !n23 || !n31 || !n14 || !n24 || !n34)
3645     return 0;
3646   if(hasConstructionFaces()) {
3647     // creation quadratic faces - not implemented
3648     return 0;
3649   }
3650   // --- retrieve nodes ID
3651   vector<vtkIdType> nodeIds;
3652   nodeIds.clear();
3653   nodeIds.push_back(n1->getVtkId());
3654   nodeIds.push_back(n3->getVtkId());
3655   nodeIds.push_back(n2->getVtkId());
3656   nodeIds.push_back(n4->getVtkId());
3657
3658   nodeIds.push_back(n31->getVtkId());
3659   nodeIds.push_back(n23->getVtkId());
3660   nodeIds.push_back(n12->getVtkId());
3661
3662   nodeIds.push_back(n14->getVtkId());
3663   nodeIds.push_back(n34->getVtkId());
3664   nodeIds.push_back(n24->getVtkId());
3665
3666   SMDS_VtkVolume *volvtk = myVolumePool->getNew();
3667   volvtk->init(nodeIds, this);
3668   if (!this->registerElement(ID,volvtk))
3669     {
3670       this->myGrid->GetCellTypesArray()->SetValue(volvtk->getVtkId(), VTK_EMPTY_CELL);
3671       myVolumePool->destroy(volvtk);
3672       return 0;
3673     }
3674   adjustmyCellsCapacity(ID);
3675   myCells[ID] = volvtk;
3676   myInfo.myNbQuadTetras++;
3677
3678 //  if (!registerElement(ID, volvtk)) {
3679 //    RemoveElement(volvtk, false);
3680 //    volvtk = NULL;
3681 //  }
3682   return volvtk;
3683 }
3684
3685
3686 //=======================================================================
3687 //function : AddVolume
3688 //purpose  :
3689 //=======================================================================
3690 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
3691                                       const SMDS_MeshNode * n2,
3692                                       const SMDS_MeshNode * n3,
3693                                       const SMDS_MeshNode * n4,
3694                                       const SMDS_MeshNode * n5,
3695                                       const SMDS_MeshNode * n12,
3696                                       const SMDS_MeshNode * n23,
3697                                       const SMDS_MeshNode * n34,
3698                                       const SMDS_MeshNode * n41,
3699                                       const SMDS_MeshNode * n15,
3700                                       const SMDS_MeshNode * n25,
3701                                       const SMDS_MeshNode * n35,
3702                                       const SMDS_MeshNode * n45)
3703 {
3704   int ID = myElementIDFactory->GetFreeID();
3705   SMDS_MeshVolume * v =
3706     SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n12, n23, n34, n41,
3707                                n15, n25, n35, n45, ID);
3708   if(v==NULL) myElementIDFactory->ReleaseID(ID);
3709   return v;
3710 }
3711
3712 //=======================================================================
3713 //function : AddVolumeWithID
3714 //purpose  :
3715 //=======================================================================
3716 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3, int n4, int n5,
3717                                             int n12,int n23,int n34,int n41,
3718                                             int n15,int n25,int n35,int n45, int ID)
3719 {
3720   return SMDS_Mesh::AddVolumeWithID
3721     ((SMDS_MeshNode*) myNodeIDFactory->MeshElement(n1) ,
3722      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n2) ,
3723      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n3) ,
3724      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n4) ,
3725      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n5) ,
3726      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n12),
3727      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n23),
3728      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n34),
3729      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n41),
3730      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n15),
3731      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n25),
3732      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n35),
3733      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n45),
3734      ID);
3735 }
3736
3737 //=======================================================================
3738 //function : AddVolumeWithID
3739 //purpose  : 2d order pyramid of 13 nodes
3740 //=======================================================================
3741 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
3742                                             const SMDS_MeshNode * n2,
3743                                             const SMDS_MeshNode * n3,
3744                                             const SMDS_MeshNode * n4,
3745                                             const SMDS_MeshNode * n5,
3746                                             const SMDS_MeshNode * n12,
3747                                             const SMDS_MeshNode * n23,
3748                                             const SMDS_MeshNode * n34,
3749                                             const SMDS_MeshNode * n41,
3750                                             const SMDS_MeshNode * n15,
3751                                             const SMDS_MeshNode * n25,
3752                                             const SMDS_MeshNode * n35,
3753                                             const SMDS_MeshNode * n45,
3754                                             int ID)
3755 {
3756   if (!n1 || !n2 || !n3 || !n4 || !n5 || !n12 || !n23 ||
3757       !n34 || !n41 || !n15 || !n25 || !n35 || !n45)
3758     return 0;
3759   if(hasConstructionFaces()) {
3760     // creation quadratic faces - not implemented
3761     return 0;
3762   }
3763   // --- retrieve nodes ID
3764   vector<vtkIdType> nodeIds;
3765   nodeIds.clear();
3766   nodeIds.push_back(n1->getVtkId());
3767   nodeIds.push_back(n4->getVtkId());
3768   nodeIds.push_back(n3->getVtkId());
3769   nodeIds.push_back(n2->getVtkId());
3770   nodeIds.push_back(n5->getVtkId());
3771
3772   nodeIds.push_back(n41->getVtkId());
3773   nodeIds.push_back(n34->getVtkId());
3774   nodeIds.push_back(n23->getVtkId());
3775   nodeIds.push_back(n12->getVtkId());
3776
3777   nodeIds.push_back(n15->getVtkId());
3778   nodeIds.push_back(n45->getVtkId());
3779   nodeIds.push_back(n35->getVtkId());
3780   nodeIds.push_back(n25->getVtkId());
3781
3782   SMDS_VtkVolume *volvtk = myVolumePool->getNew();
3783   volvtk->init(nodeIds, this);
3784   if (!this->registerElement(ID,volvtk))
3785     {
3786       this->myGrid->GetCellTypesArray()->SetValue(volvtk->getVtkId(), VTK_EMPTY_CELL);
3787       myVolumePool->destroy(volvtk);
3788       return 0;
3789     }
3790   adjustmyCellsCapacity(ID);
3791   myCells[ID] = volvtk;
3792   myInfo.myNbQuadPyramids++;
3793
3794 //  if (!registerElement(ID, volvtk)) {
3795 //    RemoveElement(volvtk, false);
3796 //    volvtk = NULL;
3797 //  }
3798   return volvtk;
3799 }
3800
3801
3802 //=======================================================================
3803 //function : AddVolume
3804 //purpose  :
3805 //=======================================================================
3806 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
3807                                       const SMDS_MeshNode * n2,
3808                                       const SMDS_MeshNode * n3,
3809                                       const SMDS_MeshNode * n4,
3810                                       const SMDS_MeshNode * n5,
3811                                       const SMDS_MeshNode * n6,
3812                                       const SMDS_MeshNode * n12,
3813                                       const SMDS_MeshNode * n23,
3814                                       const SMDS_MeshNode * n31,
3815                                       const SMDS_MeshNode * n45,
3816                                       const SMDS_MeshNode * n56,
3817                                       const SMDS_MeshNode * n64,
3818                                       const SMDS_MeshNode * n14,
3819                                       const SMDS_MeshNode * n25,
3820                                       const SMDS_MeshNode * n36)
3821 {
3822   int ID = myElementIDFactory->GetFreeID();
3823   SMDS_MeshVolume * v =
3824     SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n6, n12, n23, n31,
3825                                n45, n56, n64, n14, n25, n36, ID);
3826   if(v==NULL) myElementIDFactory->ReleaseID(ID);
3827   return v;
3828 }
3829
3830 //=======================================================================
3831 //function : AddVolumeWithID
3832 //purpose  :
3833 //=======================================================================
3834 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3,
3835                                             int n4, int n5, int n6,
3836                                             int n12,int n23,int n31,
3837                                             int n45,int n56,int n64,
3838                                             int n14,int n25,int n36, int ID)
3839 {
3840   return SMDS_Mesh::AddVolumeWithID
3841     ((SMDS_MeshNode*) myNodeIDFactory->MeshElement(n1) ,
3842      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n2) ,
3843      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n3) ,
3844      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n4) ,
3845      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n5) ,
3846      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n6) ,
3847      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n12),
3848      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n23),
3849      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n31),
3850      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n45),
3851      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n56),
3852      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n64),
3853      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n14),
3854      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n25),
3855      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n36),
3856      ID);
3857 }
3858
3859 //=======================================================================
3860 //function : AddVolumeWithID
3861 //purpose  : 2d order Pentahedron with 15 nodes
3862 //=======================================================================
3863 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
3864                                             const SMDS_MeshNode * n2,
3865                                             const SMDS_MeshNode * n3,
3866                                             const SMDS_MeshNode * n4,
3867                                             const SMDS_MeshNode * n5,
3868                                             const SMDS_MeshNode * n6,
3869                                             const SMDS_MeshNode * n12,
3870                                             const SMDS_MeshNode * n23,
3871                                             const SMDS_MeshNode * n31,
3872                                             const SMDS_MeshNode * n45,
3873                                             const SMDS_MeshNode * n56,
3874                                             const SMDS_MeshNode * n64,
3875                                             const SMDS_MeshNode * n14,
3876                                             const SMDS_MeshNode * n25,
3877                                             const SMDS_MeshNode * n36,
3878                                             int ID)
3879 {
3880   if (!n1 || !n2 || !n3 || !n4 || !n5 || !n6 || !n12 || !n23 ||
3881       !n31 || !n45 || !n56 || !n64 || !n14 || !n25 || !n36)
3882     return 0;
3883   if(hasConstructionFaces()) {
3884     // creation quadratic faces - not implemented
3885     return 0;
3886   }
3887   // --- retrieve nodes ID
3888   vector<vtkIdType> nodeIds;
3889   nodeIds.clear();
3890   nodeIds.push_back(n1->getVtkId());
3891   nodeIds.push_back(n2->getVtkId());
3892   nodeIds.push_back(n3->getVtkId());
3893
3894   nodeIds.push_back(n4->getVtkId());
3895   nodeIds.push_back(n5->getVtkId());
3896   nodeIds.push_back(n6->getVtkId());
3897
3898   nodeIds.push_back(n12->getVtkId());
3899   nodeIds.push_back(n23->getVtkId());
3900   nodeIds.push_back(n31->getVtkId());
3901
3902   nodeIds.push_back(n45->getVtkId());
3903   nodeIds.push_back(n56->getVtkId());
3904   nodeIds.push_back(n64->getVtkId());
3905
3906   nodeIds.push_back(n14->getVtkId());
3907   nodeIds.push_back(n25->getVtkId());
3908   nodeIds.push_back(n36->getVtkId());
3909
3910   SMDS_VtkVolume *volvtk = myVolumePool->getNew();
3911   volvtk->init(nodeIds, this);
3912   if (!this->registerElement(ID,volvtk))
3913     {
3914       this->myGrid->GetCellTypesArray()->SetValue(volvtk->getVtkId(), VTK_EMPTY_CELL);
3915       myVolumePool->destroy(volvtk);
3916       return 0;
3917     }
3918   adjustmyCellsCapacity(ID);
3919   myCells[ID] = volvtk;
3920   myInfo.myNbQuadPrisms++;
3921
3922 //  if (!registerElement(ID, volvtk)) {
3923 //    RemoveElement(volvtk, false);
3924 //    volvtk = NULL;
3925 //  }
3926   return volvtk;
3927 }
3928
3929
3930 //=======================================================================
3931 //function : AddVolume
3932 //purpose  :
3933 //=======================================================================
3934 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
3935                                       const SMDS_MeshNode * n2,
3936                                       const SMDS_MeshNode * n3,
3937                                       const SMDS_MeshNode * n4,
3938                                       const SMDS_MeshNode * n5,
3939                                       const SMDS_MeshNode * n6,
3940                                       const SMDS_MeshNode * n7,
3941                                       const SMDS_MeshNode * n8,
3942                                       const SMDS_MeshNode * n12,
3943                                       const SMDS_MeshNode * n23,
3944                                       const SMDS_MeshNode * n34,
3945                                       const SMDS_MeshNode * n41,
3946                                       const SMDS_MeshNode * n56,
3947                                       const SMDS_MeshNode * n67,
3948                                       const SMDS_MeshNode * n78,
3949                                       const SMDS_MeshNode * n85,
3950                                       const SMDS_MeshNode * n15,
3951                                       const SMDS_MeshNode * n26,
3952                                       const SMDS_MeshNode * n37,
3953                                       const SMDS_MeshNode * n48)
3954 {
3955   int ID = myElementIDFactory->GetFreeID();
3956   SMDS_MeshVolume * v =
3957     SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n6, n7, n8, n12, n23, n34, n41,
3958                                n56, n67, n78, n85, n15, n26, n37, n48, ID);
3959   if(v==NULL) myElementIDFactory->ReleaseID(ID);
3960   return v;
3961 }
3962
3963 //=======================================================================
3964 //function : AddVolumeWithID
3965 //purpose  :
3966 //=======================================================================
3967 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3, int n4,
3968                                             int n5, int n6, int n7, int n8,
3969                                             int n12,int n23,int n34,int n41,
3970                                             int n56,int n67,int n78,int n85,
3971                                             int n15,int n26,int n37,int n48, int ID)
3972 {
3973   return SMDS_Mesh::AddVolumeWithID
3974     ((SMDS_MeshNode*) myNodeIDFactory->MeshElement(n1),
3975      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n2),
3976      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n3),
3977      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n4),
3978      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n5),
3979      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n6),
3980      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n7),
3981      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n8),
3982      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n12),
3983      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n23),
3984      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n34),
3985      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n41),
3986      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n56),
3987      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n67),
3988      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n78),
3989      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n85),
3990      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n15),
3991      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n26),
3992      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n37),
3993      (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n48),
3994      ID);
3995 }
3996
3997 //=======================================================================
3998 //function : AddVolumeWithID
3999 //purpose  : 2d order Hexahedrons with 20 nodes
4000 //=======================================================================
4001 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
4002                                             const SMDS_MeshNode * n2,
4003                                             const SMDS_MeshNode * n3,
4004                                             const SMDS_MeshNode * n4,
4005                                             const SMDS_MeshNode * n5,
4006                                             const SMDS_MeshNode * n6,
4007                                             const SMDS_MeshNode * n7,
4008                                             const SMDS_MeshNode * n8,
4009                                             const SMDS_MeshNode * n12,
4010                                             const SMDS_MeshNode * n23,
4011                                             const SMDS_MeshNode * n34,
4012                                             const SMDS_MeshNode * n41,
4013                                             const SMDS_MeshNode * n56,
4014                                             const SMDS_MeshNode * n67,
4015                                             const SMDS_MeshNode * n78,
4016                                             const SMDS_MeshNode * n85,
4017                                             const SMDS_MeshNode * n15,
4018                                             const SMDS_MeshNode * n26,
4019                                             const SMDS_MeshNode * n37,
4020                                             const SMDS_MeshNode * n48,
4021                                             int ID)
4022 {
4023   if (!n1 || !n2 || !n3 || !n4 || !n5 || !n6 || !n7 || !n8 || !n12 || !n23 ||
4024       !n34 || !n41 || !n56 || !n67 || !n78 || !n85 || !n15 || !n26 || !n37 || !n48)
4025     return 0;
4026   if(hasConstructionFaces()) {
4027     return 0;
4028     // creation quadratic faces - not implemented
4029   }
4030   // --- retrieve nodes ID
4031   vector<vtkIdType> nodeIds;
4032   nodeIds.clear();
4033   nodeIds.push_back(n1->getVtkId());
4034   nodeIds.push_back(n4->getVtkId());
4035   nodeIds.push_back(n3->getVtkId());
4036   nodeIds.push_back(n2->getVtkId());
4037
4038   nodeIds.push_back(n5->getVtkId());
4039   nodeIds.push_back(n8->getVtkId());
4040   nodeIds.push_back(n7->getVtkId());
4041   nodeIds.push_back(n6->getVtkId());
4042
4043   nodeIds.push_back(n41->getVtkId());
4044   nodeIds.push_back(n34->getVtkId());
4045   nodeIds.push_back(n23->getVtkId());
4046   nodeIds.push_back(n12->getVtkId());
4047
4048   nodeIds.push_back(n85->getVtkId());
4049   nodeIds.push_back(n78->getVtkId());
4050   nodeIds.push_back(n67->getVtkId());
4051   nodeIds.push_back(n56->getVtkId());
4052
4053   nodeIds.push_back(n15->getVtkId());
4054   nodeIds.push_back(n48->getVtkId());
4055   nodeIds.push_back(n37->getVtkId());
4056   nodeIds.push_back(n26->getVtkId());
4057
4058   SMDS_VtkVolume *volvtk = myVolumePool->getNew();
4059   volvtk->init(nodeIds, this);
4060   if (!this->registerElement(ID,volvtk))
4061     {
4062       this->myGrid->GetCellTypesArray()->SetValue(volvtk->getVtkId(), VTK_EMPTY_CELL);
4063       myVolumePool->destroy(volvtk);
4064       return 0;
4065     }
4066   adjustmyCellsCapacity(ID);
4067   myCells[ID] = volvtk;
4068   myInfo.myNbQuadHexas++;
4069
4070 //  if (!registerElement(ID, volvtk)) {
4071 //    RemoveElement(volvtk, false);
4072 //    volvtk = NULL;
4073 //  }
4074   return volvtk;
4075 }
4076
4077 void SMDS_Mesh::updateNodeMinMax()
4078 {
4079   myNodeMin = 0;
4080   if (myNodes.size() == 0)
4081   {
4082         myNodeMax=0;
4083         return;
4084   }
4085   while (!myNodes[myNodeMin] && (myNodeMin<myNodes.size()))
4086     myNodeMin++;
4087   myNodeMax=myNodes.size()-1;
4088   while (!myNodes[myNodeMax] && (myNodeMin>=0))
4089     myNodeMin--;
4090 }
4091
4092 void SMDS_Mesh::incrementNodesCapacity(int nbNodes)
4093 {
4094 //  int val = myCellIdSmdsToVtk.size();
4095 //  MESSAGE(" ------------------- resize myCellIdSmdsToVtk " << val << " --> " << val + nbNodes);
4096 //  myCellIdSmdsToVtk.resize(val + nbNodes, -1); // fill new elements with -1
4097   int val = myNodes.size();
4098   MESSAGE(" ------------------- resize myNodes " << val << " --> " << val + nbNodes);
4099   myNodes.resize(val +nbNodes, 0);
4100 }
4101
4102 void SMDS_Mesh::incrementCellsCapacity(int nbCells)
4103 {
4104   int val = myCellIdVtkToSmds.size();
4105   MESSAGE(" ------------------- resize myCellIdVtkToSmds " << val << " --> " << val + nbCells);
4106   myCellIdVtkToSmds.resize(val + nbCells, -1); // fill new elements with -1
4107   val = myCells.size();
4108   MESSAGE(" ------------------- resize myCells " << val << " --> " << val + nbCells);
4109   myNodes.resize(val +nbCells, 0);
4110 }
4111
4112 void SMDS_Mesh::adjustStructure()
4113 {
4114   myGrid->GetPoints()->GetData()->SetNumberOfTuples(myNodeIDFactory->GetMaxID());
4115 }
4116
4117 void SMDS_Mesh::dumpGrid(string ficdump)
4118 {
4119         MESSAGE("SMDS_Mesh::dumpGrid " << ficdump);
4120 //  vtkUnstructuredGridWriter* aWriter = vtkUnstructuredGridWriter::New();
4121 //  aWriter->SetFileName(ficdump.c_str());
4122 //  aWriter->SetInput(myGrid);
4123 //  if(myGrid->GetNumberOfCells())
4124 //  {
4125 //    aWriter->Write();
4126 //  }
4127 //  aWriter->Delete();
4128   ficdump = ficdump + "_connectivity";
4129   ofstream ficcon(ficdump.c_str(), ios::out);
4130   int nbPoints = myGrid->GetNumberOfPoints();
4131   ficcon << "-------------------------------- points " <<  nbPoints << endl;
4132   for (int i=0; i<nbPoints; i++)
4133   {
4134         ficcon << i << " " << *(myGrid->GetPoint(i)) << " " << *(myGrid->GetPoint(i)+1) << " " << " " << *(myGrid->GetPoint(i)+2) << endl;
4135   }
4136   int nbCells = myGrid->GetNumberOfCells();
4137   ficcon << "-------------------------------- cells " <<  nbCells << endl;
4138   for (int i=0; i<nbCells; i++)
4139   {
4140 //      MESSAGE(i << " " << myGrid->GetCell(i));
4141 //      MESSAGE("  " << myGrid->GetCell(i)->GetCellType());
4142         ficcon << i << " - " << myGrid->GetCell(i)->GetCellType() << " -";
4143         int nbptcell = myGrid->GetCell(i)->GetNumberOfPoints();
4144         vtkIdList *listid = myGrid->GetCell(i)->GetPointIds();
4145         for (int j=0; j<nbptcell; j++)
4146         {
4147                 ficcon << " " <<  listid->GetId(j);
4148         }
4149         ficcon << endl;
4150   }
4151   ficcon << "-------------------------------- connectivity " <<  nbPoints << endl;
4152         vtkCellLinks *links = myGrid->GetCellLinks();
4153   for (int i=0; i<nbPoints; i++)
4154   {
4155         int ncells = links->GetNcells(i);
4156         vtkIdType *cells = links->GetCells(i);
4157         ficcon << i << " - " << ncells << " -";
4158         for (int j=0; j<ncells; j++)
4159         {
4160                 ficcon << " " << cells[j];
4161         }
4162         ficcon << endl;
4163   }
4164   ficcon.close();
4165
4166 }
4167
4168 void SMDS_Mesh::compactMesh()
4169 {
4170   MESSAGE("SMDS_Mesh::compactMesh do nothing!");
4171 }
4172
4173 int SMDS_Mesh::fromVtkToSmds(int vtkid)
4174 {
4175   if (vtkid >= 0 && vtkid < myCellIdVtkToSmds.size())
4176     return myCellIdVtkToSmds[vtkid];
4177   throw SALOME_Exception(LOCALIZED ("vtk id out of bounds"));
4178 }
4179
4180 void SMDS_Mesh::updateBoundingBox()
4181 {
4182   xmin = 0; xmax = 0;
4183   ymin = 0; ymax = 0;
4184   zmin = 0; zmax = 0;
4185   vtkPoints *points = myGrid->GetPoints();
4186   int myNodesSize = this->myNodes.size();
4187   for (int i = 0; i < myNodesSize; i++)
4188     {
4189       if (SMDS_MeshNode *n = myNodes[i])
4190         {
4191           double coords[3];
4192           points->GetPoint(n->myVtkID, coords);
4193           if (coords[0] < xmin) xmin = coords[0];
4194           else if (coords[0] > xmax) xmax = coords[0];
4195           if (coords[1] < ymin) ymin = coords[1];
4196           else if (coords[1] > ymax) ymax = coords[1];
4197           if (coords[2] < zmin) zmin = coords[2];
4198           else if (coords[2] > zmax) zmax = coords[2];
4199         }
4200     }
4201 }
4202
4203 double SMDS_Mesh::getMaxDim()
4204 {
4205   double dmax = 1.e-3;
4206   if ((xmax - xmin) > dmax) dmax = xmax -xmin;
4207   if ((ymax - ymin) > dmax) dmax = ymax -ymin;
4208   if ((zmax - zmin) > dmax) dmax = zmax -zmin;
4209   MESSAGE("getMaxDim " << dmax);
4210   return dmax;
4211 }
4212
4213 //! modification that needs compact structure and redraw
4214 void SMDS_Mesh::Modified()
4215 {
4216   if (this->myModified)
4217     {
4218       this->myModifTime++;
4219       MESSAGE("modified");
4220       myModified = false;
4221     }
4222 }
4223
4224 //! get last modification timeStamp
4225 unsigned long SMDS_Mesh::GetMTime()
4226 {
4227   return this->myModifTime;
4228 }
4229
4230 bool SMDS_Mesh::isCompacted()
4231 {
4232   if (this->myModifTime > this->myCompactTime)
4233     {
4234       MESSAGE(" *** isCompacted " << myCompactTime << " < " << myModifTime);
4235       this->myCompactTime = this->myModifTime;
4236       return false;
4237     }
4238   return true;
4239 }