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