1 // Copyright (C) 2007-2008 CEA/DEN, EDF R&D, OPEN CASCADE
3 // Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
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.
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.
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
20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
22 // SMESH SMDS : implementaion of Salome mesh data structure
25 #pragma warning(disable:4786)
28 #include "utilities.h"
29 #include "SMDS_Mesh.hxx"
30 #include "SMDS_VolumeOfNodes.hxx"
31 #include "SMDS_VolumeOfFaces.hxx"
32 #include "SMDS_FaceOfNodes.hxx"
33 #include "SMDS_FaceOfEdges.hxx"
34 #include "SMDS_PolyhedralVolumeOfNodes.hxx"
35 #include "SMDS_PolygonalFaceOfNodes.hxx"
36 #include "SMDS_QuadraticEdge.hxx"
37 #include "SMDS_QuadraticFaceOfNodes.hxx"
38 #include "SMDS_QuadraticVolumeOfNodes.hxx"
40 #include <vtkUnstructuredGrid.h>
47 #include <sys/sysinfo.h>
50 // number of added entitis to check memory after
51 #define CHECKMEMORY_INTERVAL 1000
53 vector<SMDS_Mesh*> SMDS_Mesh::_meshList = vector<SMDS_Mesh*>();
54 int SMDS_Mesh::chunkSize = 1000;
56 //================================================================================
58 * \brief Raise an exception if free memory (ram+swap) too low
59 * \param doNotRaise - if true, suppres exception, just return free memory size
60 * \retval int - amount of available memory in MB or negative number in failure case
62 //================================================================================
64 int SMDS_Mesh::CheckMemory(const bool doNotRaise) throw (std::bad_alloc)
68 int err = sysinfo( &si );
72 static int limit = -1;
74 int status = system("SMDS_MemoryLimit"); // it returns lower limit of free RAM
76 limit = WEXITSTATUS(status);
81 limit = int( limit * 1.5 );
83 MESSAGE ( "SMDS_Mesh::CheckMemory() memory limit = " << limit << " MB" );
87 const unsigned long Mbyte = 1024 * 1024;
88 // compute separately to avoid overflow
90 ( si.freeram * si.mem_unit ) / Mbyte +
91 ( si.freeswap * si.mem_unit ) / Mbyte;
94 return freeMb - limit;
99 MESSAGE ("SMDS_Mesh::CheckMemory() throws as free memory too low: " << freeMb <<" MB" );
101 throw std::bad_alloc();
107 ///////////////////////////////////////////////////////////////////////////////
108 /// Create a new mesh object
109 ///////////////////////////////////////////////////////////////////////////////
110 SMDS_Mesh::SMDS_Mesh()
112 myNodeIDFactory(new SMDS_MeshNodeIDFactory()),
113 myElementIDFactory(new SMDS_MeshElementIDFactory()),
114 myHasConstructionEdges(false), myHasConstructionFaces(false),
115 myHasInverseElements(true),
116 myNodeMin(0), myNodeMax(0), myCellLinksSize(0)
118 myMeshId = _meshList.size(); // --- index of the mesh to push back in the vector
119 MESSAGE("myMeshId=" << myMeshId);
120 myNodeIDFactory->SetMesh(this);
121 myElementIDFactory->SetMesh(this);
122 _meshList.push_back(this);
125 myGrid = vtkUnstructuredGrid::New();
126 myGrid->Initialize();
128 vtkPoints* points = vtkPoints::New();
129 points->SetNumberOfPoints(SMDS_Mesh::chunkSize);
130 myGrid->SetPoints( points );
132 myGrid->BuildLinks();
135 ///////////////////////////////////////////////////////////////////////////////
136 /// Create a new child mesh
137 /// Note that the tree structure of SMDS_Mesh seems to be unused in this version
138 /// (2003-09-08) of SMESH
139 ///////////////////////////////////////////////////////////////////////////////
140 SMDS_Mesh::SMDS_Mesh(SMDS_Mesh * parent)
141 :myParent(parent), myNodeIDFactory(parent->myNodeIDFactory),
142 myElementIDFactory(parent->myElementIDFactory),
143 myHasConstructionEdges(false), myHasConstructionFaces(false),
144 myHasInverseElements(true)
148 ///////////////////////////////////////////////////////////////////////////////
149 ///Create a submesh and add it to the current mesh
150 ///////////////////////////////////////////////////////////////////////////////
152 SMDS_Mesh *SMDS_Mesh::AddSubMesh()
154 SMDS_Mesh *submesh = new SMDS_Mesh(this);
155 myChildren.insert(myChildren.end(), submesh);
159 ///////////////////////////////////////////////////////////////////////////////
160 ///create a MeshNode and add it to the current Mesh
161 ///An ID is automatically assigned to the node.
162 ///@return : The created node
163 ///////////////////////////////////////////////////////////////////////////////
165 SMDS_MeshNode * SMDS_Mesh::AddNode(double x, double y, double z)
167 return SMDS_Mesh::AddNodeWithID(x,y,z,myNodeIDFactory->GetFreeID());
170 ///////////////////////////////////////////////////////////////////////////////
171 ///create a MeshNode and add it to the current Mesh
172 ///@param ID : The ID of the MeshNode to create
173 ///@return : The created node or NULL if a node with this ID already exists
174 ///////////////////////////////////////////////////////////////////////////////
175 SMDS_MeshNode * SMDS_Mesh::AddNodeWithID(double x, double y, double z, int ID)
177 // find the MeshNode corresponding to ID
178 const SMDS_MeshElement *node = myNodeIDFactory->MeshElement(ID);
180 //if ( myNodes.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
181 SMDS_MeshNode * node=new SMDS_MeshNode(ID, myMeshId, -1, x, y, z);
182 if (ID >= myNodes.size())
184 myNodes.resize(ID+SMDS_Mesh::chunkSize, 0);
185 //MESSAGE(" ------------------ myNodes resize " << ID << " --> " << ID+SMDS_Mesh::chunkSize);
188 myNodeIDFactory->BindID(ID,node);
195 ///////////////////////////////////////////////////////////////////////////////
196 /// create a Mesh0DElement and add it to the current Mesh
197 /// @return : The created Mesh0DElement
198 ///////////////////////////////////////////////////////////////////////////////
199 SMDS_Mesh0DElement* SMDS_Mesh::Add0DElementWithID(int idnode, int ID)
201 SMDS_MeshNode * node = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode);
202 if (!node) return NULL;
203 return SMDS_Mesh::Add0DElementWithID(node, ID);
206 ///////////////////////////////////////////////////////////////////////////////
207 /// create a Mesh0DElement and add it to the current Mesh
208 /// @return : The created Mesh0DElement
209 ///////////////////////////////////////////////////////////////////////////////
210 SMDS_Mesh0DElement* SMDS_Mesh::Add0DElement(const SMDS_MeshNode * node)
212 return SMDS_Mesh::Add0DElementWithID(node, myElementIDFactory->GetFreeID());
215 ///////////////////////////////////////////////////////////////////////////////
216 /// Create a new Mesh0DElement and at it to the mesh
217 /// @param idnode ID of the node
218 /// @param ID ID of the 0D element to create
219 /// @return The created 0D element or NULL if an element with this
220 /// ID already exists or if input node is not found.
221 ///////////////////////////////////////////////////////////////////////////////
222 SMDS_Mesh0DElement* SMDS_Mesh::Add0DElementWithID(const SMDS_MeshNode * n, int ID)
226 //if (my0DElements.Extent() % CHECKMEMORY_INTERVAL == 0) CheckMemory();
227 //MESSAGE("Add0DElementWithID" << ID)
228 SMDS_Mesh0DElement * el0d = new SMDS_Mesh0DElement(n);
229 if (myElementIDFactory->BindID(ID, el0d)) {
230 SMDS_MeshNode *node = const_cast<SMDS_MeshNode*>(n);
231 //node->AddInverseElement(el0d);// --- fait avec BindID
232 adjustmyCellsCapacity(ID);
234 myInfo.myNb0DElements++;
242 ///////////////////////////////////////////////////////////////////////////////
243 /// create a MeshEdge and add it to the current Mesh
244 /// @return : The created MeshEdge
245 ///////////////////////////////////////////////////////////////////////////////
247 SMDS_MeshEdge* SMDS_Mesh::AddEdgeWithID(int idnode1, int idnode2, int ID)
249 SMDS_MeshNode * node1 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode1);
250 SMDS_MeshNode * node2 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode2);
251 if(!node1 || !node2) return NULL;
252 return SMDS_Mesh::AddEdgeWithID(node1, node2, ID);
255 ///////////////////////////////////////////////////////////////////////////////
256 /// create a MeshEdge and add it to the current Mesh
257 /// @return : The created MeshEdge
258 ///////////////////////////////////////////////////////////////////////////////
260 SMDS_MeshEdge* SMDS_Mesh::AddEdge(const SMDS_MeshNode * node1,
261 const SMDS_MeshNode * node2)
263 return SMDS_Mesh::AddEdgeWithID(node1, node2, myElementIDFactory->GetFreeID());
266 ///////////////////////////////////////////////////////////////////////////////
267 /// Create a new edge and at it to the mesh
268 /// @param idnode1 ID of the first node
269 /// @param idnode2 ID of the second node
270 /// @param ID ID of the edge to create
271 /// @return The created edge or NULL if an element with this ID already exists or
272 /// if input nodes are not found.
273 ///////////////////////////////////////////////////////////////////////////////
275 SMDS_MeshEdge* SMDS_Mesh::AddEdgeWithID(const SMDS_MeshNode * n1,
276 const SMDS_MeshNode * n2,
279 if ( !n1 || !n2 ) return 0;
281 //if ( myEdges.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
282 //MESSAGE("AddEdgeWithID " << ID)
283 SMDS_MeshEdge * edge=new SMDS_MeshEdge(n1,n2);
284 adjustmyCellsCapacity(ID);
288 if (edge && !registerElement(ID, edge))
290 RemoveElement(edge, false);
296 ///////////////////////////////////////////////////////////////////////////////
297 /// Add a triangle defined by its nodes. An ID is automatically affected to the
299 ///////////////////////////////////////////////////////////////////////////////
301 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshNode * n1,
302 const SMDS_MeshNode * n2,
303 const SMDS_MeshNode * n3)
305 return SMDS_Mesh::AddFaceWithID(n1,n2,n3, myElementIDFactory->GetFreeID());
308 ///////////////////////////////////////////////////////////////////////////////
309 /// Add a triangle defined by its nodes IDs
310 ///////////////////////////////////////////////////////////////////////////////
312 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(int idnode1, int idnode2, int idnode3, int ID)
314 SMDS_MeshNode * node1 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode1);
315 SMDS_MeshNode * node2 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode2);
316 SMDS_MeshNode * node3 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode3);
317 if(!node1 || !node2 || !node3) return NULL;
318 return SMDS_Mesh::AddFaceWithID(node1, node2, node3, ID);
321 ///////////////////////////////////////////////////////////////////////////////
322 /// Add a triangle defined by its nodes
323 ///////////////////////////////////////////////////////////////////////////////
325 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshNode * n1,
326 const SMDS_MeshNode * n2,
327 const SMDS_MeshNode * n3,
330 //MESSAGE("AddFaceWithID " << ID)
331 SMDS_MeshFace * face=createTriangle(n1, n2, n3);
333 if (face && !registerElement(ID, face)) {
334 RemoveElement(face, false);
340 ///////////////////////////////////////////////////////////////////////////////
341 /// Add a quadrangle defined by its nodes. An ID is automatically affected to the
343 ///////////////////////////////////////////////////////////////////////////////
345 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshNode * n1,
346 const SMDS_MeshNode * n2,
347 const SMDS_MeshNode * n3,
348 const SMDS_MeshNode * n4)
350 return SMDS_Mesh::AddFaceWithID(n1,n2,n3, n4, myElementIDFactory->GetFreeID());
353 ///////////////////////////////////////////////////////////////////////////////
354 /// Add a quadrangle defined by its nodes IDs
355 ///////////////////////////////////////////////////////////////////////////////
357 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(int idnode1,
363 SMDS_MeshNode *node1, *node2, *node3, *node4;
364 node1 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode1);
365 node2 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode2);
366 node3 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode3);
367 node4 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode4);
368 if(!node1 || !node2 || !node3 || !node4) return NULL;
369 return SMDS_Mesh::AddFaceWithID(node1, node2, node3, node4, ID);
372 ///////////////////////////////////////////////////////////////////////////////
373 /// Add a quadrangle defined by its nodes
374 ///////////////////////////////////////////////////////////////////////////////
376 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshNode * n1,
377 const SMDS_MeshNode * n2,
378 const SMDS_MeshNode * n3,
379 const SMDS_MeshNode * n4,
382 //MESSAGE("AddFaceWithID " << ID);
383 SMDS_MeshFace * face=createQuadrangle(n1, n2, n3, n4, ID);
385 if (face && !registerElement(ID, face)) {
386 RemoveElement(face, false);
392 ///////////////////////////////////////////////////////////////////////////////
393 /// Add a triangle defined by its edges. An ID is automatically assigned to the
395 ///////////////////////////////////////////////////////////////////////////////
397 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshEdge * e1,
398 const SMDS_MeshEdge * e2,
399 const SMDS_MeshEdge * e3)
401 if (!hasConstructionEdges())
403 //MESSAGE("AddFaceWithID");
404 return AddFaceWithID(e1,e2,e3, myElementIDFactory->GetFreeID());
407 ///////////////////////////////////////////////////////////////////////////////
408 /// Add a triangle defined by its edges
409 ///////////////////////////////////////////////////////////////////////////////
411 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshEdge * e1,
412 const SMDS_MeshEdge * e2,
413 const SMDS_MeshEdge * e3,
416 if (!hasConstructionEdges())
418 if ( !e1 || !e2 || !e3 ) return 0;
420 //if ( myFaces.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
421 //MESSAGE("AddFaceWithID" << ID);
423 SMDS_MeshFace * face = new SMDS_FaceOfEdges(e1,e2,e3);
424 adjustmyCellsCapacity(ID);
426 myInfo.myNbTriangles++;
428 if (!registerElement(ID, face)) {
429 RemoveElement(face, false);
435 ///////////////////////////////////////////////////////////////////////////////
436 /// Add a quadrangle defined by its edges. An ID is automatically assigned to the
438 ///////////////////////////////////////////////////////////////////////////////
440 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshEdge * e1,
441 const SMDS_MeshEdge * e2,
442 const SMDS_MeshEdge * e3,
443 const SMDS_MeshEdge * e4)
445 if (!hasConstructionEdges())
447 //MESSAGE("AddFaceWithID" );
448 return AddFaceWithID(e1,e2,e3,e4, myElementIDFactory->GetFreeID());
451 ///////////////////////////////////////////////////////////////////////////////
452 /// Add a quadrangle defined by its edges
453 ///////////////////////////////////////////////////////////////////////////////
455 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshEdge * e1,
456 const SMDS_MeshEdge * e2,
457 const SMDS_MeshEdge * e3,
458 const SMDS_MeshEdge * e4,
461 if (!hasConstructionEdges())
463 //MESSAGE("AddFaceWithID" << ID);
464 if ( !e1 || !e2 || !e3 || !e4 ) return 0;
465 //if ( myFaces.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
466 SMDS_MeshFace * face = new SMDS_FaceOfEdges(e1,e2,e3,e4);
467 adjustmyCellsCapacity(ID);
469 myInfo.myNbQuadrangles++;
471 if (!registerElement(ID, face))
473 RemoveElement(face, false);
479 ///////////////////////////////////////////////////////////////////////////////
480 ///Create a new tetrahedron and add it to the mesh.
481 ///@return The created tetrahedron
482 ///////////////////////////////////////////////////////////////////////////////
484 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
485 const SMDS_MeshNode * n2,
486 const SMDS_MeshNode * n3,
487 const SMDS_MeshNode * n4)
489 int ID = myElementIDFactory->GetFreeID();
490 //MESSAGE("AddVolumeWithID " << ID);
491 SMDS_MeshVolume * v = SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, ID);
492 if(v==NULL) myElementIDFactory->ReleaseID(ID);
496 ///////////////////////////////////////////////////////////////////////////////
497 ///Create a new tetrahedron and add it to the mesh.
498 ///@param ID The ID of the new volume
499 ///@return The created tetrahedron or NULL if an element with this ID already exists
500 ///or if input nodes are not found.
501 ///////////////////////////////////////////////////////////////////////////////
503 SMDS_MeshVolume * SMDS_Mesh::AddVolumeWithID(int idnode1,
509 //MESSAGE("AddVolumeWithID" << ID);
510 SMDS_MeshNode *node1, *node2, *node3, *node4;
511 node1 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode1);
512 node2 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode2);
513 node3 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode3);
514 node4 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode4);
515 if(!node1 || !node2 || !node3 || !node4) return NULL;
516 return SMDS_Mesh::AddVolumeWithID(node1, node2, node3, node4, ID);
519 ///////////////////////////////////////////////////////////////////////////////
520 ///Create a new tetrahedron and add it to the mesh.
521 ///@param ID The ID of the new volume
522 ///@return The created tetrahedron
523 ///////////////////////////////////////////////////////////////////////////////
525 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
526 const SMDS_MeshNode * n2,
527 const SMDS_MeshNode * n3,
528 const SMDS_MeshNode * n4,
531 //MESSAGE("AddVolumeWithID " << ID);
532 SMDS_MeshVolume* volume = 0;
533 if ( !n1 || !n2 || !n3 || !n4) return volume;
534 //if ( myVolumes.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
535 if(hasConstructionFaces()) {
536 SMDS_MeshFace * f1=FindFaceOrCreate(n1,n2,n3);
537 SMDS_MeshFace * f2=FindFaceOrCreate(n1,n2,n4);
538 SMDS_MeshFace * f3=FindFaceOrCreate(n1,n3,n4);
539 SMDS_MeshFace * f4=FindFaceOrCreate(n2,n3,n4);
540 volume=new SMDS_VolumeOfFaces(f1,f2,f3,f4);
541 adjustmyCellsCapacity(ID);
542 myCells[ID] = volume;
545 else if(hasConstructionEdges()) {
546 MESSAGE("Error : Not implemented");
550 volume=new SMDS_VolumeOfNodes(n1,n2,n3,n4);
551 adjustmyCellsCapacity(ID);
552 myCells[ID] = volume;
556 if (!registerElement(ID, volume)) {
557 RemoveElement(volume, false);
563 ///////////////////////////////////////////////////////////////////////////////
564 ///Create a new pyramid and add it to the mesh.
565 ///Nodes 1,2,3 and 4 define the base of the pyramid
566 ///@return The created pyramid
567 ///////////////////////////////////////////////////////////////////////////////
569 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
570 const SMDS_MeshNode * n2,
571 const SMDS_MeshNode * n3,
572 const SMDS_MeshNode * n4,
573 const SMDS_MeshNode * n5)
575 int ID = myElementIDFactory->GetFreeID();
576 //MESSAGE("AddVolumeWithID " << ID);
577 SMDS_MeshVolume * v = SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, ID);
578 if(v==NULL) myElementIDFactory->ReleaseID(ID);
582 ///////////////////////////////////////////////////////////////////////////////
583 ///Create a new pyramid and add it to the mesh.
584 ///Nodes 1,2,3 and 4 define the base of the pyramid
585 ///@param ID The ID of the new volume
586 ///@return The created pyramid or NULL if an element with this ID already exists
587 ///or if input nodes are not found.
588 ///////////////////////////////////////////////////////////////////////////////
590 SMDS_MeshVolume * SMDS_Mesh::AddVolumeWithID(int idnode1,
597 //MESSAGE("AddVolumeWithID " << ID);
598 SMDS_MeshNode *node1, *node2, *node3, *node4, *node5;
599 node1 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode1);
600 node2 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode2);
601 node3 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode3);
602 node4 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode4);
603 node5 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode5);
604 if(!node1 || !node2 || !node3 || !node4 || !node5) return NULL;
605 return SMDS_Mesh::AddVolumeWithID(node1, node2, node3, node4, node5, ID);
608 ///////////////////////////////////////////////////////////////////////////////
609 ///Create a new pyramid and add it to the mesh.
610 ///Nodes 1,2,3 and 4 define the base of the pyramid
611 ///@param ID The ID of the new volume
612 ///@return The created pyramid
613 ///////////////////////////////////////////////////////////////////////////////
615 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
616 const SMDS_MeshNode * n2,
617 const SMDS_MeshNode * n3,
618 const SMDS_MeshNode * n4,
619 const SMDS_MeshNode * n5,
622 //MESSAGE("AddVolumeWithID " << ID);
623 SMDS_MeshVolume* volume = 0;
624 if ( !n1 || !n2 || !n3 || !n4 || !n5) return volume;
625 //if ( myVolumes.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
626 if(hasConstructionFaces()) {
627 SMDS_MeshFace * f1=FindFaceOrCreate(n1,n2,n3,n4);
628 SMDS_MeshFace * f2=FindFaceOrCreate(n1,n2,n5);
629 SMDS_MeshFace * f3=FindFaceOrCreate(n2,n3,n5);
630 SMDS_MeshFace * f4=FindFaceOrCreate(n3,n4,n5);
631 volume=new SMDS_VolumeOfFaces(f1,f2,f3,f4);
632 adjustmyCellsCapacity(ID);
633 myCells[ID] = volume;
634 myInfo.myNbPyramids++;
636 else if(hasConstructionEdges()) {
637 MESSAGE("Error : Not implemented");
641 volume=new SMDS_VolumeOfNodes(n1,n2,n3,n4,n5);
642 adjustmyCellsCapacity(ID);
643 myCells[ID] = volume;
644 myInfo.myNbPyramids++;
647 if (!registerElement(ID, volume)) {
648 RemoveElement(volume, false);
654 ///////////////////////////////////////////////////////////////////////////////
655 ///Create a new prism and add it to the mesh.
656 ///Nodes 1,2,3 is a triangle and 1,2,5,4 a quadrangle.
657 ///@return The created prism
658 ///////////////////////////////////////////////////////////////////////////////
660 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
661 const SMDS_MeshNode * n2,
662 const SMDS_MeshNode * n3,
663 const SMDS_MeshNode * n4,
664 const SMDS_MeshNode * n5,
665 const SMDS_MeshNode * n6)
667 int ID = myElementIDFactory->GetFreeID();
668 //MESSAGE("AddVolumeWithID " << ID);
669 SMDS_MeshVolume * v = SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n6, ID);
670 if(v==NULL) myElementIDFactory->ReleaseID(ID);
674 ///////////////////////////////////////////////////////////////////////////////
675 ///Create a new prism and add it to the mesh.
676 ///Nodes 1,2,3 is a triangle and 1,2,5,4 a quadrangle.
677 ///@param ID The ID of the new volume
678 ///@return The created prism or NULL if an element with this ID already exists
679 ///or if input nodes are not found.
680 ///////////////////////////////////////////////////////////////////////////////
682 SMDS_MeshVolume * SMDS_Mesh::AddVolumeWithID(int idnode1,
690 //MESSAGE("AddVolumeWithID " << ID);
691 SMDS_MeshNode *node1, *node2, *node3, *node4, *node5, *node6;
692 node1 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode1);
693 node2 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode2);
694 node3 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode3);
695 node4 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode4);
696 node5 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode5);
697 node6 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode6);
698 if(!node1 || !node2 || !node3 || !node4 || !node5 || !node6) return NULL;
699 return SMDS_Mesh::AddVolumeWithID(node1, node2, node3, node4, node5, node6, ID);
702 ///////////////////////////////////////////////////////////////////////////////
703 ///Create a new prism and add it to the mesh.
704 ///Nodes 1,2,3 is a triangle and 1,2,5,4 a quadrangle.
705 ///@param ID The ID of the new volume
706 ///@return The created prism
707 ///////////////////////////////////////////////////////////////////////////////
709 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
710 const SMDS_MeshNode * n2,
711 const SMDS_MeshNode * n3,
712 const SMDS_MeshNode * n4,
713 const SMDS_MeshNode * n5,
714 const SMDS_MeshNode * n6,
717 //MESSAGE("AddVolumeWithID " << ID);
718 SMDS_MeshVolume* volume = 0;
719 if ( !n1 || !n2 || !n3 || !n4 || !n5 || !n6) return volume;
720 //if ( myVolumes.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
721 if(hasConstructionFaces()) {
722 SMDS_MeshFace * f1=FindFaceOrCreate(n1,n2,n3);
723 SMDS_MeshFace * f2=FindFaceOrCreate(n4,n5,n6);
724 SMDS_MeshFace * f3=FindFaceOrCreate(n1,n4,n5,n2);
725 SMDS_MeshFace * f4=FindFaceOrCreate(n2,n5,n6,n3);
726 SMDS_MeshFace * f5=FindFaceOrCreate(n3,n6,n4,n1);
727 volume=new SMDS_VolumeOfFaces(f1,f2,f3,f4,f5);
728 adjustmyCellsCapacity(ID);
729 myCells[ID] = volume;
732 else if(hasConstructionEdges()) {
733 MESSAGE("Error : Not implemented");
737 volume=new SMDS_VolumeOfNodes(n1,n2,n3,n4,n5,n6);
738 adjustmyCellsCapacity(ID);
739 myCells[ID] = volume;
743 if (!registerElement(ID, volume)) {
744 RemoveElement(volume, false);
750 ///////////////////////////////////////////////////////////////////////////////
751 ///Create a new hexahedron and add it to the mesh.
752 ///Nodes 1,2,3,4 and 5,6,7,8 are quadrangle and 5,1 and 7,3 are an edges.
753 ///@return The created hexahedron
754 ///////////////////////////////////////////////////////////////////////////////
756 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
757 const SMDS_MeshNode * n2,
758 const SMDS_MeshNode * n3,
759 const SMDS_MeshNode * n4,
760 const SMDS_MeshNode * n5,
761 const SMDS_MeshNode * n6,
762 const SMDS_MeshNode * n7,
763 const SMDS_MeshNode * n8)
765 int ID = myElementIDFactory->GetFreeID();
766 //MESSAGE("AddVolumeWithID " << ID);
767 SMDS_MeshVolume * v = SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n6, n7, n8, ID);
768 if(v==NULL) myElementIDFactory->ReleaseID(ID);
772 ///////////////////////////////////////////////////////////////////////////////
773 ///Create a new hexahedron and add it to the mesh.
774 ///Nodes 1,2,3,4 and 5,6,7,8 are quadrangle and 5,1 and 7,3 are an edges.
775 ///@param ID The ID of the new volume
776 ///@return The created hexahedron or NULL if an element with this ID already
777 ///exists or if input nodes are not found.
778 ///////////////////////////////////////////////////////////////////////////////
780 SMDS_MeshVolume * SMDS_Mesh::AddVolumeWithID(int idnode1,
790 //MESSAGE("AddVolumeWithID " << ID);
791 SMDS_MeshNode *node1, *node2, *node3, *node4, *node5, *node6, *node7, *node8;
792 node1 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode1);
793 node2 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode2);
794 node3 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode3);
795 node4 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode4);
796 node5 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode5);
797 node6 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode6);
798 node7 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode7);
799 node8 = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(idnode8);
800 if(!node1 || !node2 || !node3 || !node4 || !node5 || !node6 || !node7 || !node8)
802 return SMDS_Mesh::AddVolumeWithID(node1, node2, node3, node4, node5, node6,
806 ///////////////////////////////////////////////////////////////////////////////
807 ///Create a new hexahedron and add it to the mesh.
808 ///Nodes 1,2,3,4 and 5,6,7,8 are quadrangle and 5,1 and 7,3 are an edges.
809 ///@param ID The ID of the new volume
810 ///@return The created prism or NULL if an element with this ID already exists
811 ///or if input nodes are not found.
812 ///////////////////////////////////////////////////////////////////////////////
814 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
815 const SMDS_MeshNode * n2,
816 const SMDS_MeshNode * n3,
817 const SMDS_MeshNode * n4,
818 const SMDS_MeshNode * n5,
819 const SMDS_MeshNode * n6,
820 const SMDS_MeshNode * n7,
821 const SMDS_MeshNode * n8,
824 //MESSAGE("AddVolumeWithID " << ID);
825 SMDS_MeshVolume* volume = 0;
826 if ( !n1 || !n2 || !n3 || !n4 || !n5 || !n6 || !n7 || !n8) return volume;
827 //if ( myVolumes.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
828 if(hasConstructionFaces()) {
829 SMDS_MeshFace * f1=FindFaceOrCreate(n1,n2,n3,n4);
830 SMDS_MeshFace * f2=FindFaceOrCreate(n5,n6,n7,n8);
831 SMDS_MeshFace * f3=FindFaceOrCreate(n1,n4,n8,n5);
832 SMDS_MeshFace * f4=FindFaceOrCreate(n1,n2,n6,n5);
833 SMDS_MeshFace * f5=FindFaceOrCreate(n2,n3,n7,n6);
834 SMDS_MeshFace * f6=FindFaceOrCreate(n3,n4,n8,n7);
835 volume=new SMDS_VolumeOfFaces(f1,f2,f3,f4,f5,f6);
836 adjustmyCellsCapacity(ID);
837 myCells[ID] = volume;
840 else if(hasConstructionEdges()) {
841 MESSAGE("Error : Not implemented");
845 // volume=new SMDS_HexahedronOfNodes(n1,n2,n3,n4,n5,n6,n7,n8);
846 volume=new SMDS_VolumeOfNodes(n1,n2,n3,n4,n5,n6,n7,n8);
847 adjustmyCellsCapacity(ID);
848 myCells[ID] = volume;
852 if (!registerElement(ID, volume)) {
853 RemoveElement(volume, false);
859 ///////////////////////////////////////////////////////////////////////////////
860 ///Create a new tetrahedron defined by its faces and add it to the mesh.
861 ///@return The created tetrahedron
862 ///////////////////////////////////////////////////////////////////////////////
864 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshFace * f1,
865 const SMDS_MeshFace * f2,
866 const SMDS_MeshFace * f3,
867 const SMDS_MeshFace * f4)
869 //MESSAGE("AddVolumeWithID");
870 if (!hasConstructionFaces())
872 return AddVolumeWithID(f1,f2,f3,f4, myElementIDFactory->GetFreeID());
875 ///////////////////////////////////////////////////////////////////////////////
876 ///Create a new tetrahedron defined by its faces and add it to the mesh.
877 ///@param ID The ID of the new volume
878 ///@return The created tetrahedron
879 ///////////////////////////////////////////////////////////////////////////////
881 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshFace * f1,
882 const SMDS_MeshFace * f2,
883 const SMDS_MeshFace * f3,
884 const SMDS_MeshFace * f4,
887 //MESSAGE("AddVolumeWithID" << ID);
888 if (!hasConstructionFaces())
890 if ( !f1 || !f2 || !f3 || !f4) return 0;
891 //if ( myVolumes.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
892 SMDS_MeshVolume * volume = new SMDS_VolumeOfFaces(f1,f2,f3,f4);
893 adjustmyCellsCapacity(ID);
894 myCells[ID] = volume;
897 if (!registerElement(ID, volume)) {
898 RemoveElement(volume, false);
904 ///////////////////////////////////////////////////////////////////////////////
905 ///Create a new pyramid defined by its faces and add it to the mesh.
906 ///@return The created pyramid
907 ///////////////////////////////////////////////////////////////////////////////
909 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshFace * f1,
910 const SMDS_MeshFace * f2,
911 const SMDS_MeshFace * f3,
912 const SMDS_MeshFace * f4,
913 const SMDS_MeshFace * f5)
915 //MESSAGE("AddVolumeWithID");
916 if (!hasConstructionFaces())
918 return AddVolumeWithID(f1,f2,f3,f4,f5, myElementIDFactory->GetFreeID());
921 ///////////////////////////////////////////////////////////////////////////////
922 ///Create a new pyramid defined by its faces and add it to the mesh.
923 ///@param ID The ID of the new volume
924 ///@return The created pyramid
925 ///////////////////////////////////////////////////////////////////////////////
927 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshFace * f1,
928 const SMDS_MeshFace * f2,
929 const SMDS_MeshFace * f3,
930 const SMDS_MeshFace * f4,
931 const SMDS_MeshFace * f5,
934 //MESSAGE("AddVolumeWithID" << ID);
935 if (!hasConstructionFaces())
937 if ( !f1 || !f2 || !f3 || !f4 || !f5) return 0;
938 //if ( myVolumes.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
939 SMDS_MeshVolume * volume = new SMDS_VolumeOfFaces(f1,f2,f3,f4,f5);
940 adjustmyCellsCapacity(ID);
941 myCells[ID] = volume;
942 myInfo.myNbPyramids++;
944 if (!registerElement(ID, volume)) {
945 RemoveElement(volume, false);
951 ///////////////////////////////////////////////////////////////////////////////
952 ///Create a new prism defined by its faces and add it to the mesh.
953 ///@return The created prism
954 ///////////////////////////////////////////////////////////////////////////////
956 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshFace * f1,
957 const SMDS_MeshFace * f2,
958 const SMDS_MeshFace * f3,
959 const SMDS_MeshFace * f4,
960 const SMDS_MeshFace * f5,
961 const SMDS_MeshFace * f6)
963 //MESSAGE("AddVolumeWithID" );
964 if (!hasConstructionFaces())
966 return AddVolumeWithID(f1,f2,f3,f4,f5,f6, myElementIDFactory->GetFreeID());
969 ///////////////////////////////////////////////////////////////////////////////
970 ///Create a new prism defined by its faces and add it to the mesh.
971 ///@param ID The ID of the new volume
972 ///@return The created prism
973 ///////////////////////////////////////////////////////////////////////////////
975 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshFace * f1,
976 const SMDS_MeshFace * f2,
977 const SMDS_MeshFace * f3,
978 const SMDS_MeshFace * f4,
979 const SMDS_MeshFace * f5,
980 const SMDS_MeshFace * f6,
983 //MESSAGE("AddVolumeWithID" << ID);
984 if (!hasConstructionFaces())
986 if ( !f1 || !f2 || !f3 || !f4 || !f5 || !f6) return 0;
987 //if ( myVolumes.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
988 SMDS_MeshVolume * volume = new SMDS_VolumeOfFaces(f1,f2,f3,f4,f5,f6);
989 adjustmyCellsCapacity(ID);
990 myCells[ID] = volume;
993 if (!registerElement(ID, volume)) {
994 RemoveElement(volume, false);
1000 ///////////////////////////////////////////////////////////////////////////////
1001 /// Add a polygon defined by its nodes IDs
1002 ///////////////////////////////////////////////////////////////////////////////
1004 SMDS_MeshFace* SMDS_Mesh::AddPolygonalFaceWithID (std::vector<int> nodes_ids,
1007 int nbNodes = nodes_ids.size();
1008 std::vector<const SMDS_MeshNode*> nodes (nbNodes);
1009 for (int i = 0; i < nbNodes; i++) {
1010 nodes[i] = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(nodes_ids[i]);
1011 if (!nodes[i]) return NULL;
1013 return SMDS_Mesh::AddPolygonalFaceWithID(nodes, ID);
1016 ///////////////////////////////////////////////////////////////////////////////
1017 /// Add a polygon defined by its nodes
1018 ///////////////////////////////////////////////////////////////////////////////
1020 SMDS_MeshFace* SMDS_Mesh::AddPolygonalFaceWithID
1021 (std::vector<const SMDS_MeshNode*> nodes,
1024 SMDS_MeshFace * face;
1026 //if ( myFaces.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
1027 if (hasConstructionEdges())
1029 MESSAGE("Error : Not implemented");
1034 for ( int i = 0; i < nodes.size(); ++i )
1035 if ( !nodes[ i ] ) return 0;
1036 face = new SMDS_PolygonalFaceOfNodes(nodes);
1037 adjustmyCellsCapacity(ID);
1039 myInfo.myNbPolygons++;
1042 if (!registerElement(ID, face)) {
1043 RemoveElement(face, false);
1049 ///////////////////////////////////////////////////////////////////////////////
1050 /// Add a polygon defined by its nodes.
1051 /// An ID is automatically affected to the created face.
1052 ///////////////////////////////////////////////////////////////////////////////
1054 SMDS_MeshFace* SMDS_Mesh::AddPolygonalFace (std::vector<const SMDS_MeshNode*> nodes)
1056 return SMDS_Mesh::AddPolygonalFaceWithID(nodes, myElementIDFactory->GetFreeID());
1059 ///////////////////////////////////////////////////////////////////////////////
1060 /// Create a new polyhedral volume and add it to the mesh.
1061 /// @param ID The ID of the new volume
1062 /// @return The created volume or NULL if an element with this ID already exists
1063 /// or if input nodes are not found.
1064 ///////////////////////////////////////////////////////////////////////////////
1066 SMDS_MeshVolume * SMDS_Mesh::AddPolyhedralVolumeWithID
1067 (std::vector<int> nodes_ids,
1068 std::vector<int> quantities,
1071 int nbNodes = nodes_ids.size();
1072 std::vector<const SMDS_MeshNode*> nodes (nbNodes);
1073 for (int i = 0; i < nbNodes; i++) {
1074 nodes[i] = (SMDS_MeshNode *)myNodeIDFactory->MeshElement(nodes_ids[i]);
1075 if (!nodes[i]) return NULL;
1077 return SMDS_Mesh::AddPolyhedralVolumeWithID(nodes, quantities, ID);
1080 ///////////////////////////////////////////////////////////////////////////////
1081 /// Create a new polyhedral volume and add it to the mesh.
1082 /// @param ID The ID of the new volume
1083 /// @return The created volume
1084 ///////////////////////////////////////////////////////////////////////////////
1086 SMDS_MeshVolume* SMDS_Mesh::AddPolyhedralVolumeWithID
1087 (std::vector<const SMDS_MeshNode*> nodes,
1088 std::vector<int> quantities,
1091 SMDS_MeshVolume* volume;
1092 //if ( myVolumes.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
1093 if (hasConstructionFaces()) {
1094 MESSAGE("Error : Not implemented");
1096 } else if (hasConstructionEdges()) {
1097 MESSAGE("Error : Not implemented");
1100 for ( int i = 0; i < nodes.size(); ++i )
1101 if ( !nodes[ i ] ) return 0;
1102 volume = new SMDS_PolyhedralVolumeOfNodes(nodes, quantities);
1103 adjustmyCellsCapacity(ID);
1104 myCells[ID] = volume;
1105 myInfo.myNbPolyhedrons++;
1108 if (!registerElement(ID, volume)) {
1109 RemoveElement(volume, false);
1115 ///////////////////////////////////////////////////////////////////////////////
1116 /// Create a new polyhedral volume and add it to the mesh.
1117 /// @return The created volume
1118 ///////////////////////////////////////////////////////////////////////////////
1120 SMDS_MeshVolume* SMDS_Mesh::AddPolyhedralVolume
1121 (std::vector<const SMDS_MeshNode*> nodes,
1122 std::vector<int> quantities)
1124 int ID = myElementIDFactory->GetFreeID();
1125 SMDS_MeshVolume * v = SMDS_Mesh::AddPolyhedralVolumeWithID(nodes, quantities, ID);
1126 if (v == NULL) myElementIDFactory->ReleaseID(ID);
1130 ///////////////////////////////////////////////////////////////////////////////
1131 /// Registers element with the given ID, maintains inverse connections
1132 ///////////////////////////////////////////////////////////////////////////////
1133 bool SMDS_Mesh::registerElement(int ID, SMDS_MeshElement * element)
1135 //MESSAGE("registerElement " << ID)
1136 if (myElementIDFactory->BindID(ID, element)) {
1139 MESSAGE("BindID " << ID << " false!---------------");
1143 ///////////////////////////////////////////////////////////////////////////////
1144 /// Return the node whose ID is 'ID'.
1145 ///////////////////////////////////////////////////////////////////////////////
1146 const SMDS_MeshNode * SMDS_Mesh::FindNode(int ID) const
1148 if (ID <0 || ID >= myNodes.size())
1150 return (const SMDS_MeshNode *)myNodes[ID];
1153 ///////////////////////////////////////////////////////////////////////////////
1154 ///Create a triangle and add it to the current mesh. This methode do not bind a
1155 ///ID to the create triangle.
1156 ///////////////////////////////////////////////////////////////////////////////
1157 SMDS_MeshFace * SMDS_Mesh::createTriangle(const SMDS_MeshNode * node1,
1158 const SMDS_MeshNode * node2,
1159 const SMDS_MeshNode * node3)
1161 if ( !node1 || !node2 || !node3) return 0;
1162 // if ( myFaces.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
1163 if(hasConstructionEdges())
1165 SMDS_MeshEdge *edge1, *edge2, *edge3;
1166 edge1=FindEdgeOrCreate(node1,node2);
1167 edge2=FindEdgeOrCreate(node2,node3);
1168 edge3=FindEdgeOrCreate(node3,node1);
1170 int ID = myElementIDFactory->GetFreeID(); // -PR- voir si on range cet element
1171 SMDS_MeshFace * face = new SMDS_FaceOfEdges(edge1,edge2,edge3);
1172 adjustmyCellsCapacity(ID);
1174 myInfo.myNbTriangles++;
1179 int ID = myElementIDFactory->GetFreeID(); // -PR- voir si on range cet element
1180 SMDS_MeshFace * face = new SMDS_FaceOfNodes(node1,node2,node3);
1181 adjustmyCellsCapacity(ID);
1183 myInfo.myNbTriangles++;
1188 ///////////////////////////////////////////////////////////////////////////////
1189 ///Create a quadrangle and add it to the current mesh. This methode do not bind
1190 ///a ID to the create triangle.
1191 ///////////////////////////////////////////////////////////////////////////////
1192 SMDS_MeshFace * SMDS_Mesh::createQuadrangle(const SMDS_MeshNode * node1,
1193 const SMDS_MeshNode * node2,
1194 const SMDS_MeshNode * node3,
1195 const SMDS_MeshNode * node4,
1198 if ( !node1 || !node2 || !node3 || !node4 ) return 0;
1199 // if ( myFaces.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
1200 if(hasConstructionEdges())
1202 //MESSAGE("createQuadrangle hasConstructionEdges "<< ID);
1203 SMDS_MeshEdge *edge1, *edge2, *edge3, *edge4;
1204 edge1=FindEdgeOrCreate(node1,node2);
1205 edge2=FindEdgeOrCreate(node2,node3);
1206 edge3=FindEdgeOrCreate(node3,node4);
1207 edge4=FindEdgeOrCreate(node4,node1);
1209 SMDS_MeshFace * face = new SMDS_FaceOfEdges(edge1,edge2,edge3,edge4);
1210 adjustmyCellsCapacity(ID);
1212 myInfo.myNbQuadrangles++;
1217 //MESSAGE("createQuadrangle " << ID);
1218 SMDS_MeshFace * face = new SMDS_FaceOfNodes(node1,node2,node3,node4);
1219 adjustmyCellsCapacity(ID);
1221 myInfo.myNbQuadrangles++;
1226 ///////////////////////////////////////////////////////////////////////////////
1227 /// Remove a node and all the elements which own this node
1228 ///////////////////////////////////////////////////////////////////////////////
1230 void SMDS_Mesh::RemoveNode(const SMDS_MeshNode * node)
1232 MESSAGE("RemoveNode");
1233 RemoveElement(node, true);
1236 ///////////////////////////////////////////////////////////////////////////////
1237 /// Remove an edge and all the elements which own this edge
1238 ///////////////////////////////////////////////////////////////////////////////
1240 void SMDS_Mesh::Remove0DElement(const SMDS_Mesh0DElement * elem0d)
1242 MESSAGE("Remove0DElement");
1243 RemoveElement(elem0d,true);
1246 ///////////////////////////////////////////////////////////////////////////////
1247 /// Remove an edge and all the elements which own this edge
1248 ///////////////////////////////////////////////////////////////////////////////
1250 void SMDS_Mesh::RemoveEdge(const SMDS_MeshEdge * edge)
1252 MESSAGE("RemoveEdge");
1253 RemoveElement(edge,true);
1256 ///////////////////////////////////////////////////////////////////////////////
1257 /// Remove an face and all the elements which own this face
1258 ///////////////////////////////////////////////////////////////////////////////
1260 void SMDS_Mesh::RemoveFace(const SMDS_MeshFace * face)
1262 MESSAGE("RemoveFace");
1263 RemoveElement(face, true);
1266 ///////////////////////////////////////////////////////////////////////////////
1268 ///////////////////////////////////////////////////////////////////////////////
1270 void SMDS_Mesh::RemoveVolume(const SMDS_MeshVolume * volume)
1272 MESSAGE("RemoveVolume");
1273 RemoveElement(volume, true);
1276 //=======================================================================
1277 //function : RemoveFromParent
1279 //=======================================================================
1281 bool SMDS_Mesh::RemoveFromParent()
1283 if (myParent==NULL) return false;
1284 else return (myParent->RemoveSubMesh(this));
1287 //=======================================================================
1288 //function : RemoveSubMesh
1290 //=======================================================================
1292 bool SMDS_Mesh::RemoveSubMesh(const SMDS_Mesh * aMesh)
1296 list<SMDS_Mesh *>::iterator itmsh=myChildren.begin();
1297 for (; itmsh!=myChildren.end() && !found; itmsh++)
1299 SMDS_Mesh * submesh = *itmsh;
1300 if (submesh == aMesh)
1303 myChildren.erase(itmsh);
1310 //=======================================================================
1311 //function : ChangeElementNodes
1313 //=======================================================================
1315 bool SMDS_Mesh::ChangeElementNodes(const SMDS_MeshElement * element,
1316 const SMDS_MeshNode * nodes[],
1319 // keep current nodes of elem
1320 set<const SMDS_MeshElement*> oldNodes;
1321 SMDS_ElemIteratorPtr itn = element->nodesIterator();
1323 oldNodes.insert( itn->next() );
1325 if ( !element->IsPoly() )
1326 myInfo.remove( element ); // element may change type
1330 SMDS_MeshElement* elem = const_cast<SMDS_MeshElement*>(element);
1331 switch ( elem->GetType() )
1333 case SMDSAbs_0DElement: {
1334 if ( SMDS_Mesh0DElement* elem0d = dynamic_cast<SMDS_Mesh0DElement*>( elem ))
1335 Ok = elem0d->ChangeNode( nodes[0] );
1338 case SMDSAbs_Edge: {
1339 if ( nbnodes == 2 ) {
1340 if ( SMDS_MeshEdge* edge = dynamic_cast<SMDS_MeshEdge*>( elem ))
1341 Ok = edge->ChangeNodes( nodes[0], nodes[1] );
1343 else if ( nbnodes == 3 ) {
1344 if ( SMDS_QuadraticEdge* edge = dynamic_cast<SMDS_QuadraticEdge*>( elem ))
1345 Ok = edge->ChangeNodes( nodes[0], nodes[1], nodes[2] );
1349 case SMDSAbs_Face: {
1350 if ( SMDS_FaceOfNodes* face = dynamic_cast<SMDS_FaceOfNodes*>( elem ))
1351 Ok = face->ChangeNodes( nodes, nbnodes );
1353 if ( SMDS_QuadraticFaceOfNodes* QF = dynamic_cast<SMDS_QuadraticFaceOfNodes*>( elem ))
1354 Ok = QF->ChangeNodes( nodes, nbnodes );
1356 if (SMDS_PolygonalFaceOfNodes* face = dynamic_cast<SMDS_PolygonalFaceOfNodes*>(elem))
1357 Ok = face->ChangeNodes(nodes, nbnodes);
1360 case SMDSAbs_Volume: {
1361 if ( SMDS_VolumeOfNodes* vol = dynamic_cast<SMDS_VolumeOfNodes*>( elem ))
1362 Ok = vol->ChangeNodes( nodes, nbnodes );
1364 if ( SMDS_QuadraticVolumeOfNodes* QV = dynamic_cast<SMDS_QuadraticVolumeOfNodes*>( elem ))
1365 Ok = QV->ChangeNodes( nodes, nbnodes );
1369 MESSAGE ( "WRONG ELEM TYPE");
1372 if ( Ok ) { // update InverseElements
1374 set<const SMDS_MeshElement*>::iterator it;
1376 // AddInverseElement to new nodes
1377 for ( int i = 0; i < nbnodes; i++ ) {
1378 it = oldNodes.find( nodes[i] );
1379 if ( it == oldNodes.end() )
1381 const_cast<SMDS_MeshNode*>( nodes[i] )->AddInverseElement( elem );
1383 // remove from oldNodes a node that remains in elem
1384 oldNodes.erase( it );
1386 // RemoveInverseElement from the nodes removed from elem
1387 for ( it = oldNodes.begin(); it != oldNodes.end(); it++ )
1389 SMDS_MeshNode * n = static_cast<SMDS_MeshNode *>
1390 (const_cast<SMDS_MeshElement *>( *it ));
1391 n->RemoveInverseElement( elem );
1395 if ( !element->IsPoly() )
1396 myInfo.add( element ); // element may change type
1401 //=======================================================================
1402 //function : ChangePolyhedronNodes
1403 //purpose : to change nodes of polyhedral volume
1404 //=======================================================================
1405 bool SMDS_Mesh::ChangePolyhedronNodes (const SMDS_MeshElement * elem,
1406 const vector<const SMDS_MeshNode*>& nodes,
1407 const vector<int> & quantities)
1409 if (elem->GetType() != SMDSAbs_Volume) {
1410 MESSAGE("WRONG ELEM TYPE");
1414 const SMDS_PolyhedralVolumeOfNodes* vol = dynamic_cast<const SMDS_PolyhedralVolumeOfNodes*>(elem);
1419 // keep current nodes of elem
1420 set<const SMDS_MeshElement*> oldNodes;
1421 SMDS_ElemIteratorPtr itn = elem->nodesIterator();
1422 while (itn->more()) {
1423 oldNodes.insert(itn->next());
1427 bool Ok = const_cast<SMDS_PolyhedralVolumeOfNodes*>(vol)->ChangeNodes(nodes, quantities);
1432 // update InverseElements
1434 // AddInverseElement to new nodes
1435 int nbnodes = nodes.size();
1436 set<const SMDS_MeshElement*>::iterator it;
1437 for (int i = 0; i < nbnodes; i++) {
1438 it = oldNodes.find(nodes[i]);
1439 if (it == oldNodes.end()) {
1441 const_cast<SMDS_MeshNode*>(nodes[i])->AddInverseElement(elem);
1443 // remove from oldNodes a node that remains in elem
1448 // RemoveInverseElement from the nodes removed from elem
1449 for (it = oldNodes.begin(); it != oldNodes.end(); it++) {
1450 SMDS_MeshNode * n = static_cast<SMDS_MeshNode *>
1451 (const_cast<SMDS_MeshElement *>( *it ));
1452 n->RemoveInverseElement(elem);
1459 //=======================================================================
1460 //function : Find0DElement
1462 //=======================================================================
1463 const SMDS_Mesh0DElement* SMDS_Mesh::Find0DElement(int idnode) const
1465 const SMDS_MeshNode * node = FindNode(idnode);
1466 if(node == NULL) return NULL;
1467 return Find0DElement(node);
1470 const SMDS_Mesh0DElement* SMDS_Mesh::Find0DElement(const SMDS_MeshNode * node)
1472 if (!node) return 0;
1473 const SMDS_Mesh0DElement* toReturn = NULL;
1474 SMDS_ElemIteratorPtr it1 = node->GetInverseElementIterator(SMDSAbs_0DElement);
1475 while (it1->more() && (toReturn == NULL)) {
1476 const SMDS_MeshElement* e = it1->next();
1477 if (e->NbNodes() == 1) {
1478 toReturn = static_cast<const SMDS_Mesh0DElement*>(e);
1484 //=======================================================================
1485 //function : Find0DElementOrCreate
1487 //=======================================================================
1488 //SMDS_Mesh0DElement* SMDS_Mesh::Find0DElementOrCreate(const SMDS_MeshNode * node)
1490 // if (!node) return 0;
1491 // SMDS_Mesh0DElement * toReturn = NULL;
1492 // toReturn = const_cast<SMDS_Mesh0DElement*>(Find0DElement(node));
1493 // if (toReturn == NULL) {
1494 // //if (my0DElements.Extent() % CHECKMEMORY_INTERVAL == 0) CheckMemory();
1495 // toReturn = new SMDS_Mesh0DElement(node);
1496 // my0DElements.Add(toReturn);
1497 // myInfo.myNb0DElements++;
1503 //=======================================================================
1504 //function : FindEdge
1506 //=======================================================================
1508 const SMDS_MeshEdge* SMDS_Mesh::FindEdge(int idnode1, int idnode2) const
1510 const SMDS_MeshNode * node1=FindNode(idnode1);
1511 const SMDS_MeshNode * node2=FindNode(idnode2);
1512 if((node1==NULL)||(node2==NULL)) return NULL;
1513 return FindEdge(node1,node2);
1516 //#include "Profiler.h"
1517 const SMDS_MeshEdge* SMDS_Mesh::FindEdge(const SMDS_MeshNode * node1,
1518 const SMDS_MeshNode * node2)
1520 if ( !node1 ) return 0;
1521 const SMDS_MeshEdge * toReturn=NULL;
1524 SMDS_ElemIteratorPtr it1=node1->GetInverseElementIterator(SMDSAbs_Edge);
1527 while(it1->more()) {
1528 const SMDS_MeshElement * e = it1->next();
1529 if ( e->NbNodes() == 2 && e->GetNodeIndex( node2 ) >= 0 ) {
1530 toReturn = static_cast<const SMDS_MeshEdge*>( e );
1539 //=======================================================================
1540 //function : FindEdgeOrCreate
1542 //=======================================================================
1544 SMDS_MeshEdge* SMDS_Mesh::FindEdgeOrCreate(const SMDS_MeshNode * node1,
1545 const SMDS_MeshNode * node2)
1547 if ( !node1 || !node2) return 0;
1548 SMDS_MeshEdge * toReturn=NULL;
1549 toReturn=const_cast<SMDS_MeshEdge*>(FindEdge(node1,node2));
1550 if(toReturn==NULL) {
1551 //if ( myEdges.Extent() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
1552 int ID = myElementIDFactory->GetFreeID(); // -PR- voir si on range cet element
1553 adjustmyCellsCapacity(ID);
1554 toReturn=new SMDS_MeshEdge(node1,node2);
1555 myCells[ID] = toReturn;
1562 //=======================================================================
1563 //function : FindEdge
1565 //=======================================================================
1567 const SMDS_MeshEdge* SMDS_Mesh::FindEdge(int idnode1, int idnode2,
1570 const SMDS_MeshNode * node1=FindNode(idnode1);
1571 const SMDS_MeshNode * node2=FindNode(idnode2);
1572 const SMDS_MeshNode * node3=FindNode(idnode3);
1573 return FindEdge(node1,node2,node3);
1576 const SMDS_MeshEdge* SMDS_Mesh::FindEdge(const SMDS_MeshNode * node1,
1577 const SMDS_MeshNode * node2,
1578 const SMDS_MeshNode * node3)
1580 if ( !node1 ) return 0;
1581 SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Edge);
1582 while(it1->more()) {
1583 const SMDS_MeshElement * e = it1->next();
1584 if ( e->NbNodes() == 3 ) {
1585 SMDS_ElemIteratorPtr it2 = e->nodesIterator();
1586 while(it2->more()) {
1587 const SMDS_MeshElement* n = it2->next();
1597 return static_cast<const SMDS_MeshEdge *> (e);
1604 //=======================================================================
1605 //function : FindFace
1607 //=======================================================================
1609 const SMDS_MeshFace* SMDS_Mesh::FindFace(int idnode1, int idnode2,
1612 const SMDS_MeshNode * node1=FindNode(idnode1);
1613 const SMDS_MeshNode * node2=FindNode(idnode2);
1614 const SMDS_MeshNode * node3=FindNode(idnode3);
1615 return FindFace(node1, node2, node3);
1618 const SMDS_MeshFace* SMDS_Mesh::FindFace(const SMDS_MeshNode *node1,
1619 const SMDS_MeshNode *node2,
1620 const SMDS_MeshNode *node3)
1622 if ( !node1 ) return 0;
1623 SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Face);
1624 while(it1->more()) {
1625 const SMDS_MeshElement * e = it1->next();
1626 if ( e->NbNodes() == 3 ) {
1627 SMDS_ElemIteratorPtr it2 = e->nodesIterator();
1628 while(it2->more()) {
1629 const SMDS_MeshElement* n = it2->next();
1639 return static_cast<const SMDS_MeshFace *> (e);
1645 SMDS_MeshFace* SMDS_Mesh::FindFaceOrCreate(const SMDS_MeshNode *node1,
1646 const SMDS_MeshNode *node2,
1647 const SMDS_MeshNode *node3)
1649 SMDS_MeshFace * toReturn=NULL;
1650 toReturn = const_cast<SMDS_MeshFace*>(FindFace(node1,node2,node3));
1651 if(toReturn==NULL) {
1652 toReturn = createTriangle(node1,node2,node3);
1658 //=======================================================================
1659 //function : FindFace
1661 //=======================================================================
1663 const SMDS_MeshFace* SMDS_Mesh::FindFace(int idnode1, int idnode2,
1664 int idnode3, int idnode4) const
1666 const SMDS_MeshNode * node1=FindNode(idnode1);
1667 const SMDS_MeshNode * node2=FindNode(idnode2);
1668 const SMDS_MeshNode * node3=FindNode(idnode3);
1669 const SMDS_MeshNode * node4=FindNode(idnode4);
1670 return FindFace(node1, node2, node3, node4);
1673 const SMDS_MeshFace* SMDS_Mesh::FindFace(const SMDS_MeshNode *node1,
1674 const SMDS_MeshNode *node2,
1675 const SMDS_MeshNode *node3,
1676 const SMDS_MeshNode *node4)
1678 if ( !node1 ) return 0;
1679 SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Face);
1680 while(it1->more()) {
1681 const SMDS_MeshElement * e = it1->next();
1682 if ( e->NbNodes() == 4 ) {
1683 SMDS_ElemIteratorPtr it2 = e->nodesIterator();
1684 while(it2->more()) {
1685 const SMDS_MeshElement* n = it2->next();
1696 return static_cast<const SMDS_MeshFace *> (e);
1702 SMDS_MeshFace* SMDS_Mesh::FindFaceOrCreate(const SMDS_MeshNode *node1,
1703 const SMDS_MeshNode *node2,
1704 const SMDS_MeshNode *node3,
1705 const SMDS_MeshNode *node4)
1707 SMDS_MeshFace * toReturn=NULL;
1708 toReturn=const_cast<SMDS_MeshFace*>(FindFace(node1,node2,node3,node4));
1709 if(toReturn==NULL) {
1710 int ID = myElementIDFactory->GetFreeID();
1711 toReturn=createQuadrangle(node1,node2,node3,node4,ID);
1717 //=======================================================================
1718 //function : FindFace
1719 //purpose :quadratic triangle
1720 //=======================================================================
1722 const SMDS_MeshFace* SMDS_Mesh::FindFace(int idnode1, int idnode2,
1723 int idnode3, int idnode4,
1724 int idnode5, int idnode6) const
1726 const SMDS_MeshNode * node1 = FindNode(idnode1);
1727 const SMDS_MeshNode * node2 = FindNode(idnode2);
1728 const SMDS_MeshNode * node3 = FindNode(idnode3);
1729 const SMDS_MeshNode * node4 = FindNode(idnode4);
1730 const SMDS_MeshNode * node5 = FindNode(idnode5);
1731 const SMDS_MeshNode * node6 = FindNode(idnode6);
1732 return FindFace(node1, node2, node3, node4, node5, node6);
1735 const SMDS_MeshFace* SMDS_Mesh::FindFace(const SMDS_MeshNode *node1,
1736 const SMDS_MeshNode *node2,
1737 const SMDS_MeshNode *node3,
1738 const SMDS_MeshNode *node4,
1739 const SMDS_MeshNode *node5,
1740 const SMDS_MeshNode *node6)
1742 if ( !node1 ) return 0;
1743 SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Face);
1744 while(it1->more()) {
1745 const SMDS_MeshElement * e = it1->next();
1746 if ( e->NbNodes() == 6 ) {
1747 SMDS_ElemIteratorPtr it2 = e->nodesIterator();
1748 while(it2->more()) {
1749 const SMDS_MeshElement* n = it2->next();
1762 return static_cast<const SMDS_MeshFace *> (e);
1769 //=======================================================================
1770 //function : FindFace
1771 //purpose : quadratic quadrangle
1772 //=======================================================================
1774 const SMDS_MeshFace* SMDS_Mesh::FindFace(int idnode1, int idnode2,
1775 int idnode3, int idnode4,
1776 int idnode5, int idnode6,
1777 int idnode7, int idnode8) const
1779 const SMDS_MeshNode * node1 = FindNode(idnode1);
1780 const SMDS_MeshNode * node2 = FindNode(idnode2);
1781 const SMDS_MeshNode * node3 = FindNode(idnode3);
1782 const SMDS_MeshNode * node4 = FindNode(idnode4);
1783 const SMDS_MeshNode * node5 = FindNode(idnode5);
1784 const SMDS_MeshNode * node6 = FindNode(idnode6);
1785 const SMDS_MeshNode * node7 = FindNode(idnode7);
1786 const SMDS_MeshNode * node8 = FindNode(idnode8);
1787 return FindFace(node1, node2, node3, node4, node5, node6, node7, node8);
1790 const SMDS_MeshFace* SMDS_Mesh::FindFace(const SMDS_MeshNode *node1,
1791 const SMDS_MeshNode *node2,
1792 const SMDS_MeshNode *node3,
1793 const SMDS_MeshNode *node4,
1794 const SMDS_MeshNode *node5,
1795 const SMDS_MeshNode *node6,
1796 const SMDS_MeshNode *node7,
1797 const SMDS_MeshNode *node8)
1799 if ( !node1 ) return 0;
1800 SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Face);
1801 while(it1->more()) {
1802 const SMDS_MeshElement * e = it1->next();
1803 if ( e->NbNodes() == 8 ) {
1804 SMDS_ElemIteratorPtr it2 = e->nodesIterator();
1805 while(it2->more()) {
1806 const SMDS_MeshElement* n = it2->next();
1821 return static_cast<const SMDS_MeshFace *> (e);
1828 //=======================================================================
1829 //function : FindElement
1831 //=======================================================================
1833 const SMDS_MeshElement* SMDS_Mesh::FindElement(int IDelem) const
1835 if ((IDelem < 0) || IDelem >= myCells.size())
1837 return myCells[IDelem];
1840 //=======================================================================
1841 //function : FindFace
1842 //purpose : find polygon
1843 //=======================================================================
1845 const SMDS_MeshFace* SMDS_Mesh::FindFace (std::vector<int> nodes_ids) const
1847 int nbnodes = nodes_ids.size();
1848 std::vector<const SMDS_MeshNode *> poly_nodes (nbnodes);
1849 for (int inode = 0; inode < nbnodes; inode++) {
1850 const SMDS_MeshNode * node = FindNode(nodes_ids[inode]);
1851 if (node == NULL) return NULL;
1852 poly_nodes[inode] = node;
1854 return FindFace(poly_nodes);
1857 const SMDS_MeshFace* SMDS_Mesh::FindFace (std::vector<const SMDS_MeshNode *> nodes)
1859 if ( nodes.size() > 2 && nodes[0] ) {
1860 SMDS_ElemIteratorPtr itF = nodes[0]->GetInverseElementIterator(SMDSAbs_Face);
1861 while (itF->more()) {
1862 const SMDS_MeshElement* f = itF->next();
1863 if ( f->NbNodes() == nodes.size() ) {
1864 SMDS_ElemIteratorPtr it2 = f->nodesIterator();
1865 while(it2->more()) {
1866 if ( find( nodes.begin(), nodes.end(), it2->next() ) == nodes.end() ) {
1872 return static_cast<const SMDS_MeshFace *> (f);
1879 //=======================================================================
1880 //function : DumpNodes
1882 //=======================================================================
1884 void SMDS_Mesh::DumpNodes() const
1886 MESSAGE("dump nodes of mesh : ");
1887 SMDS_NodeIteratorPtr itnode=nodesIterator();
1888 while(itnode->more()) ; //MESSAGE(itnode->next());
1891 //=======================================================================
1892 //function : Dump0DElements
1894 //=======================================================================
1895 void SMDS_Mesh::Dump0DElements() const
1897 MESSAGE("dump 0D elements of mesh : ");
1898 SMDS_0DElementIteratorPtr it0d = elements0dIterator();
1899 while(it0d->more()) ; //MESSAGE(it0d->next());
1902 //=======================================================================
1903 //function : DumpEdges
1905 //=======================================================================
1907 void SMDS_Mesh::DumpEdges() const
1909 MESSAGE("dump edges of mesh : ");
1910 SMDS_EdgeIteratorPtr itedge=edgesIterator();
1911 while(itedge->more()) ; //MESSAGE(itedge->next());
1914 //=======================================================================
1915 //function : DumpFaces
1917 //=======================================================================
1919 void SMDS_Mesh::DumpFaces() const
1921 MESSAGE("dump faces of mesh : ");
1922 SMDS_FaceIteratorPtr itface=facesIterator();
1923 while(itface->more()) ; //MESSAGE(itface->next());
1926 //=======================================================================
1927 //function : DumpVolumes
1929 //=======================================================================
1931 void SMDS_Mesh::DumpVolumes() const
1933 MESSAGE("dump volumes of mesh : ");
1934 SMDS_VolumeIteratorPtr itvol=volumesIterator();
1935 while(itvol->more()) ; //MESSAGE(itvol->next());
1938 //=======================================================================
1939 //function : DebugStats
1941 //=======================================================================
1943 void SMDS_Mesh::DebugStats() const
1945 MESSAGE("Debug stats of mesh : ");
1947 MESSAGE("===== NODES ====="<<NbNodes());
1948 MESSAGE("===== 0DELEMS ====="<<Nb0DElements());
1949 MESSAGE("===== EDGES ====="<<NbEdges());
1950 MESSAGE("===== FACES ====="<<NbFaces());
1951 MESSAGE("===== VOLUMES ====="<<NbVolumes());
1953 MESSAGE("End Debug stats of mesh ");
1957 SMDS_NodeIteratorPtr itnode=nodesIterator();
1958 int sizeofnodes = 0;
1959 int sizeoffaces = 0;
1961 while(itnode->more())
1963 const SMDS_MeshNode *node = itnode->next();
1965 sizeofnodes += sizeof(*node);
1967 SMDS_ElemIteratorPtr it = node->GetInverseElementIterator();
1970 const SMDS_MeshElement *me = it->next();
1971 sizeofnodes += sizeof(me);
1975 SMDS_FaceIteratorPtr itface=facesIterator();
1976 while(itface->more())
1978 const SMDS_MeshElement *face = itface->next();
1979 sizeoffaces += sizeof(*face);
1982 MESSAGE("total size of node elements = " << sizeofnodes);;
1983 MESSAGE("total size of face elements = " << sizeoffaces);;
1988 ///////////////////////////////////////////////////////////////////////////////
1989 /// Return the number of nodes
1990 ///////////////////////////////////////////////////////////////////////////////
1991 int SMDS_Mesh::NbNodes() const
1993 return myNodes.size();
1996 ///////////////////////////////////////////////////////////////////////////////
1997 /// Return the number of 0D elements
1998 ///////////////////////////////////////////////////////////////////////////////
1999 int SMDS_Mesh::Nb0DElements() const
2001 return myInfo.Nb0DElements(); // -PR- a verfier
2004 ///////////////////////////////////////////////////////////////////////////////
2005 /// Return the number of edges (including construction edges)
2006 ///////////////////////////////////////////////////////////////////////////////
2007 int SMDS_Mesh::NbEdges() const
2009 return myInfo.NbEdges(); // -PR- a verfier
2012 ///////////////////////////////////////////////////////////////////////////////
2013 /// Return the number of faces (including construction faces)
2014 ///////////////////////////////////////////////////////////////////////////////
2015 int SMDS_Mesh::NbFaces() const
2017 return myInfo.NbFaces(); // -PR- a verfier
2020 ///////////////////////////////////////////////////////////////////////////////
2021 /// Return the number of volumes
2022 ///////////////////////////////////////////////////////////////////////////////
2023 int SMDS_Mesh::NbVolumes() const
2025 return myInfo.NbVolumes(); // -PR- a verfier
2028 ///////////////////////////////////////////////////////////////////////////////
2029 /// Return the number of child mesh of this mesh.
2030 /// Note that the tree structure of SMDS_Mesh seems to be unused in this version
2031 /// (2003-09-08) of SMESH
2032 ///////////////////////////////////////////////////////////////////////////////
2033 int SMDS_Mesh::NbSubMesh() const
2035 return myChildren.size();
2038 ///////////////////////////////////////////////////////////////////////////////
2039 /// Destroy the mesh and all its elements
2040 /// All pointer on elements owned by this mesh become illegals.
2041 ///////////////////////////////////////////////////////////////////////////////
2042 SMDS_Mesh::~SMDS_Mesh()
2044 list<SMDS_Mesh*>::iterator itc=myChildren.begin();
2045 while(itc!=myChildren.end())
2053 delete myNodeIDFactory;
2054 delete myElementIDFactory;
2058 SMDS_ElemIteratorPtr eIt = elementsIterator();
2059 while ( eIt->more() )
2060 myElementIDFactory->ReleaseID(eIt->next()->GetID());
2061 SMDS_NodeIteratorPtr itn = nodesIterator();
2063 myNodeIDFactory->ReleaseID(itn->next()->GetID());
2066 // SetOfNodes::Iterator itn(myNodes);
2067 // for (; itn.More(); itn.Next())
2068 // delete itn.Value();
2070 // SetOf0DElements::Iterator it0d (my0DElements);
2071 // for (; it0d.More(); it0d.Next())
2073 // SMDS_MeshElement* elem = it0d.Value();
2077 // SetOfEdges::Iterator ite(myEdges);
2078 // for (; ite.More(); ite.Next())
2080 // SMDS_MeshElement* elem = ite.Value();
2084 // SetOfFaces::Iterator itf(myFaces);
2085 // for (; itf.More(); itf.Next())
2087 // SMDS_MeshElement* elem = itf.Value();
2091 // SetOfVolumes::Iterator itv(myVolumes);
2092 // for (; itv.More(); itv.Next())
2094 // SMDS_MeshElement* elem = itv.Value();
2099 //================================================================================
2101 * \brief Clear all data
2103 //================================================================================
2105 void SMDS_Mesh::Clear()
2107 if (myParent!=NULL) {
2108 SMDS_ElemIteratorPtr eIt = elementsIterator();
2109 while ( eIt->more() )
2110 myElementIDFactory->ReleaseID(eIt->next()->GetID());
2111 SMDS_NodeIteratorPtr itn = nodesIterator();
2113 myNodeIDFactory->ReleaseID(itn->next()->GetID());
2116 myNodeIDFactory->Clear();
2117 myElementIDFactory->Clear();
2120 SMDS_ElemIteratorPtr itv = elementsIterator();
2125 // SMDS_VolumeIteratorPtr itv = volumesIterator();
2126 // while (itv->more())
2127 // delete itv->next();
2128 // myVolumes.Clear();
2130 // SMDS_FaceIteratorPtr itf = facesIterator();
2131 // while (itf->more())
2132 // delete itf->next();
2135 // SMDS_EdgeIteratorPtr ite = edgesIterator();
2136 // while (ite->more())
2137 // delete ite->next();
2140 // SMDS_0DElementIteratorPtr it0d = elements0dIterator();
2141 // while (it0d->more())
2142 // delete it0d->next();
2143 // my0DElements.Clear();
2145 SMDS_NodeIteratorPtr itn = nodesIterator();
2150 list<SMDS_Mesh*>::iterator itc=myChildren.begin();
2151 while(itc!=myChildren.end())
2157 ///////////////////////////////////////////////////////////////////////////////
2158 /// Return true if this mesh create faces with edges.
2159 /// A false returned value mean that faces are created with nodes. A concequence
2160 /// is, iteration on edges (SMDS_Element::edgesIterator) will be unavailable.
2161 ///////////////////////////////////////////////////////////////////////////////
2162 bool SMDS_Mesh::hasConstructionEdges()
2164 return myHasConstructionEdges;
2167 ///////////////////////////////////////////////////////////////////////////////
2168 /// Return true if this mesh create volumes with faces
2169 /// A false returned value mean that volumes are created with nodes or edges.
2170 /// (see hasConstructionEdges)
2171 /// A concequence is, iteration on faces (SMDS_Element::facesIterator) will be
2173 ///////////////////////////////////////////////////////////////////////////////
2174 bool SMDS_Mesh::hasConstructionFaces()
2176 return myHasConstructionFaces;
2179 ///////////////////////////////////////////////////////////////////////////////
2180 /// Return true if nodes are linked to the finit elements, they are belonging to.
2181 /// Currently, It always return true.
2182 ///////////////////////////////////////////////////////////////////////////////
2183 bool SMDS_Mesh::hasInverseElements()
2185 return myHasInverseElements;
2188 ///////////////////////////////////////////////////////////////////////////////
2189 /// Make this mesh creating construction edges (see hasConstructionEdges)
2190 /// @param b true to have construction edges, else false.
2191 ///////////////////////////////////////////////////////////////////////////////
2192 void SMDS_Mesh::setConstructionEdges(bool b)
2194 myHasConstructionEdges=b;
2197 ///////////////////////////////////////////////////////////////////////////////
2198 /// Make this mesh creating construction faces (see hasConstructionFaces)
2199 /// @param b true to have construction faces, else false.
2200 ///////////////////////////////////////////////////////////////////////////////
2201 void SMDS_Mesh::setConstructionFaces(bool b)
2203 myHasConstructionFaces=b;
2206 ///////////////////////////////////////////////////////////////////////////////
2207 /// Make this mesh creating link from nodes to elements (see hasInverseElements)
2208 /// @param b true to link nodes to elements, else false.
2209 ///////////////////////////////////////////////////////////////////////////////
2210 void SMDS_Mesh::setInverseElements(bool b)
2212 if(!b) MESSAGE("Error : inverseElement=false not implemented");
2213 myHasInverseElements=b;
2216 ///////////////////////////////////////////////////////////////////////////////
2217 ///Iterator on NCollection_Map
2218 ///////////////////////////////////////////////////////////////////////////////
2219 template <class MAP, typename ELEM=const SMDS_MeshElement*, class FATHER=SMDS_ElemIterator>
2220 struct MYNode_Map_Iterator: public FATHER
2224 MYNode_Map_Iterator(const MAP& map): _map(map) // map is a std::vector<ELEM>
2231 while (_ctr < _map.size())
2242 ELEM current = _map[_ctr];
2248 template <class MAP, typename ELEM=const SMDS_MeshElement*, class FATHER=SMDS_ElemIterator>
2249 struct MYElem_Map_Iterator: public FATHER
2254 MYElem_Map_Iterator(const MAP& map, int typ): _map(map) // map is a std::vector<ELEM>
2262 while (_ctr < _map.size())
2265 if ( (_type == SMDSAbs_All) || (_map[_ctr]->GetType() == _type))
2274 ELEM current = dynamic_cast<ELEM> (_map[_ctr]);
2280 ///////////////////////////////////////////////////////////////////////////////
2281 /// Return an iterator on nodes of the current mesh factory
2282 ///////////////////////////////////////////////////////////////////////////////
2284 SMDS_NodeIteratorPtr SMDS_Mesh::nodesIterator() const
2286 //return SMDS_NodeIteratorPtr
2287 // (new SMDS_Mesh_MyNodeIterator(myNodeIDFactory->elementsIterator()));
2288 typedef MYNode_Map_Iterator
2289 < SetOfNodes, const SMDS_MeshNode*, SMDS_NodeIterator > TIterator;
2290 return SMDS_NodeIteratorPtr(new TIterator(myNodes));
2293 ///////////////////////////////////////////////////////////////////////////////
2294 ///Return an iterator on 0D elements of the current mesh.
2295 ///////////////////////////////////////////////////////////////////////////////
2297 SMDS_0DElementIteratorPtr SMDS_Mesh::elements0dIterator() const
2299 typedef MYElem_Map_Iterator
2300 < SetOfCells, const SMDS_Mesh0DElement*, SMDS_0DElementIterator > TIterator;
2301 return SMDS_0DElementIteratorPtr(new TIterator(myCells, SMDSAbs_0DElement));
2304 ///////////////////////////////////////////////////////////////////////////////
2305 ///Return an iterator on edges of the current mesh.
2306 ///////////////////////////////////////////////////////////////////////////////
2308 SMDS_EdgeIteratorPtr SMDS_Mesh::edgesIterator() const
2310 typedef MYElem_Map_Iterator
2311 < SetOfCells, const SMDS_MeshEdge*, SMDS_EdgeIterator > TIterator;
2312 return SMDS_EdgeIteratorPtr(new TIterator(myCells, SMDSAbs_Edge));
2315 ///////////////////////////////////////////////////////////////////////////////
2316 ///Return an iterator on faces of the current mesh.
2317 ///////////////////////////////////////////////////////////////////////////////
2319 SMDS_FaceIteratorPtr SMDS_Mesh::facesIterator() const
2321 typedef MYElem_Map_Iterator
2322 < SetOfCells, const SMDS_MeshFace*, SMDS_FaceIterator > TIterator;
2323 return SMDS_FaceIteratorPtr(new TIterator(myCells, SMDSAbs_Face));
2326 ///////////////////////////////////////////////////////////////////////////////
2327 ///Return an iterator on volumes of the current mesh.
2328 ///////////////////////////////////////////////////////////////////////////////
2330 SMDS_VolumeIteratorPtr SMDS_Mesh::volumesIterator() const
2332 typedef MYElem_Map_Iterator
2333 < SetOfCells, const SMDS_MeshVolume*, SMDS_VolumeIterator > TIterator;
2334 return SMDS_VolumeIteratorPtr(new TIterator(myCells, SMDSAbs_Volume));
2337 ///////////////////////////////////////////////////////////////////////////////
2338 /// Return an iterator on elements of the current mesh factory
2339 ///////////////////////////////////////////////////////////////////////////////
2340 SMDS_ElemIteratorPtr SMDS_Mesh::elementsIterator(SMDSAbs_ElementType type) const
2344 return SMDS_ElemIteratorPtr (new MYElem_Map_Iterator< SetOfCells >(myCells, SMDSAbs_All));
2346 case SMDSAbs_Volume:
2347 return SMDS_ElemIteratorPtr (new MYElem_Map_Iterator< SetOfCells >(myCells, SMDSAbs_Volume));
2349 return SMDS_ElemIteratorPtr (new MYElem_Map_Iterator< SetOfCells >(myCells, SMDSAbs_Face));
2351 return SMDS_ElemIteratorPtr (new MYElem_Map_Iterator< SetOfCells >(myCells, SMDSAbs_Edge));
2352 case SMDSAbs_0DElement:
2353 return SMDS_ElemIteratorPtr (new MYElem_Map_Iterator< SetOfCells >(myCells, SMDSAbs_0DElement));
2355 return SMDS_ElemIteratorPtr (new MYElem_Map_Iterator< SetOfNodes >(myNodes, SMDSAbs_All));
2356 //return myNodeIDFactory->elementsIterator();
2359 return myElementIDFactory->elementsIterator();
2362 ///////////////////////////////////////////////////////////////////////////////
2363 /// Do intersection of sets (more than 2)
2364 ///////////////////////////////////////////////////////////////////////////////
2365 static set<const SMDS_MeshElement*> * intersectionOfSets(
2366 set<const SMDS_MeshElement*> vs[], int numberOfSets)
2368 set<const SMDS_MeshElement*>* rsetA=new set<const SMDS_MeshElement*>(vs[0]);
2369 set<const SMDS_MeshElement*>* rsetB;
2371 for(int i=0; i<numberOfSets-1; i++)
2373 rsetB=new set<const SMDS_MeshElement*>();
2375 rsetA->begin(), rsetA->end(),
2376 vs[i+1].begin(), vs[i+1].end(),
2377 inserter(*rsetB, rsetB->begin()));
2384 ///////////////////////////////////////////////////////////////////////////////
2385 /// Return the list of finit elements owning the given element
2386 ///////////////////////////////////////////////////////////////////////////////
2387 static set<const SMDS_MeshElement*> * getFinitElements(const SMDS_MeshElement * element)
2389 int numberOfSets=element->NbNodes();
2390 set<const SMDS_MeshElement*> *initSet = new set<const SMDS_MeshElement*>[numberOfSets];
2392 SMDS_ElemIteratorPtr itNodes=element->nodesIterator();
2395 while(itNodes->more())
2397 const SMDS_MeshNode * n=static_cast<const SMDS_MeshNode*>(itNodes->next());
2398 SMDS_ElemIteratorPtr itFe = n->GetInverseElementIterator();
2400 //initSet[i]=set<const SMDS_MeshElement*>();
2402 initSet[i].insert(itFe->next());
2406 set<const SMDS_MeshElement*> *retSet=intersectionOfSets(initSet, numberOfSets);
2411 ///////////////////////////////////////////////////////////////////////////////
2412 /// Return the list of nodes used only by the given elements
2413 ///////////////////////////////////////////////////////////////////////////////
2414 static set<const SMDS_MeshElement*> * getExclusiveNodes(
2415 set<const SMDS_MeshElement*>& elements)
2417 set<const SMDS_MeshElement*> * toReturn=new set<const SMDS_MeshElement*>();
2418 set<const SMDS_MeshElement*>::iterator itElements=elements.begin();
2420 while(itElements!=elements.end())
2422 SMDS_ElemIteratorPtr itNodes = (*itElements)->nodesIterator();
2425 while(itNodes->more())
2427 const SMDS_MeshNode * n=static_cast<const SMDS_MeshNode*>(itNodes->next());
2428 SMDS_ElemIteratorPtr itFe = n->GetInverseElementIterator();
2429 set<const SMDS_MeshElement*> s;
2431 s.insert(itFe->next());
2432 if(s==elements) toReturn->insert(n);
2438 ///////////////////////////////////////////////////////////////////////////////
2439 ///Find the children of an element that are made of given nodes
2440 ///@param setOfChildren The set in which matching children will be inserted
2441 ///@param element The element were to search matching children
2442 ///@param nodes The nodes that the children must have to be selected
2443 ///////////////////////////////////////////////////////////////////////////////
2444 void SMDS_Mesh::addChildrenWithNodes(set<const SMDS_MeshElement*>& setOfChildren,
2445 const SMDS_MeshElement * element,
2446 set<const SMDS_MeshElement*>& nodes)
2448 switch(element->GetType())
2451 MESSAGE("Internal Error: This should not happend");
2453 case SMDSAbs_0DElement:
2459 SMDS_ElemIteratorPtr itn=element->nodesIterator();
2462 const SMDS_MeshElement * e=itn->next();
2463 if(nodes.find(e)!=nodes.end())
2465 setOfChildren.insert(element);
2472 SMDS_ElemIteratorPtr itn=element->nodesIterator();
2475 const SMDS_MeshElement * e=itn->next();
2476 if(nodes.find(e)!=nodes.end())
2478 setOfChildren.insert(element);
2482 if(hasConstructionEdges())
2484 SMDS_ElemIteratorPtr ite=element->edgesIterator();
2486 addChildrenWithNodes(setOfChildren, ite->next(), nodes);
2489 case SMDSAbs_Volume:
2491 if(hasConstructionFaces())
2493 SMDS_ElemIteratorPtr ite=element->facesIterator();
2495 addChildrenWithNodes(setOfChildren, ite->next(), nodes);
2497 else if(hasConstructionEdges())
2499 SMDS_ElemIteratorPtr ite=element->edgesIterator();
2501 addChildrenWithNodes(setOfChildren, ite->next(), nodes);
2507 ///////////////////////////////////////////////////////////////////////////////
2508 ///@param elem The element to delete
2509 ///@param removenodes if true remaining nodes will be removed
2510 ///////////////////////////////////////////////////////////////////////////////
2511 void SMDS_Mesh::RemoveElement(const SMDS_MeshElement * elem,
2512 const bool removenodes)
2514 list<const SMDS_MeshElement *> removedElems;
2515 list<const SMDS_MeshElement *> removedNodes;
2516 RemoveElement( elem, removedElems, removedNodes, removenodes );
2519 ///////////////////////////////////////////////////////////////////////////////
2520 ///@param elem The element to delete
2521 ///@param removedElems contains all removed elements
2522 ///@param removedNodes contains all removed nodes
2523 ///@param removenodes if true remaining nodes will be removed
2524 ///////////////////////////////////////////////////////////////////////////////
2525 void SMDS_Mesh::RemoveElement(const SMDS_MeshElement * elem,
2526 list<const SMDS_MeshElement *>& removedElems,
2527 list<const SMDS_MeshElement *>& removedNodes,
2530 // get finite elements built on elem
2531 set<const SMDS_MeshElement*> * s1;
2532 if (elem->GetType() == SMDSAbs_0DElement ||
2533 elem->GetType() == SMDSAbs_Edge && !hasConstructionEdges() ||
2534 elem->GetType() == SMDSAbs_Face && !hasConstructionFaces() ||
2535 elem->GetType() == SMDSAbs_Volume)
2537 s1 = new set<const SMDS_MeshElement*>();
2541 s1 = getFinitElements(elem);
2543 // get exclusive nodes (which would become free afterwards)
2544 set<const SMDS_MeshElement*> * s2;
2545 if (elem->GetType() == SMDSAbs_Node) // a node is removed
2547 // do not remove nodes except elem
2548 s2 = new set<const SMDS_MeshElement*>();
2553 s2 = getExclusiveNodes(*s1);
2555 // form the set of finite and construction elements to remove
2556 set<const SMDS_MeshElement*> s3;
2557 set<const SMDS_MeshElement*>::iterator it=s1->begin();
2558 while(it!=s1->end())
2560 addChildrenWithNodes(s3, *it ,*s2);
2564 if(elem->GetType()!=SMDSAbs_Node) s3.insert(elem);
2566 // remove finite and construction elements
2570 // Remove element from <InverseElements> of its nodes
2571 SMDS_ElemIteratorPtr itn=(*it)->nodesIterator();
2574 SMDS_MeshNode * n = static_cast<SMDS_MeshNode *>
2575 (const_cast<SMDS_MeshElement *>(itn->next()));
2576 n->RemoveInverseElement( (*it) );
2579 switch((*it)->GetType())
2582 MESSAGE("Internal Error: This should not happen");
2584 case SMDSAbs_0DElement:
2585 myCells[(*it)->GetID()] = 0; // -PR- ici ou dans myElementIDFactory->ReleaseID ?
2589 myCells[(*it)->GetID()] = 0;
2590 myInfo.RemoveEdge(*it);
2593 myCells[(*it)->GetID()] = 0;
2594 myInfo.RemoveFace(*it);
2596 case SMDSAbs_Volume:
2597 myCells[(*it)->GetID()] = 0;
2598 myInfo.RemoveVolume(*it);
2601 //MESSAGE( "SMDS: RM elem " << (*it)->GetID() );
2602 removedElems.push_back( (*it) );
2603 myElementIDFactory->ReleaseID((*it)->GetID());
2608 // remove exclusive (free) nodes
2612 while(it!=s2->end())
2614 //MESSAGE( "SMDS: RM node " << (*it)->GetID() );
2615 myNodes[(*it)->GetID()] = 0;
2617 myNodeIDFactory->ReleaseID((*it)->GetID());
2618 removedNodes.push_back( (*it) );
2629 ///////////////////////////////////////////////////////////////////////////////
2630 ///@param elem The element to delete
2631 ///////////////////////////////////////////////////////////////////////////////
2632 void SMDS_Mesh::RemoveFreeElement(const SMDS_MeshElement * elem)
2634 SMDSAbs_ElementType aType = elem->GetType();
2635 if (aType == SMDSAbs_Node) {
2636 // only free node can be removed by this method
2637 const SMDS_MeshNode* n = static_cast<const SMDS_MeshNode*>(elem);
2638 SMDS_ElemIteratorPtr itFe = n->GetInverseElementIterator();
2639 if (!itFe->more()) { // free node
2640 myNodes[elem->GetID()] = 0;
2642 myNodeIDFactory->ReleaseID(elem->GetID());
2646 if (hasConstructionEdges() || hasConstructionFaces())
2647 // this methods is only for meshes without descendants
2650 // Remove element from <InverseElements> of its nodes
2651 SMDS_ElemIteratorPtr itn = elem->nodesIterator();
2652 while (itn->more()) {
2653 SMDS_MeshNode * n = static_cast<SMDS_MeshNode *>
2654 (const_cast<SMDS_MeshElement *>(itn->next()));
2655 n->RemoveInverseElement(elem);
2658 // in meshes without descendants elements are always free
2660 case SMDSAbs_0DElement:
2661 myCells[elem->GetID()] = 0;
2662 myInfo.remove(elem);
2665 myCells[elem->GetID()] = 0;
2666 myInfo.RemoveEdge(elem);
2669 myCells[elem->GetID()] = 0;
2670 myInfo.RemoveFace(elem);
2672 case SMDSAbs_Volume:
2673 myCells[elem->GetID()] = 0;
2674 myInfo.RemoveVolume(elem);
2679 myElementIDFactory->ReleaseID(elem->GetID());
2685 * Checks if the element is present in mesh.
2686 * Useful to determine dead pointers.
2688 bool SMDS_Mesh::Contains (const SMDS_MeshElement* elem) const
2690 // we should not imply on validity of *elem, so iterate on containers
2691 // of all types in the hope of finding <elem> somewhere there
2692 SMDS_NodeIteratorPtr itn = nodesIterator();
2694 if (elem == itn->next())
2696 SMDS_0DElementIteratorPtr it0d = elements0dIterator();
2697 while (it0d->more())
2698 if (elem == it0d->next())
2700 SMDS_EdgeIteratorPtr ite = edgesIterator();
2702 if (elem == ite->next())
2704 SMDS_FaceIteratorPtr itf = facesIterator();
2706 if (elem == itf->next())
2708 SMDS_VolumeIteratorPtr itv = volumesIterator();
2710 if (elem == itv->next())
2715 //=======================================================================
2716 //function : MaxNodeID
2718 //=======================================================================
2720 int SMDS_Mesh::MaxNodeID() const
2725 //=======================================================================
2726 //function : MinNodeID
2728 //=======================================================================
2730 int SMDS_Mesh::MinNodeID() const
2735 //=======================================================================
2736 //function : MaxElementID
2738 //=======================================================================
2740 int SMDS_Mesh::MaxElementID() const
2742 return myElementIDFactory->GetMaxID();
2745 //=======================================================================
2746 //function : MinElementID
2748 //=======================================================================
2750 int SMDS_Mesh::MinElementID() const
2752 return myElementIDFactory->GetMinID();
2755 //=======================================================================
2756 //function : Renumber
2757 //purpose : Renumber all nodes or elements.
2758 //=======================================================================
2760 void SMDS_Mesh::Renumber (const bool isNodes, const int startID, const int deltaID)
2762 MESSAGE("Renumber");
2766 SMDS_MeshNodeIDFactory * idFactory =
2767 isNodes ? myNodeIDFactory : myElementIDFactory;
2769 // get existing elements in the order of ID increasing
2770 map<int,SMDS_MeshElement*> elemMap;
2771 SMDS_ElemIteratorPtr idElemIt = idFactory->elementsIterator();
2772 while ( idElemIt->more() ) {
2773 SMDS_MeshElement* elem = const_cast<SMDS_MeshElement*>(idElemIt->next());
2774 int id = elem->GetID();
2775 elemMap.insert(map<int,SMDS_MeshElement*>::value_type(id, elem));
2777 // release their ids
2778 map<int,SMDS_MeshElement*>::iterator elemIt = elemMap.begin();
2780 // for ( ; elemIt != elemMap.end(); elemIt++ )
2782 // int id = (*elemIt).first;
2783 // idFactory->ReleaseID( id );
2787 elemIt = elemMap.begin();
2788 for ( ; elemIt != elemMap.end(); elemIt++ )
2790 idFactory->BindID( ID, (*elemIt).second );
2795 //=======================================================================
2796 //function : GetElementType
2797 //purpose : Return type of element or node with id
2798 //=======================================================================
2800 SMDSAbs_ElementType SMDS_Mesh::GetElementType( const int id, const bool iselem ) const
2802 SMDS_MeshElement* elem = 0;
2804 elem = myElementIDFactory->MeshElement( id );
2806 elem = myNodeIDFactory->MeshElement( id );
2810 //throw SALOME_Exception(LOCALIZED ("this element isn't exist"));
2814 return elem->GetType();
2819 //********************************************************************
2820 //********************************************************************
2821 //******** *********
2822 //***** Methods for addition of quadratic elements ******
2823 //******** *********
2824 //********************************************************************
2825 //********************************************************************
2827 //=======================================================================
2828 //function : AddEdgeWithID
2830 //=======================================================================
2831 SMDS_MeshEdge* SMDS_Mesh::AddEdgeWithID(int n1, int n2, int n12, int ID)
2833 return SMDS_Mesh::AddEdgeWithID
2834 ((SMDS_MeshNode*) myNodeIDFactory->MeshElement(n1),
2835 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n2),
2836 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n12),
2840 //=======================================================================
2841 //function : AddEdge
2843 //=======================================================================
2844 SMDS_MeshEdge* SMDS_Mesh::AddEdge(const SMDS_MeshNode* n1,
2845 const SMDS_MeshNode* n2,
2846 const SMDS_MeshNode* n12)
2848 return SMDS_Mesh::AddEdgeWithID(n1, n2, n12, myElementIDFactory->GetFreeID());
2851 //=======================================================================
2852 //function : AddEdgeWithID
2854 //=======================================================================
2855 SMDS_MeshEdge* SMDS_Mesh::AddEdgeWithID(const SMDS_MeshNode * n1,
2856 const SMDS_MeshNode * n2,
2857 const SMDS_MeshNode * n12,
2860 if ( !n1 || !n2 || !n12 ) return 0;
2861 SMDS_QuadraticEdge* edge = new SMDS_QuadraticEdge(n1,n2,n12);
2862 if(myElementIDFactory->BindID(ID, edge)) {
2863 SMDS_MeshNode *node1,*node2, *node12;
2864 //node1 = const_cast<SMDS_MeshNode*>(n1);
2865 //node2 = const_cast<SMDS_MeshNode*>(n2);
2866 //node12 = const_cast<SMDS_MeshNode*>(n12);
2867 //node1->AddInverseElement(edge); // --- fait avec BindID
2868 //node2->AddInverseElement(edge);
2869 //node12->AddInverseElement(edge);
2870 adjustmyCellsCapacity(ID);
2872 myInfo.myNbQuadEdges++;
2882 //=======================================================================
2883 //function : AddFace
2885 //=======================================================================
2886 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshNode * n1,
2887 const SMDS_MeshNode * n2,
2888 const SMDS_MeshNode * n3,
2889 const SMDS_MeshNode * n12,
2890 const SMDS_MeshNode * n23,
2891 const SMDS_MeshNode * n31)
2893 return SMDS_Mesh::AddFaceWithID(n1,n2,n3,n12,n23,n31,
2894 myElementIDFactory->GetFreeID());
2897 //=======================================================================
2898 //function : AddFaceWithID
2900 //=======================================================================
2901 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(int n1, int n2, int n3,
2902 int n12,int n23,int n31, int ID)
2904 return SMDS_Mesh::AddFaceWithID
2905 ((SMDS_MeshNode *)myNodeIDFactory->MeshElement(n1) ,
2906 (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n2) ,
2907 (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n3) ,
2908 (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n12),
2909 (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n23),
2910 (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n31),
2914 //=======================================================================
2915 //function : AddFaceWithID
2917 //=======================================================================
2918 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshNode * n1,
2919 const SMDS_MeshNode * n2,
2920 const SMDS_MeshNode * n3,
2921 const SMDS_MeshNode * n12,
2922 const SMDS_MeshNode * n23,
2923 const SMDS_MeshNode * n31,
2926 if ( !n1 || !n2 || !n3 || !n12 || !n23 || !n31) return 0;
2927 if(hasConstructionEdges()) {
2928 // creation quadratic edges - not implemented
2931 SMDS_QuadraticFaceOfNodes* face =
2932 new SMDS_QuadraticFaceOfNodes(n1,n2,n3,n12,n23,n31);
2933 adjustmyCellsCapacity(ID);
2935 myInfo.myNbQuadTriangles++;
2937 if (!registerElement(ID, face)) {
2938 RemoveElement(face, false);
2945 //=======================================================================
2946 //function : AddFace
2948 //=======================================================================
2949 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshNode * n1,
2950 const SMDS_MeshNode * n2,
2951 const SMDS_MeshNode * n3,
2952 const SMDS_MeshNode * n4,
2953 const SMDS_MeshNode * n12,
2954 const SMDS_MeshNode * n23,
2955 const SMDS_MeshNode * n34,
2956 const SMDS_MeshNode * n41)
2958 return SMDS_Mesh::AddFaceWithID(n1,n2,n3,n4,n12,n23,n34,n41,
2959 myElementIDFactory->GetFreeID());
2962 //=======================================================================
2963 //function : AddFaceWithID
2965 //=======================================================================
2966 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(int n1, int n2, int n3, int n4,
2967 int n12,int n23,int n34,int n41, int ID)
2969 return SMDS_Mesh::AddFaceWithID
2970 ((SMDS_MeshNode *)myNodeIDFactory->MeshElement(n1) ,
2971 (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n2) ,
2972 (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n3) ,
2973 (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n4) ,
2974 (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n12),
2975 (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n23),
2976 (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n34),
2977 (SMDS_MeshNode *)myNodeIDFactory->MeshElement(n41),
2981 //=======================================================================
2982 //function : AddFaceWithID
2984 //=======================================================================
2985 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshNode * n1,
2986 const SMDS_MeshNode * n2,
2987 const SMDS_MeshNode * n3,
2988 const SMDS_MeshNode * n4,
2989 const SMDS_MeshNode * n12,
2990 const SMDS_MeshNode * n23,
2991 const SMDS_MeshNode * n34,
2992 const SMDS_MeshNode * n41,
2995 if ( !n1 || !n2 || !n3 || !n4 || !n12 || !n23 || !n34 || !n41) return 0;
2996 if(hasConstructionEdges()) {
2997 // creation quadratic edges - not implemented
2999 SMDS_QuadraticFaceOfNodes* face =
3000 new SMDS_QuadraticFaceOfNodes(n1,n2,n3,n4,n12,n23,n34,n41);
3001 adjustmyCellsCapacity(ID);
3003 myInfo.myNbQuadQuadrangles++;
3005 if (!registerElement(ID, face)) {
3006 RemoveElement(face, false);
3013 //=======================================================================
3014 //function : AddVolume
3016 //=======================================================================
3017 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
3018 const SMDS_MeshNode * n2,
3019 const SMDS_MeshNode * n3,
3020 const SMDS_MeshNode * n4,
3021 const SMDS_MeshNode * n12,
3022 const SMDS_MeshNode * n23,
3023 const SMDS_MeshNode * n31,
3024 const SMDS_MeshNode * n14,
3025 const SMDS_MeshNode * n24,
3026 const SMDS_MeshNode * n34)
3028 int ID = myElementIDFactory->GetFreeID();
3029 SMDS_MeshVolume * v = SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n12, n23,
3030 n31, n14, n24, n34, ID);
3031 if(v==NULL) myElementIDFactory->ReleaseID(ID);
3035 //=======================================================================
3036 //function : AddVolumeWithID
3038 //=======================================================================
3039 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3, int n4,
3040 int n12,int n23,int n31,
3041 int n14,int n24,int n34, int ID)
3043 return SMDS_Mesh::AddVolumeWithID
3044 ((SMDS_MeshNode*) myNodeIDFactory->MeshElement(n1) ,
3045 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n2) ,
3046 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n3) ,
3047 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n4) ,
3048 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n12),
3049 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n23),
3050 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n31),
3051 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n14),
3052 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n24),
3053 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n34),
3057 //=======================================================================
3058 //function : AddVolumeWithID
3059 //purpose : 2d order tetrahedron of 10 nodes
3060 //=======================================================================
3061 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
3062 const SMDS_MeshNode * n2,
3063 const SMDS_MeshNode * n3,
3064 const SMDS_MeshNode * n4,
3065 const SMDS_MeshNode * n12,
3066 const SMDS_MeshNode * n23,
3067 const SMDS_MeshNode * n31,
3068 const SMDS_MeshNode * n14,
3069 const SMDS_MeshNode * n24,
3070 const SMDS_MeshNode * n34,
3073 if ( !n1 || !n2 || !n3 || !n4 || !n12 || !n23 || !n31 || !n14 || !n24 || !n34)
3075 if(hasConstructionFaces()) {
3076 // creation quadratic faces - not implemented
3079 SMDS_QuadraticVolumeOfNodes * volume =
3080 new SMDS_QuadraticVolumeOfNodes(n1,n2,n3,n4,n12,n23,n31,n14,n24,n34);
3081 adjustmyCellsCapacity(ID);
3082 myCells[ID] = volume;
3083 myInfo.myNbQuadTetras++;
3085 if (!registerElement(ID, volume)) {
3086 RemoveElement(volume, false);
3093 //=======================================================================
3094 //function : AddVolume
3096 //=======================================================================
3097 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
3098 const SMDS_MeshNode * n2,
3099 const SMDS_MeshNode * n3,
3100 const SMDS_MeshNode * n4,
3101 const SMDS_MeshNode * n5,
3102 const SMDS_MeshNode * n12,
3103 const SMDS_MeshNode * n23,
3104 const SMDS_MeshNode * n34,
3105 const SMDS_MeshNode * n41,
3106 const SMDS_MeshNode * n15,
3107 const SMDS_MeshNode * n25,
3108 const SMDS_MeshNode * n35,
3109 const SMDS_MeshNode * n45)
3111 int ID = myElementIDFactory->GetFreeID();
3112 SMDS_MeshVolume * v =
3113 SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n12, n23, n34, n41,
3114 n15, n25, n35, n45, ID);
3115 if(v==NULL) myElementIDFactory->ReleaseID(ID);
3119 //=======================================================================
3120 //function : AddVolumeWithID
3122 //=======================================================================
3123 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3, int n4, int n5,
3124 int n12,int n23,int n34,int n41,
3125 int n15,int n25,int n35,int n45, int ID)
3127 return SMDS_Mesh::AddVolumeWithID
3128 ((SMDS_MeshNode*) myNodeIDFactory->MeshElement(n1) ,
3129 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n2) ,
3130 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n3) ,
3131 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n4) ,
3132 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n5) ,
3133 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n12),
3134 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n23),
3135 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n34),
3136 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n41),
3137 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n15),
3138 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n25),
3139 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n35),
3140 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n45),
3144 //=======================================================================
3145 //function : AddVolumeWithID
3146 //purpose : 2d order pyramid of 13 nodes
3147 //=======================================================================
3148 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
3149 const SMDS_MeshNode * n2,
3150 const SMDS_MeshNode * n3,
3151 const SMDS_MeshNode * n4,
3152 const SMDS_MeshNode * n5,
3153 const SMDS_MeshNode * n12,
3154 const SMDS_MeshNode * n23,
3155 const SMDS_MeshNode * n34,
3156 const SMDS_MeshNode * n41,
3157 const SMDS_MeshNode * n15,
3158 const SMDS_MeshNode * n25,
3159 const SMDS_MeshNode * n35,
3160 const SMDS_MeshNode * n45,
3163 if (!n1 || !n2 || !n3 || !n4 || !n5 || !n12 || !n23 ||
3164 !n34 || !n41 || !n15 || !n25 || !n35 || !n45)
3166 if(hasConstructionFaces()) {
3167 // creation quadratic faces - not implemented
3170 SMDS_QuadraticVolumeOfNodes * volume =
3171 new SMDS_QuadraticVolumeOfNodes(n1,n2,n3,n4,n5,n12,n23,
3172 n34,n41,n15,n25,n35,n45);
3173 adjustmyCellsCapacity(ID);
3174 myCells[ID] = volume;
3175 myInfo.myNbQuadPyramids++;
3177 if (!registerElement(ID, volume)) {
3178 RemoveElement(volume, false);
3185 //=======================================================================
3186 //function : AddVolume
3188 //=======================================================================
3189 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
3190 const SMDS_MeshNode * n2,
3191 const SMDS_MeshNode * n3,
3192 const SMDS_MeshNode * n4,
3193 const SMDS_MeshNode * n5,
3194 const SMDS_MeshNode * n6,
3195 const SMDS_MeshNode * n12,
3196 const SMDS_MeshNode * n23,
3197 const SMDS_MeshNode * n31,
3198 const SMDS_MeshNode * n45,
3199 const SMDS_MeshNode * n56,
3200 const SMDS_MeshNode * n64,
3201 const SMDS_MeshNode * n14,
3202 const SMDS_MeshNode * n25,
3203 const SMDS_MeshNode * n36)
3205 int ID = myElementIDFactory->GetFreeID();
3206 SMDS_MeshVolume * v =
3207 SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n6, n12, n23, n31,
3208 n45, n56, n64, n14, n25, n36, ID);
3209 if(v==NULL) myElementIDFactory->ReleaseID(ID);
3213 //=======================================================================
3214 //function : AddVolumeWithID
3216 //=======================================================================
3217 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3,
3218 int n4, int n5, int n6,
3219 int n12,int n23,int n31,
3220 int n45,int n56,int n64,
3221 int n14,int n25,int n36, int ID)
3223 return SMDS_Mesh::AddVolumeWithID
3224 ((SMDS_MeshNode*) myNodeIDFactory->MeshElement(n1) ,
3225 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n2) ,
3226 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n3) ,
3227 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n4) ,
3228 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n5) ,
3229 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n6) ,
3230 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n12),
3231 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n23),
3232 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n31),
3233 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n45),
3234 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n56),
3235 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n64),
3236 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n14),
3237 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n25),
3238 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n36),
3242 //=======================================================================
3243 //function : AddVolumeWithID
3244 //purpose : 2d order Pentahedron with 15 nodes
3245 //=======================================================================
3246 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
3247 const SMDS_MeshNode * n2,
3248 const SMDS_MeshNode * n3,
3249 const SMDS_MeshNode * n4,
3250 const SMDS_MeshNode * n5,
3251 const SMDS_MeshNode * n6,
3252 const SMDS_MeshNode * n12,
3253 const SMDS_MeshNode * n23,
3254 const SMDS_MeshNode * n31,
3255 const SMDS_MeshNode * n45,
3256 const SMDS_MeshNode * n56,
3257 const SMDS_MeshNode * n64,
3258 const SMDS_MeshNode * n14,
3259 const SMDS_MeshNode * n25,
3260 const SMDS_MeshNode * n36,
3263 if (!n1 || !n2 || !n3 || !n4 || !n5 || !n6 || !n12 || !n23 ||
3264 !n31 || !n45 || !n56 || !n64 || !n14 || !n25 || !n36)
3266 if(hasConstructionFaces()) {
3267 // creation quadratic faces - not implemented
3270 SMDS_QuadraticVolumeOfNodes * volume =
3271 new SMDS_QuadraticVolumeOfNodes(n1,n2,n3,n4,n5,n6,n12,n23,n31,
3272 n45,n56,n64,n14,n25,n36);
3273 adjustmyCellsCapacity(ID);
3274 myCells[ID] = volume;
3275 myInfo.myNbQuadPrisms++;
3277 if (!registerElement(ID, volume)) {
3278 RemoveElement(volume, false);
3285 //=======================================================================
3286 //function : AddVolume
3288 //=======================================================================
3289 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
3290 const SMDS_MeshNode * n2,
3291 const SMDS_MeshNode * n3,
3292 const SMDS_MeshNode * n4,
3293 const SMDS_MeshNode * n5,
3294 const SMDS_MeshNode * n6,
3295 const SMDS_MeshNode * n7,
3296 const SMDS_MeshNode * n8,
3297 const SMDS_MeshNode * n12,
3298 const SMDS_MeshNode * n23,
3299 const SMDS_MeshNode * n34,
3300 const SMDS_MeshNode * n41,
3301 const SMDS_MeshNode * n56,
3302 const SMDS_MeshNode * n67,
3303 const SMDS_MeshNode * n78,
3304 const SMDS_MeshNode * n85,
3305 const SMDS_MeshNode * n15,
3306 const SMDS_MeshNode * n26,
3307 const SMDS_MeshNode * n37,
3308 const SMDS_MeshNode * n48)
3310 int ID = myElementIDFactory->GetFreeID();
3311 SMDS_MeshVolume * v =
3312 SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n6, n7, n8, n12, n23, n34, n41,
3313 n56, n67, n78, n85, n15, n26, n37, n48, ID);
3314 if(v==NULL) myElementIDFactory->ReleaseID(ID);
3318 //=======================================================================
3319 //function : AddVolumeWithID
3321 //=======================================================================
3322 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3, int n4,
3323 int n5, int n6, int n7, int n8,
3324 int n12,int n23,int n34,int n41,
3325 int n56,int n67,int n78,int n85,
3326 int n15,int n26,int n37,int n48, int ID)
3328 return SMDS_Mesh::AddVolumeWithID
3329 ((SMDS_MeshNode*) myNodeIDFactory->MeshElement(n1),
3330 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n2),
3331 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n3),
3332 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n4),
3333 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n5),
3334 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n6),
3335 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n7),
3336 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n8),
3337 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n12),
3338 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n23),
3339 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n34),
3340 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n41),
3341 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n56),
3342 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n67),
3343 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n78),
3344 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n85),
3345 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n15),
3346 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n26),
3347 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n37),
3348 (SMDS_MeshNode*) myNodeIDFactory->MeshElement(n48),
3352 //=======================================================================
3353 //function : AddVolumeWithID
3354 //purpose : 2d order Hexahedrons with 20 nodes
3355 //=======================================================================
3356 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
3357 const SMDS_MeshNode * n2,
3358 const SMDS_MeshNode * n3,
3359 const SMDS_MeshNode * n4,
3360 const SMDS_MeshNode * n5,
3361 const SMDS_MeshNode * n6,
3362 const SMDS_MeshNode * n7,
3363 const SMDS_MeshNode * n8,
3364 const SMDS_MeshNode * n12,
3365 const SMDS_MeshNode * n23,
3366 const SMDS_MeshNode * n34,
3367 const SMDS_MeshNode * n41,
3368 const SMDS_MeshNode * n56,
3369 const SMDS_MeshNode * n67,
3370 const SMDS_MeshNode * n78,
3371 const SMDS_MeshNode * n85,
3372 const SMDS_MeshNode * n15,
3373 const SMDS_MeshNode * n26,
3374 const SMDS_MeshNode * n37,
3375 const SMDS_MeshNode * n48,
3378 if (!n1 || !n2 || !n3 || !n4 || !n5 || !n6 || !n7 || !n8 || !n12 || !n23 ||
3379 !n34 || !n41 || !n56 || !n67 || !n78 || !n85 || !n15 || !n26 || !n37 || !n48)
3381 if(hasConstructionFaces()) {
3383 // creation quadratic faces - not implemented
3385 SMDS_QuadraticVolumeOfNodes * volume =
3386 new SMDS_QuadraticVolumeOfNodes(n1,n2,n3,n4,n5,n6,n7,n8,n12,n23,n34,n41,
3387 n56,n67,n78,n85,n15,n26,n37,n48);
3388 adjustmyCellsCapacity(ID);
3389 myCells[ID] = volume;
3390 myInfo.myNbQuadHexas++;
3392 if (!registerElement(ID, volume)) {
3393 RemoveElement(volume, false);
3399 void SMDS_Mesh::updateNodeMinMax()
3402 while (!myNodes[myNodeMin] && (myNodeMin<myNodes.size()))
3404 myNodeMax=myNodes.size()-1;
3405 while (!myNodes[myNodeMax] && (myNodeMin>=0))