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