1 // Copyright (C) 2007-2020 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, or (at your option) any later version.
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
23 // SMESH SMDS : implementation of Salome mesh data structure
26 #pragma warning(disable:4786)
29 #include "SMDS_Mesh.hxx"
31 #include "SMDS_ElementFactory.hxx"
32 #include "SMDS_ElementHolder.hxx"
33 #include "SMDS_SetIterator.hxx"
34 #include "SMDS_SpacePosition.hxx"
35 #include "SMDS_UnstructuredGrid.hxx"
37 #include <utilities.h>
39 #include <vtkUnstructuredGrid.h>
40 //#include <vtkUnstructuredGridWriter.h>
42 #include <vtkUnsignedCharArray.h>
43 #include <vtkCellLinks.h>
44 #include <vtkIdList.h>
50 #include <boost/make_shared.hpp>
52 #if !defined WIN32 && !defined __APPLE__
53 #include <sys/sysinfo.h>
56 // number of added entities to check memory after
57 #define CHECKMEMORY_INTERVAL 100000
59 #define MYASSERT(val) if (!(val)) throw SALOME_Exception(LOCALIZED("assertion not verified"));
61 int SMDS_Mesh::chunkSize = 1024;
63 //================================================================================
65 * \brief Raise an exception if free memory (ram+swap) too low
66 * \param doNotRaise - if true, suppress exception, just return free memory size
67 * \retval int - amount of available memory in MB or negative number in failure case
69 //================================================================================
71 int SMDS_Mesh::CheckMemory(const bool doNotRaise)
74 #if !defined WIN32 && !defined __APPLE__
76 int err = sysinfo( &si );
80 const unsigned long Mbyte = 1024 * 1024;
82 static int limit = -1;
84 if ( si.totalswap == 0 )
86 int status = system("SMDS_MemoryLimit"); // it returns lower limit of free RAM
88 limit = WEXITSTATUS(status);
91 double factor = ( si.totalswap == 0 ) ? 0.1 : 0.2;
92 limit = int(( factor * si.totalram * si.mem_unit ) / Mbyte );
98 limit = int ( limit * 1.5 );
99 MESSAGE ( "SMDS_Mesh::CheckMemory() memory limit = " << limit << " MB" );
102 // compute separately to avoid overflow
104 ( si.freeram * si.mem_unit ) / Mbyte +
105 ( si.freeswap * si.mem_unit ) / Mbyte;
107 if ( freeMb > limit )
108 return freeMb - limit;
113 MESSAGE ("SMDS_Mesh::CheckMemory() throws as free memory too low: " << freeMb <<" MB" );
114 throw std::bad_alloc();
120 ///////////////////////////////////////////////////////////////////////////////
121 /// Create a new mesh object
122 ///////////////////////////////////////////////////////////////////////////////
123 SMDS_Mesh::SMDS_Mesh():
124 myNodeFactory( new SMDS_NodeFactory( this )),
125 myCellFactory( new SMDS_ElementFactory( this )),
127 myModified(false), myModifTime(0), myCompactTime(0),
128 xmin(0), xmax(0), ymin(0), ymax(0), zmin(0), zmax(0)
130 myGrid = SMDS_UnstructuredGrid::New();
131 myGrid->setSMDS_mesh(this);
132 myGrid->Initialize();
134 vtkPoints* points = vtkPoints::New();
135 // bug "21125: EDF 1233 SMESH: Degrardation of precision in a test case for quadratic conversion"
136 // Use double type for storing coordinates of nodes instead of float.
137 points->SetDataType(VTK_DOUBLE);
138 points->SetNumberOfPoints( 0 );
139 myGrid->SetPoints( points );
143 // initialize static maps in SMDS_MeshCell, to be thread-safe
144 SMDS_MeshCell::InitStaticMembers();
147 ///////////////////////////////////////////////////////////////////////////////
148 /// Create a new child mesh
149 /// Note that the tree structure of SMDS_Mesh seems to be unused in this version
150 /// (2003-09-08) of SMESH
151 ///////////////////////////////////////////////////////////////////////////////
152 SMDS_Mesh::SMDS_Mesh(SMDS_Mesh * parent):
153 myNodeFactory( new SMDS_NodeFactory( this )),
154 myCellFactory( new SMDS_ElementFactory( this )),
159 ///////////////////////////////////////////////////////////////////////////////
160 ///Create a submesh and add it to the current mesh
161 ///////////////////////////////////////////////////////////////////////////////
163 SMDS_Mesh *SMDS_Mesh::AddSubMesh()
165 SMDS_Mesh *submesh = new SMDS_Mesh(this);
166 myChildren.insert(myChildren.end(), submesh);
170 ///////////////////////////////////////////////////////////////////////////////
171 ///create a MeshNode and add it to the current Mesh
172 ///An ID is automatically assigned to the node.
173 ///@return : The created node
174 ///////////////////////////////////////////////////////////////////////////////
176 SMDS_MeshNode * SMDS_Mesh::AddNode(double x, double y, double z)
178 return SMDS_Mesh::AddNodeWithID( x,y,z, myNodeFactory->GetFreeID() );
181 ///////////////////////////////////////////////////////////////////////////////
182 ///create a MeshNode and add it to the current Mesh
183 ///@param ID : The ID of the MeshNode to create
184 ///@return : The created node or NULL if a node with this ID already exists
185 ///////////////////////////////////////////////////////////////////////////////
186 SMDS_MeshNode * SMDS_Mesh::AddNodeWithID( double x, double y, double z, int ID )
188 // find the MeshNode corresponding to ID
189 SMDS_MeshNode *node = myNodeFactory->NewNode( ID );
192 node->init( x, y, z );
195 this->adjustBoundingBox(x, y, z);
200 ///////////////////////////////////////////////////////////////////////////////
201 /// create a Mesh0DElement and add it to the current Mesh
202 /// @return : The created Mesh0DElement
203 ///////////////////////////////////////////////////////////////////////////////
204 SMDS_Mesh0DElement* SMDS_Mesh::Add0DElementWithID(int idnode, int ID)
206 const SMDS_MeshNode * node = myNodeFactory->FindNode(idnode);
207 if (!node) return NULL;
208 return SMDS_Mesh::Add0DElementWithID(node, ID);
211 ///////////////////////////////////////////////////////////////////////////////
212 /// create a Mesh0DElement and add it to the current Mesh
213 /// @return : The created Mesh0DElement
214 ///////////////////////////////////////////////////////////////////////////////
215 SMDS_Mesh0DElement* SMDS_Mesh::Add0DElement(const SMDS_MeshNode * node)
217 return SMDS_Mesh::Add0DElementWithID( node, myCellFactory->GetFreeID() );
220 ///////////////////////////////////////////////////////////////////////////////
221 /// Create a new Mesh0DElement and at it to the mesh
222 /// @param idnode ID of the node
223 /// @param ID ID of the 0D element to create
224 /// @return The created 0D element or NULL if an element with this
225 /// ID already exists or if input node is not found.
226 ///////////////////////////////////////////////////////////////////////////////
227 SMDS_Mesh0DElement* SMDS_Mesh::Add0DElementWithID(const SMDS_MeshNode * n, int ID)
231 if (Nb0DElements() % CHECKMEMORY_INTERVAL == 0) CheckMemory();
233 if ( SMDS_MeshCell * cell = myCellFactory->NewCell( ID ))
235 cell->init( SMDSEntity_0D, /*nbNodes=*/1, n );
236 myInfo.myNb0DElements++;
237 return static_cast< SMDS_Mesh0DElement*> ( cell );
243 ///////////////////////////////////////////////////////////////////////////////
244 /// create a Ball and add it to the current Mesh
245 /// @return : The created Ball
246 ///////////////////////////////////////////////////////////////////////////////
247 SMDS_BallElement* SMDS_Mesh::AddBallWithID( int idnode, double diameter, int ID )
249 const SMDS_MeshNode * node = myNodeFactory->FindNode( idnode );
250 if (!node) return NULL;
251 return SMDS_Mesh::AddBallWithID( node, diameter, ID );
254 ///////////////////////////////////////////////////////////////////////////////
255 /// create a Ball and add it to the current Mesh
256 /// @return : The created Ball
257 ///////////////////////////////////////////////////////////////////////////////
258 SMDS_BallElement* SMDS_Mesh::AddBall(const SMDS_MeshNode * node, double diameter)
260 return SMDS_Mesh::AddBallWithID(node, diameter, myCellFactory->GetFreeID());
263 ///////////////////////////////////////////////////////////////////////////////
264 /// Create a new Ball and at it to the mesh
265 /// @param idnode ID of the node
266 // @param diameter ball diameter
267 /// @param ID ID of the 0D element to create
268 /// @return The created 0D element or NULL if an element with this
269 /// ID already exists or if input node is not found.
270 ///////////////////////////////////////////////////////////////////////////////
271 SMDS_BallElement* SMDS_Mesh::AddBallWithID(const SMDS_MeshNode * n, double diameter, int ID)
275 if (NbBalls() % CHECKMEMORY_INTERVAL == 0) CheckMemory();
277 SMDS_BallElement* ball = static_cast< SMDS_BallElement*>( myCellFactory->NewElement( ID ));
280 ball->init( n, diameter );
286 ///////////////////////////////////////////////////////////////////////////////
287 /// create a MeshEdge and add it to the current Mesh
288 /// @return : The created MeshEdge
289 ///////////////////////////////////////////////////////////////////////////////
291 SMDS_MeshEdge* SMDS_Mesh::AddEdgeWithID(int idnode1, int idnode2, int ID)
293 const SMDS_MeshNode * node1 = myNodeFactory->FindNode(idnode1);
294 const SMDS_MeshNode * node2 = myNodeFactory->FindNode(idnode2);
295 if(!node1 || !node2) return NULL;
296 return SMDS_Mesh::AddEdgeWithID(node1, node2, ID);
299 ///////////////////////////////////////////////////////////////////////////////
300 /// create a MeshEdge and add it to the current Mesh
301 /// @return : The created MeshEdge
302 ///////////////////////////////////////////////////////////////////////////////
304 SMDS_MeshEdge* SMDS_Mesh::AddEdge(const SMDS_MeshNode * node1,
305 const SMDS_MeshNode * node2)
307 return SMDS_Mesh::AddEdgeWithID(node1, node2, myCellFactory->GetFreeID());
310 ///////////////////////////////////////////////////////////////////////////////
311 /// Create a new edge and at it to the mesh
312 /// @param idnode1 ID of the first node
313 /// @param idnode2 ID of the second node
314 /// @param ID ID of the edge to create
315 /// @return The created edge or NULL if an element with this ID already exists or
316 /// if input nodes are not found.
317 ///////////////////////////////////////////////////////////////////////////////
319 SMDS_MeshEdge* SMDS_Mesh::AddEdgeWithID(const SMDS_MeshNode * n1,
320 const SMDS_MeshNode * n2,
323 if ( !n1 || !n2 ) return 0;
325 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
327 cell->init( SMDSEntity_Edge, /*nbNodes=*/2, n1, n2 );
329 return static_cast<SMDS_MeshEdge*>( cell );
334 ///////////////////////////////////////////////////////////////////////////////
335 /// Add a triangle defined by its nodes. An ID is automatically affected to the
337 ///////////////////////////////////////////////////////////////////////////////
339 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshNode * n1,
340 const SMDS_MeshNode * n2,
341 const SMDS_MeshNode * n3)
343 return SMDS_Mesh::AddFaceWithID(n1,n2,n3, myCellFactory->GetFreeID());
346 ///////////////////////////////////////////////////////////////////////////////
347 /// Add a triangle defined by its nodes IDs
348 ///////////////////////////////////////////////////////////////////////////////
350 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(int idnode1, int idnode2, int idnode3, int ID)
352 const SMDS_MeshNode * node1 = myNodeFactory->FindNode(idnode1);
353 const SMDS_MeshNode * node2 = myNodeFactory->FindNode(idnode2);
354 const SMDS_MeshNode * node3 = myNodeFactory->FindNode(idnode3);
355 if(!node1 || !node2 || !node3) return NULL;
356 return SMDS_Mesh::AddFaceWithID(node1, node2, node3, ID);
359 ///////////////////////////////////////////////////////////////////////////////
360 /// Add a triangle defined by its nodes
361 ///////////////////////////////////////////////////////////////////////////////
363 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshNode * n1,
364 const SMDS_MeshNode * n2,
365 const SMDS_MeshNode * n3,
368 if ( !n1 || !n2 || !n3 ) return 0;
369 if ( NbFaces() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
371 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
373 cell->init( SMDSEntity_Triangle, /*nbNodes=*/3, n1, n2, n3 );
374 myInfo.myNbTriangles++;
375 return static_cast<SMDS_MeshFace*>( cell );
380 ///////////////////////////////////////////////////////////////////////////////
381 /// Add a quadrangle defined by its nodes. An ID is automatically affected to the
383 ///////////////////////////////////////////////////////////////////////////////
385 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshNode * n1,
386 const SMDS_MeshNode * n2,
387 const SMDS_MeshNode * n3,
388 const SMDS_MeshNode * n4)
390 return SMDS_Mesh::AddFaceWithID(n1,n2,n3, n4, myCellFactory->GetFreeID());
393 ///////////////////////////////////////////////////////////////////////////////
394 /// Add a quadrangle defined by its nodes IDs
395 ///////////////////////////////////////////////////////////////////////////////
397 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(int idnode1,
403 const SMDS_MeshNode *node1, *node2, *node3, *node4;
404 node1 = myNodeFactory->FindNode(idnode1);
405 node2 = myNodeFactory->FindNode(idnode2);
406 node3 = myNodeFactory->FindNode(idnode3);
407 node4 = myNodeFactory->FindNode(idnode4);
408 if ( !node1 || !node2 || !node3 || !node4 ) return NULL;
409 return SMDS_Mesh::AddFaceWithID(node1, node2, node3, node4, ID);
412 ///////////////////////////////////////////////////////////////////////////////
413 /// Add a quadrangle defined by its nodes
414 ///////////////////////////////////////////////////////////////////////////////
416 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshNode * n1,
417 const SMDS_MeshNode * n2,
418 const SMDS_MeshNode * n3,
419 const SMDS_MeshNode * n4,
422 if ( !n1 || !n2 || !n3 || !n4 ) return 0;
423 if ( NbFaces() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
425 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
427 cell->init( SMDSEntity_Quadrangle, /*nbNodes=*/4, n1, n2, n3, n4 );
428 myInfo.myNbQuadrangles++;
429 return static_cast<SMDS_MeshFace*>( cell );
434 ///////////////////////////////////////////////////////////////////////////////
435 ///Create a new tetrahedron and add it to the mesh.
436 ///@return The created tetrahedron
437 ///////////////////////////////////////////////////////////////////////////////
439 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
440 const SMDS_MeshNode * n2,
441 const SMDS_MeshNode * n3,
442 const SMDS_MeshNode * n4)
444 return SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, myCellFactory->GetFreeID() );
447 ///////////////////////////////////////////////////////////////////////////////
448 ///Create a new tetrahedron and add it to the mesh.
449 ///@param ID The ID of the new volume
450 ///@return The created tetrahedron or NULL if an element with this ID already exists
451 ///or if input nodes are not found.
452 ///////////////////////////////////////////////////////////////////////////////
454 SMDS_MeshVolume * SMDS_Mesh::AddVolumeWithID(int idnode1,
460 const SMDS_MeshNode *node1, *node2, *node3, *node4;
461 node1 = myNodeFactory->FindNode(idnode1);
462 node2 = myNodeFactory->FindNode(idnode2);
463 node3 = myNodeFactory->FindNode(idnode3);
464 node4 = myNodeFactory->FindNode(idnode4);
465 if(!node1 || !node2 || !node3 || !node4) return NULL;
466 return SMDS_Mesh::AddVolumeWithID(node1, node2, node3, node4, ID);
469 ///////////////////////////////////////////////////////////////////////////////
470 ///Create a new tetrahedron and add it to the mesh.
471 ///@param ID The ID of the new volume
472 ///@return The created tetrahedron
473 ///////////////////////////////////////////////////////////////////////////////
475 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
476 const SMDS_MeshNode * n2,
477 const SMDS_MeshNode * n3,
478 const SMDS_MeshNode * n4,
481 if ( !n1 || !n2 || !n3 || !n4 ) return 0;
482 if ( NbVolumes() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
484 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
486 cell->init( SMDSEntity_Tetra, /*nbNodes=*/4, n1, n2, n3, n4 );
488 return static_cast<SMDS_MeshVolume*>( cell );
493 ///////////////////////////////////////////////////////////////////////////////
494 ///Create a new pyramid and add it to the mesh.
495 ///Nodes 1,2,3 and 4 define the base of the pyramid
496 ///@return The created pyramid
497 ///////////////////////////////////////////////////////////////////////////////
499 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
500 const SMDS_MeshNode * n2,
501 const SMDS_MeshNode * n3,
502 const SMDS_MeshNode * n4,
503 const SMDS_MeshNode * n5)
505 return SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, myCellFactory->GetFreeID() );
508 ///////////////////////////////////////////////////////////////////////////////
509 ///Create a new pyramid and add it to the mesh.
510 ///Nodes 1,2,3 and 4 define the base of the pyramid
511 ///@param ID The ID of the new volume
512 ///@return The created pyramid or NULL if an element with this ID already exists
513 ///or if input nodes are not found.
514 ///////////////////////////////////////////////////////////////////////////////
516 SMDS_MeshVolume * SMDS_Mesh::AddVolumeWithID(int idnode1,
523 const SMDS_MeshNode *node1, *node2, *node3, *node4, *node5;
524 node1 = myNodeFactory->FindNode(idnode1);
525 node2 = myNodeFactory->FindNode(idnode2);
526 node3 = myNodeFactory->FindNode(idnode3);
527 node4 = myNodeFactory->FindNode(idnode4);
528 node5 = myNodeFactory->FindNode(idnode5);
529 if(!node1 || !node2 || !node3 || !node4 || !node5) return NULL;
530 return SMDS_Mesh::AddVolumeWithID(node1, node2, node3, node4, node5, ID);
533 ///////////////////////////////////////////////////////////////////////////////
534 ///Create a new pyramid and add it to the mesh.
535 ///Nodes 1,2,3 and 4 define the base of the pyramid
536 ///@param ID The ID of the new volume
537 ///@return The created pyramid
538 ///////////////////////////////////////////////////////////////////////////////
540 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
541 const SMDS_MeshNode * n2,
542 const SMDS_MeshNode * n3,
543 const SMDS_MeshNode * n4,
544 const SMDS_MeshNode * n5,
547 if ( !n1 || !n2 || !n3 || !n4 || !n5 ) return 0;
548 if ( NbVolumes() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
550 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
552 cell->init( SMDSEntity_Pyramid, /*nbNodes=*/5, n1, n2, n3, n4, n5 );
553 myInfo.myNbPyramids++;
554 return static_cast<SMDS_MeshVolume*>( cell );
559 ///////////////////////////////////////////////////////////////////////////////
560 ///Create a new prism and add it to the mesh.
561 ///Nodes 1,2,3 is a triangle and 1,2,5,4 a quadrangle.
562 ///@return The created prism
563 ///////////////////////////////////////////////////////////////////////////////
565 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
566 const SMDS_MeshNode * n2,
567 const SMDS_MeshNode * n3,
568 const SMDS_MeshNode * n4,
569 const SMDS_MeshNode * n5,
570 const SMDS_MeshNode * n6)
572 int ID = myCellFactory->GetFreeID();
573 return SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n6, ID);
576 ///////////////////////////////////////////////////////////////////////////////
577 ///Create a new prism and add it to the mesh.
578 ///Nodes 1,2,3 is a triangle and 1,2,5,4 a quadrangle.
579 ///@param ID The ID of the new volume
580 ///@return The created prism or NULL if an element with this ID already exists
581 ///or if input nodes are not found.
582 ///////////////////////////////////////////////////////////////////////////////
584 SMDS_MeshVolume * SMDS_Mesh::AddVolumeWithID(int idnode1,
592 const SMDS_MeshNode *node1, *node2, *node3, *node4, *node5, *node6;
593 node1 = myNodeFactory->FindNode(idnode1);
594 node2 = myNodeFactory->FindNode(idnode2);
595 node3 = myNodeFactory->FindNode(idnode3);
596 node4 = myNodeFactory->FindNode(idnode4);
597 node5 = myNodeFactory->FindNode(idnode5);
598 node6 = myNodeFactory->FindNode(idnode6);
599 return SMDS_Mesh::AddVolumeWithID(node1, node2, node3, node4, node5, node6, ID);
602 ///////////////////////////////////////////////////////////////////////////////
603 ///Create a new prism and add it to the mesh.
604 ///Nodes 1,2,3 is a triangle and 1,2,5,4 a quadrangle.
605 ///@param ID The ID of the new volume
606 ///@return The created prism
607 ///////////////////////////////////////////////////////////////////////////////
609 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
610 const SMDS_MeshNode * n2,
611 const SMDS_MeshNode * n3,
612 const SMDS_MeshNode * n4,
613 const SMDS_MeshNode * n5,
614 const SMDS_MeshNode * n6,
617 if ( !n1 || !n2 || !n3 || !n4 || !n5 || !n6 ) return 0;
618 if ( NbVolumes() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
620 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
622 cell->init( SMDSEntity_Penta, /*nbNodes=*/6, n1, n2, n3, n4, n5, n6 );
624 return static_cast<SMDS_MeshVolume*>( cell );
629 ///////////////////////////////////////////////////////////////////////////////
630 ///Create a new hexagonal prism and add it to the mesh.
631 ///@return The created prism
632 ///////////////////////////////////////////////////////////////////////////////
634 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
635 const SMDS_MeshNode * n2,
636 const SMDS_MeshNode * n3,
637 const SMDS_MeshNode * n4,
638 const SMDS_MeshNode * n5,
639 const SMDS_MeshNode * n6,
640 const SMDS_MeshNode * n7,
641 const SMDS_MeshNode * n8,
642 const SMDS_MeshNode * n9,
643 const SMDS_MeshNode * n10,
644 const SMDS_MeshNode * n11,
645 const SMDS_MeshNode * n12)
647 return SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n6,
648 n7, n8, n9, n10, n11, n12,
649 myCellFactory->GetFreeID() );
652 ///////////////////////////////////////////////////////////////////////////////
653 ///Create a new hexagonal prism and add it to the mesh.
654 ///@param ID The ID of the new volume
655 ///@return The created prism or NULL if an element with this ID already exists
656 ///or if input nodes are not found.
657 ///////////////////////////////////////////////////////////////////////////////
659 SMDS_MeshVolume * SMDS_Mesh::AddVolumeWithID(int idnode1,
673 const SMDS_MeshNode *node1 = myNodeFactory->FindNode(idnode1);
674 const SMDS_MeshNode *node2 = myNodeFactory->FindNode(idnode2);
675 const SMDS_MeshNode *node3 = myNodeFactory->FindNode(idnode3);
676 const SMDS_MeshNode *node4 = myNodeFactory->FindNode(idnode4);
677 const SMDS_MeshNode *node5 = myNodeFactory->FindNode(idnode5);
678 const SMDS_MeshNode *node6 = myNodeFactory->FindNode(idnode6);
679 const SMDS_MeshNode *node7 = myNodeFactory->FindNode(idnode7);
680 const SMDS_MeshNode *node8 = myNodeFactory->FindNode(idnode8);
681 const SMDS_MeshNode *node9 = myNodeFactory->FindNode(idnode9);
682 const SMDS_MeshNode *node10 = myNodeFactory->FindNode(idnode10);
683 const SMDS_MeshNode *node11 = myNodeFactory->FindNode(idnode11);
684 const SMDS_MeshNode *node12 = myNodeFactory->FindNode(idnode12);
685 return SMDS_Mesh::AddVolumeWithID(node1, node2, node3, node4, node5, node6,
686 node7, node8, node9, node10, node11, node12,
690 ///////////////////////////////////////////////////////////////////////////////
691 ///Create a new hexagonal prism and add it to the mesh.
692 ///@param ID The ID of the new volume
693 ///@return The created prism
694 ///////////////////////////////////////////////////////////////////////////////
696 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
697 const SMDS_MeshNode * n2,
698 const SMDS_MeshNode * n3,
699 const SMDS_MeshNode * n4,
700 const SMDS_MeshNode * n5,
701 const SMDS_MeshNode * n6,
702 const SMDS_MeshNode * n7,
703 const SMDS_MeshNode * n8,
704 const SMDS_MeshNode * n9,
705 const SMDS_MeshNode * n10,
706 const SMDS_MeshNode * n11,
707 const SMDS_MeshNode * n12,
710 SMDS_MeshVolume* volume = 0;
711 if(!n1 || !n2 || !n3 || !n4 || !n5 || !n6 ||
712 !n7 || !n8 || !n9 || !n10 || !n11 || !n12 )
714 if ( NbVolumes() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
716 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
718 cell->init( SMDSEntity_Hexagonal_Prism,
719 /*nbNodes=*/12, n1, n2, n3, n4, n5, n6, n7, n8, n9, n10, n11, n12 );
720 myInfo.myNbHexPrism++;
721 return static_cast<SMDS_MeshVolume*>( cell );
726 ///////////////////////////////////////////////////////////////////////////////
727 ///Create a new hexahedron and add it to the mesh.
728 ///Nodes 1,2,3,4 and 5,6,7,8 are quadrangle and 5,1 and 7,3 are an edges.
729 ///@return The created hexahedron
730 ///////////////////////////////////////////////////////////////////////////////
732 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
733 const SMDS_MeshNode * n2,
734 const SMDS_MeshNode * n3,
735 const SMDS_MeshNode * n4,
736 const SMDS_MeshNode * n5,
737 const SMDS_MeshNode * n6,
738 const SMDS_MeshNode * n7,
739 const SMDS_MeshNode * n8)
741 int ID = myCellFactory->GetFreeID();
742 return SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n6, n7, n8, ID);
745 ///////////////////////////////////////////////////////////////////////////////
746 ///Create a new hexahedron and add it to the mesh.
747 ///Nodes 1,2,3,4 and 5,6,7,8 are quadrangle and 5,1 and 7,3 are an edges.
748 ///@param ID The ID of the new volume
749 ///@return The created hexahedron or NULL if an element with this ID already
750 ///exists or if input nodes are not found.
751 ///////////////////////////////////////////////////////////////////////////////
753 SMDS_MeshVolume * SMDS_Mesh::AddVolumeWithID(int idnode1,
763 const SMDS_MeshNode *node1, *node2, *node3, *node4, *node5, *node6, *node7, *node8;
764 node1 = myNodeFactory->FindNode(idnode1);
765 node2 = myNodeFactory->FindNode(idnode2);
766 node3 = myNodeFactory->FindNode(idnode3);
767 node4 = myNodeFactory->FindNode(idnode4);
768 node5 = myNodeFactory->FindNode(idnode5);
769 node6 = myNodeFactory->FindNode(idnode6);
770 node7 = myNodeFactory->FindNode(idnode7);
771 node8 = myNodeFactory->FindNode(idnode8);
772 return SMDS_Mesh::AddVolumeWithID(node1, node2, node3, node4, node5, node6,
776 ///////////////////////////////////////////////////////////////////////////////
777 ///Create a new hexahedron and add it to the mesh.
778 ///Nodes 1,2,3,4 and 5,6,7,8 are quadrangle and 5,1 and 7,3 are an edges.
779 ///@param ID The ID of the new volume
780 ///@return The created prism or NULL if an element with this ID already exists
781 ///or if input nodes are not found.
782 ///////////////////////////////////////////////////////////////////////////////
784 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
785 const SMDS_MeshNode * n2,
786 const SMDS_MeshNode * n3,
787 const SMDS_MeshNode * n4,
788 const SMDS_MeshNode * n5,
789 const SMDS_MeshNode * n6,
790 const SMDS_MeshNode * n7,
791 const SMDS_MeshNode * n8,
794 SMDS_MeshVolume* volume = 0;
795 if ( !n1 || !n2 || !n3 || !n4 || !n5 || !n6 || !n7 || !n8) return volume;
796 if ( NbVolumes() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
798 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
800 cell->init( SMDSEntity_Hexa,
801 /*nbNodes=*/8, n1, n2, n3, n4, n5, n6, n7, n8 );
803 return static_cast<SMDS_MeshVolume*>( cell );
808 ///////////////////////////////////////////////////////////////////////////////
809 /// Add a polygon defined by its nodes IDs
810 ///////////////////////////////////////////////////////////////////////////////
812 SMDS_MeshFace* SMDS_Mesh::AddPolygonalFaceWithID (const std::vector<int> & nodes_ids,
815 int nbNodes = nodes_ids.size();
816 std::vector<const SMDS_MeshNode*> nodes (nbNodes);
817 for (int i = 0; i < nbNodes; i++) {
818 nodes[i] = myNodeFactory->FindNode( nodes_ids[i] );
819 if (!nodes[i]) return NULL;
821 return SMDS_Mesh::AddPolygonalFaceWithID(nodes, ID);
824 ///////////////////////////////////////////////////////////////////////////////
825 /// Add a polygon defined by its nodes
826 ///////////////////////////////////////////////////////////////////////////////
829 SMDS_Mesh::AddPolygonalFaceWithID (const std::vector<const SMDS_MeshNode*> & nodes,
832 if ( NbFaces() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
835 throw std::invalid_argument("Polygon without nodes is forbidden");
836 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
838 cell->init( SMDSEntity_Polygon, nodes );
839 myInfo.myNbPolygons++;
840 return static_cast<SMDS_MeshFace*>( cell );
845 ///////////////////////////////////////////////////////////////////////////////
846 /// Add a polygon defined by its nodes.
847 /// An ID is automatically affected to the created face.
848 ///////////////////////////////////////////////////////////////////////////////
850 SMDS_MeshFace* SMDS_Mesh::AddPolygonalFace (const std::vector<const SMDS_MeshNode*> & nodes)
852 return SMDS_Mesh::AddPolygonalFaceWithID(nodes, myCellFactory->GetFreeID());
855 ///////////////////////////////////////////////////////////////////////////////
856 /// Add a quadratic polygon defined by its nodes IDs
857 ///////////////////////////////////////////////////////////////////////////////
859 SMDS_MeshFace* SMDS_Mesh::AddQuadPolygonalFaceWithID (const std::vector<int> & nodes_ids,
862 std::vector<const SMDS_MeshNode*> nodes( nodes_ids.size() );
863 for ( size_t i = 0; i < nodes.size(); i++) {
864 nodes[i] = myNodeFactory->FindNode(nodes_ids[i]);
865 if (!nodes[i]) return NULL;
867 return SMDS_Mesh::AddQuadPolygonalFaceWithID(nodes, ID);
870 ///////////////////////////////////////////////////////////////////////////////
871 /// Add a quadratic polygon defined by its nodes
872 ///////////////////////////////////////////////////////////////////////////////
875 SMDS_Mesh::AddQuadPolygonalFaceWithID (const std::vector<const SMDS_MeshNode*> & nodes,
878 if ( NbFaces() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
880 throw std::invalid_argument("Polygon without nodes is forbidden");
881 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
883 cell->init( SMDSEntity_Quad_Polygon, nodes );
884 myInfo.myNbQuadPolygons++;
885 return static_cast<SMDS_MeshFace*>( cell );
890 ///////////////////////////////////////////////////////////////////////////////
891 /// Add a quadratic polygon defined by its nodes.
892 /// An ID is automatically affected to the created face.
893 ///////////////////////////////////////////////////////////////////////////////
895 SMDS_MeshFace* SMDS_Mesh::AddQuadPolygonalFace (const std::vector<const SMDS_MeshNode*> & nodes)
897 return SMDS_Mesh::AddQuadPolygonalFaceWithID(nodes, myCellFactory->GetFreeID());
900 ///////////////////////////////////////////////////////////////////////////////
901 /// Create a new polyhedral volume and add it to the mesh.
902 /// @param ID The ID of the new volume
903 /// @return The created volume or NULL if an element with this ID already exists
904 /// or if input nodes are not found.
905 ///////////////////////////////////////////////////////////////////////////////
907 SMDS_MeshVolume * SMDS_Mesh::AddPolyhedralVolumeWithID (const std::vector<int> & nodes_ids,
908 const std::vector<int> & quantities,
911 int nbNodes = nodes_ids.size();
912 std::vector<const SMDS_MeshNode*> nodes (nbNodes);
913 for (int i = 0; i < nbNodes; i++) {
914 nodes[i] = myNodeFactory->FindNode(nodes_ids[i]);
915 if (!nodes[i]) return NULL;
917 return SMDS_Mesh::AddPolyhedralVolumeWithID(nodes, quantities, ID);
920 ///////////////////////////////////////////////////////////////////////////////
921 /// Create a new polyhedral volume and add it to the mesh.
922 /// @param ID The ID of the new volume
923 /// @return The created volume
924 ///////////////////////////////////////////////////////////////////////////////
927 SMDS_Mesh::AddPolyhedralVolumeWithID (const std::vector<const SMDS_MeshNode*>& nodes,
928 const std::vector<int> & quantities,
931 if ( nodes.empty() || quantities.empty() )
933 if ( NbVolumes() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
935 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
937 SMDS_MeshVolume* volume = static_cast<SMDS_MeshVolume*>( cell );
938 volume->init( nodes, quantities );
939 myInfo.myNbPolyhedrons++;
945 ///////////////////////////////////////////////////////////////////////////////
946 /// Create a new polyhedral volume and add it to the mesh.
947 /// @return The created volume
948 ///////////////////////////////////////////////////////////////////////////////
950 SMDS_MeshVolume* SMDS_Mesh::AddPolyhedralVolume
951 (const std::vector<const SMDS_MeshNode*> & nodes,
952 const std::vector<int> & quantities)
954 int ID = myCellFactory->GetFreeID();
955 return SMDS_Mesh::AddPolyhedralVolumeWithID(nodes, quantities, ID);
958 SMDS_MeshVolume* SMDS_Mesh::AddVolumeFromVtkIds(const std::vector<vtkIdType>& vtkNodeIds)
960 SMDS_MeshCell* cell = myCellFactory->NewCell( myCellFactory->GetFreeID() );
961 SMDS_MeshVolume * vol = static_cast<SMDS_MeshVolume*>( cell );
962 vol->init( vtkNodeIds );
967 SMDS_MeshFace* SMDS_Mesh::AddFaceFromVtkIds(const std::vector<vtkIdType>& vtkNodeIds)
969 SMDS_MeshCell* cell = myCellFactory->NewCell( myCellFactory->GetFreeID() );
970 SMDS_MeshFace * f = static_cast<SMDS_MeshFace*>( cell );
971 f->init( vtkNodeIds );
976 //=======================================================================
977 //function : MoveNode
979 //=======================================================================
981 void SMDS_Mesh::MoveNode(const SMDS_MeshNode *n, double x, double y, double z)
983 SMDS_MeshNode * node=const_cast<SMDS_MeshNode*>(n);
987 ///////////////////////////////////////////////////////////////////////////////
988 /// Return the node whose SMDS ID is 'ID'.
989 ///////////////////////////////////////////////////////////////////////////////
990 const SMDS_MeshNode * SMDS_Mesh::FindNode(int ID) const
992 return myNodeFactory->FindNode( ID );
995 ///////////////////////////////////////////////////////////////////////////////
996 /// Return the node whose VTK ID is 'vtkId'.
997 ///////////////////////////////////////////////////////////////////////////////
998 const SMDS_MeshNode * SMDS_Mesh::FindNodeVtk(int vtkId) const
1000 return myNodeFactory->FindNode( vtkId + 1 );
1003 const SMDS_MeshElement * SMDS_Mesh::FindElementVtk(int IDelem) const
1005 return myCellFactory->FindElement( FromVtkToSmds( IDelem ));
1008 ///////////////////////////////////////////////////////////////////////////////
1009 /// Remove a node and all the elements which own this node
1010 ///////////////////////////////////////////////////////////////////////////////
1012 void SMDS_Mesh::RemoveNode(const SMDS_MeshNode * node)
1014 RemoveElement(node, true);
1017 //=======================================================================
1018 //function : RemoveFromParent
1020 //=======================================================================
1022 bool SMDS_Mesh::RemoveFromParent()
1024 if (myParent==NULL) return false;
1025 else return (myParent->RemoveSubMesh(this));
1028 //=======================================================================
1029 //function : RemoveSubMesh
1031 //=======================================================================
1033 bool SMDS_Mesh::RemoveSubMesh(const SMDS_Mesh * aMesh)
1037 std::list<SMDS_Mesh *>::iterator itmsh=myChildren.begin();
1038 for (; itmsh!=myChildren.end() && !found; itmsh++)
1040 SMDS_Mesh * submesh = *itmsh;
1041 if (submesh == aMesh)
1044 myChildren.erase(itmsh);
1051 //=======================================================================
1052 //function : ChangePolyhedronNodes
1054 //=======================================================================
1056 bool SMDS_Mesh::ChangePolyhedronNodes(const SMDS_MeshElement * element,
1057 const std::vector<const SMDS_MeshNode*>& nodes,
1058 const std::vector<int>& quantities)
1060 // keep current nodes of element
1061 std::set<const SMDS_MeshNode*> oldNodes( element->begin_nodes(), element->end_nodes() );
1065 if ( const SMDS_MeshVolume* vol = DownCast<SMDS_MeshVolume>( element ))
1066 Ok = vol->ChangeNodes( nodes, quantities );
1071 updateInverseElements( element, &nodes[0], nodes.size(), oldNodes );
1076 //=======================================================================
1077 //function : ChangeElementNodes
1079 //=======================================================================
1081 bool SMDS_Mesh::ChangeElementNodes(const SMDS_MeshElement * element,
1082 const SMDS_MeshNode * nodes[],
1085 // keep current nodes of element
1086 std::set<const SMDS_MeshNode*> oldNodes( element->begin_nodes(), element->end_nodes() );
1090 if ( SMDS_MeshCell* cell = dynamic_cast<SMDS_MeshCell*>((SMDS_MeshElement*) element))
1091 Ok = cell->ChangeNodes(nodes, nbnodes);
1096 updateInverseElements( element, nodes, nbnodes, oldNodes );
1101 //=======================================================================
1102 //function : updateInverseElements
1103 //purpose : update InverseElements when element changes node
1104 //=======================================================================
1106 void SMDS_Mesh::updateInverseElements( const SMDS_MeshElement * element,
1107 const SMDS_MeshNode* const* nodes,
1109 std::set<const SMDS_MeshNode*>& oldNodes )
1111 if ( GetGrid()->HasLinks() ) // update InverseElements
1113 std::set<const SMDS_MeshNode*>::iterator it;
1115 // AddInverseElement to new nodes
1116 for ( int i = 0; i < nbnodes; i++ )
1118 it = oldNodes.find( nodes[i] );
1119 if ( it == oldNodes.end() )
1121 const_cast<SMDS_MeshNode*>( nodes[i] )->AddInverseElement( element );
1123 // remove from oldNodes a node that remains in elem
1124 oldNodes.erase( it );
1126 // RemoveInverseElement from the nodes removed from elem
1127 for ( it = oldNodes.begin(); it != oldNodes.end(); it++ )
1129 SMDS_MeshNode * n = const_cast<SMDS_MeshNode *>( *it );
1130 n->RemoveInverseElement( element );
1136 const SMDS_Mesh0DElement* SMDS_Mesh::Find0DElement(const SMDS_MeshNode * node)
1138 if (!node) return 0;
1139 const SMDS_Mesh0DElement* toReturn = NULL;
1140 SMDS_ElemIteratorPtr it1 = node->GetInverseElementIterator(SMDSAbs_0DElement);
1141 while (it1->more() && (toReturn == NULL)) {
1142 const SMDS_MeshElement* e = it1->next();
1143 if (e->NbNodes() == 1) {
1144 toReturn = static_cast<const SMDS_Mesh0DElement*>(e);
1150 const SMDS_BallElement* SMDS_Mesh::FindBall(const SMDS_MeshNode * node)
1152 if (!node) return 0;
1153 const SMDS_BallElement* toReturn = NULL;
1154 SMDS_ElemIteratorPtr it1 = node->GetInverseElementIterator(SMDSAbs_Ball);
1155 while (it1->more() && (toReturn == NULL)) {
1156 const SMDS_MeshElement* e = it1->next();
1157 if (e->GetGeomType() == SMDSGeom_BALL)
1158 toReturn = static_cast<const SMDS_BallElement*>(e);
1163 const SMDS_MeshEdge* SMDS_Mesh::FindEdge(const SMDS_MeshNode * node1,
1164 const SMDS_MeshNode * node2)
1166 if ( !node1 ) return 0;
1167 const SMDS_MeshEdge * toReturn=NULL;
1168 SMDS_ElemIteratorPtr it1=node1->GetInverseElementIterator(SMDSAbs_Edge);
1169 while(it1->more()) {
1170 const SMDS_MeshElement * e = it1->next();
1171 if ( e->NbNodes() == 2 && e->GetNodeIndex( node2 ) >= 0 ) {
1172 toReturn = static_cast<const SMDS_MeshEdge*>( e );
1179 const SMDS_MeshEdge* SMDS_Mesh::FindEdge(const SMDS_MeshNode * node1,
1180 const SMDS_MeshNode * node2,
1181 const SMDS_MeshNode * node3)
1183 if ( !node1 ) return 0;
1184 SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Edge);
1185 while(it1->more()) {
1186 const SMDS_MeshElement * e = it1->next();
1187 if ( e->NbNodes() == 3 ) {
1188 SMDS_ElemIteratorPtr it2 = e->nodesIterator();
1189 while(it2->more()) {
1190 const SMDS_MeshElement* n = it2->next();
1200 return static_cast<const SMDS_MeshEdge *> (e);
1206 //=======================================================================
1207 //function : FindFace
1209 //=======================================================================
1211 const SMDS_MeshFace* SMDS_Mesh::FindFace(const SMDS_MeshNode *node1,
1212 const SMDS_MeshNode *node2,
1213 const SMDS_MeshNode *node3)
1215 if ( !node1 ) return 0;
1216 SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Face);
1217 while(it1->more()) {
1218 const SMDS_MeshElement * e = it1->next();
1219 if ( e->NbNodes() == 3 ) {
1220 SMDS_ElemIteratorPtr it2 = e->nodesIterator();
1221 while(it2->more()) {
1222 const SMDS_MeshElement* n = it2->next();
1232 return static_cast<const SMDS_MeshFace *> (e);
1238 //=======================================================================
1239 //function : FindFace
1241 //=======================================================================
1243 const SMDS_MeshFace* SMDS_Mesh::FindFace(const SMDS_MeshNode *node1,
1244 const SMDS_MeshNode *node2,
1245 const SMDS_MeshNode *node3,
1246 const SMDS_MeshNode *node4)
1248 if ( !node1 ) return 0;
1249 SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Face);
1250 while(it1->more()) {
1251 const SMDS_MeshElement * e = it1->next();
1252 if ( e->NbNodes() == 4 ) {
1253 SMDS_ElemIteratorPtr it2 = e->nodesIterator();
1254 while(it2->more()) {
1255 const SMDS_MeshElement* n = it2->next();
1266 return static_cast<const SMDS_MeshFace *> (e);
1272 //=======================================================================
1273 //function : FindFace
1274 //purpose :quadratic triangle
1275 //=======================================================================
1277 const SMDS_MeshFace* SMDS_Mesh::FindFace(const SMDS_MeshNode *node1,
1278 const SMDS_MeshNode *node2,
1279 const SMDS_MeshNode *node3,
1280 const SMDS_MeshNode *node4,
1281 const SMDS_MeshNode *node5,
1282 const SMDS_MeshNode *node6)
1284 if ( !node1 ) return 0;
1285 SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Face);
1286 while(it1->more()) {
1287 const SMDS_MeshElement * e = it1->next();
1288 if ( e->NbNodes() == 6 ) {
1289 SMDS_ElemIteratorPtr it2 = e->nodesIterator();
1290 while(it2->more()) {
1291 const SMDS_MeshElement* n = it2->next();
1304 return static_cast<const SMDS_MeshFace *> (e);
1311 //=======================================================================
1312 //function : FindFace
1313 //purpose : quadratic quadrangle
1314 //=======================================================================
1316 const SMDS_MeshFace* SMDS_Mesh::FindFace(const SMDS_MeshNode *node1,
1317 const SMDS_MeshNode *node2,
1318 const SMDS_MeshNode *node3,
1319 const SMDS_MeshNode *node4,
1320 const SMDS_MeshNode *node5,
1321 const SMDS_MeshNode *node6,
1322 const SMDS_MeshNode *node7,
1323 const SMDS_MeshNode *node8)
1325 if ( !node1 ) return 0;
1326 SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Face);
1327 while(it1->more()) {
1328 const SMDS_MeshElement * e = it1->next();
1329 if ( e->NbNodes() == 8 ) {
1330 SMDS_ElemIteratorPtr it2 = e->nodesIterator();
1331 while(it2->more()) {
1332 const SMDS_MeshElement* n = it2->next();
1347 return static_cast<const SMDS_MeshFace *> (e);
1354 //=======================================================================
1355 //function : FindElement
1357 //=======================================================================
1359 const SMDS_MeshElement* SMDS_Mesh::FindElement(int IDelem) const
1361 return myCellFactory->FindElement( IDelem );
1364 //=======================================================================
1365 //function : FindFace
1366 //purpose : find polygon
1367 //=======================================================================
1370 const SMDS_MeshFace* SMDS_Mesh::FindFace (const std::vector<const SMDS_MeshNode *>& nodes)
1372 return (const SMDS_MeshFace*) FindElement( nodes, SMDSAbs_Face );
1376 //================================================================================
1378 * \brief Return element based on all given nodes
1379 * \param nodes - node of element
1380 * \param type - type of element
1381 * \param noMedium - true if medium nodes of quadratic element are not included in <nodes>
1382 * \retval const SMDS_MeshElement* - found element or NULL
1384 //================================================================================
1386 const SMDS_MeshElement* SMDS_Mesh::FindElement (const std::vector<const SMDS_MeshNode *>& nodes,
1387 const SMDSAbs_ElementType type,
1388 const bool noMedium)
1390 if ( nodes.size() > 0 && nodes[0] )
1392 SMDS_ElemIteratorPtr itF = nodes[0]->GetInverseElementIterator(type);
1395 const SMDS_MeshElement* e = itF->next();
1396 int nbNodesToCheck = noMedium ? e->NbCornerNodes() : e->NbNodes();
1397 if ( nbNodesToCheck == (int)nodes.size() )
1399 for ( size_t i = 1; e && i < nodes.size(); ++i )
1401 int nodeIndex = e->GetNodeIndex( nodes[ i ]);
1402 if ( nodeIndex < 0 || nodeIndex >= nbNodesToCheck )
1413 //================================================================================
1415 * \brief Return elements including all given nodes
1416 * \param [in] nodes - nodes to find elements around
1417 * \param [out] foundElems - the found elements
1418 * \param [in] type - type of elements to find
1419 * \return int - a number of found elements
1421 //================================================================================
1423 int SMDS_Mesh::GetElementsByNodes(const std::vector<const SMDS_MeshNode *>& nodes,
1424 std::vector<const SMDS_MeshElement *>& foundElems,
1425 const SMDSAbs_ElementType type)
1427 // chose a node with minimal number of inverse elements
1428 const SMDS_MeshNode* n0 = nodes[0];
1429 int minNbInverse = n0 ? n0->NbInverseElements( type ) : 1000;
1430 for ( size_t i = 1; i < nodes.size(); ++i )
1431 if ( nodes[i] && nodes[i]->NbInverseElements( type ) < minNbInverse )
1434 minNbInverse = n0->NbInverseElements( type );
1440 foundElems.reserve( minNbInverse );
1441 SMDS_ElemIteratorPtr eIt = n0->GetInverseElementIterator( type );
1442 while ( eIt->more() )
1444 const SMDS_MeshElement* e = eIt->next();
1445 bool includeAll = true;
1446 for ( size_t i = 0; i < nodes.size() && includeAll; ++i )
1447 if ( nodes[i] != n0 && e->GetNodeIndex( nodes[i] ) < 0 )
1450 foundElems.push_back( e );
1453 return foundElems.size();
1456 ///////////////////////////////////////////////////////////////////////////////
1457 /// Return the number of nodes
1458 ///////////////////////////////////////////////////////////////////////////////
1459 int SMDS_Mesh::NbNodes() const
1461 return myInfo.NbNodes();
1464 ///////////////////////////////////////////////////////////////////////////////
1465 /// Return the number of elements
1466 ///////////////////////////////////////////////////////////////////////////////
1467 int SMDS_Mesh::NbElements() const
1469 return myInfo.NbElements();
1471 ///////////////////////////////////////////////////////////////////////////////
1472 /// Return the number of 0D elements
1473 ///////////////////////////////////////////////////////////////////////////////
1474 int SMDS_Mesh::Nb0DElements() const
1476 return myInfo.Nb0DElements();
1479 ///////////////////////////////////////////////////////////////////////////////
1480 /// Return the number of 0D elements
1481 ///////////////////////////////////////////////////////////////////////////////
1482 int SMDS_Mesh::NbBalls() const
1484 return myInfo.NbBalls();
1487 ///////////////////////////////////////////////////////////////////////////////
1488 /// Return the number of edges (including construction edges)
1489 ///////////////////////////////////////////////////////////////////////////////
1490 int SMDS_Mesh::NbEdges() const
1492 return myInfo.NbEdges();
1495 ///////////////////////////////////////////////////////////////////////////////
1496 /// Return the number of faces (including construction faces)
1497 ///////////////////////////////////////////////////////////////////////////////
1498 int SMDS_Mesh::NbFaces() const
1500 return myInfo.NbFaces();
1503 ///////////////////////////////////////////////////////////////////////////////
1504 /// Return the number of volumes
1505 ///////////////////////////////////////////////////////////////////////////////
1506 int SMDS_Mesh::NbVolumes() const
1508 return myInfo.NbVolumes();
1511 ///////////////////////////////////////////////////////////////////////////////
1512 /// Return the number of child mesh of this mesh.
1513 /// Note that the tree structure of SMDS_Mesh is unused in SMESH
1514 ///////////////////////////////////////////////////////////////////////////////
1515 int SMDS_Mesh::NbSubMesh() const
1517 return myChildren.size();
1520 ///////////////////////////////////////////////////////////////////////////////
1521 /// Destroy the mesh and all its elements
1522 /// All pointer on elements owned by this mesh become illegals.
1523 ///////////////////////////////////////////////////////////////////////////////
1524 SMDS_Mesh::~SMDS_Mesh()
1526 std::list<SMDS_Mesh*>::iterator itc=myChildren.begin();
1527 while(itc!=myChildren.end())
1533 delete myNodeFactory;
1534 delete myCellFactory;
1539 //================================================================================
1541 * \brief Clear all data
1543 //================================================================================
1545 void SMDS_Mesh::Clear()
1547 std::set< SMDS_ElementHolder* >::iterator holder = myElemHolders.begin();
1548 for ( ; holder != myElemHolders.end(); ++holder )
1551 myNodeFactory->Clear();
1552 myCellFactory->Clear();
1554 std::list<SMDS_Mesh*>::iterator itc=myChildren.begin();
1555 while(itc!=myChildren.end())
1566 myGrid->Initialize();
1568 vtkPoints* points = vtkPoints::New();
1569 // rnv: to fix bug "21125: EDF 1233 SMESH: Degrardation of precision in a test case for quadratic conversion"
1570 // using double type for storing coordinates of nodes instead float.
1571 points->SetDataType(VTK_DOUBLE);
1572 points->SetNumberOfPoints( 0 );
1573 myGrid->SetPoints( points );
1575 myGrid->DeleteLinks();
1578 ///////////////////////////////////////////////////////////////////////////////
1579 /// Return an iterator on nodes of the current mesh factory
1580 ///////////////////////////////////////////////////////////////////////////////
1582 SMDS_NodeIteratorPtr SMDS_Mesh::nodesIterator() const
1584 return myNodeFactory->GetIterator< SMDS_NodeIterator >( new SMDS_MeshElement::NonNullFilter );
1587 SMDS_ElemIteratorPtr SMDS_Mesh::elementGeomIterator(SMDSAbs_GeometryType type) const
1589 int nbElems = myCellFactory->CompactChangePointers() ? -1 : myInfo.NbElements( type );
1590 return myCellFactory->GetIterator< SMDS_ElemIterator >( new SMDS_MeshElement::GeomFilter( type ),
1594 SMDS_ElemIteratorPtr SMDS_Mesh::elementEntityIterator(SMDSAbs_EntityType type) const
1596 if ( type == SMDSEntity_Node )
1598 return myNodeFactory->GetIterator< SMDS_ElemIterator >( new SMDS_MeshElement::NonNullFilter );
1600 int nbElems = myCellFactory->CompactChangePointers() ? -1 : myInfo.NbElements( type );
1601 return myCellFactory->GetIterator<SMDS_ElemIterator>( new SMDS_MeshElement::EntityFilter( type ),
1605 ///////////////////////////////////////////////////////////////////////////////
1606 /// Return an iterator on elements of the current mesh factory
1607 ///////////////////////////////////////////////////////////////////////////////
1608 SMDS_ElemIteratorPtr SMDS_Mesh::elementsIterator(SMDSAbs_ElementType type) const
1610 typedef SMDS_ElemIterator TIterator;
1614 return myCellFactory->GetIterator< TIterator >( new SMDS_MeshElement::NonNullFilter );
1617 return myNodeFactory->GetIterator< TIterator >( new SMDS_MeshElement::NonNullFilter );
1620 int nbElems = myCellFactory->CompactChangePointers() ? -1 : myInfo.NbElements( type );
1621 return myCellFactory->GetIterator< TIterator >( new SMDS_MeshElement::TypeFilter( type ),
1624 return SMDS_ElemIteratorPtr();
1627 ///////////////////////////////////////////////////////////////////////////////
1628 ///Return an iterator on edges of the current mesh.
1629 ///////////////////////////////////////////////////////////////////////////////
1631 SMDS_EdgeIteratorPtr SMDS_Mesh::edgesIterator() const
1633 typedef SMDS_EdgeIterator TIterator;
1634 int nbElems = myCellFactory->CompactChangePointers() ? -1 : myInfo.NbEdges();
1635 return myCellFactory->GetIterator< TIterator >( new SMDS_MeshElement::TypeFilter( SMDSAbs_Edge ),
1639 ///////////////////////////////////////////////////////////////////////////////
1640 ///Return an iterator on faces of the current mesh.
1641 ///////////////////////////////////////////////////////////////////////////////
1643 SMDS_FaceIteratorPtr SMDS_Mesh::facesIterator() const
1645 typedef SMDS_FaceIterator TIterator;
1646 int nbElems = myCellFactory->CompactChangePointers() ? -1 : myInfo.NbFaces();
1647 return myCellFactory->GetIterator< TIterator >( new SMDS_MeshElement::TypeFilter( SMDSAbs_Face ),
1651 ///////////////////////////////////////////////////////////////////////////////
1652 ///Return an iterator on volumes of the current mesh.
1653 ///////////////////////////////////////////////////////////////////////////////
1655 SMDS_VolumeIteratorPtr SMDS_Mesh::volumesIterator() const
1657 typedef SMDS_VolumeIterator TIterator;
1658 int nbElems = myCellFactory->CompactChangePointers() ? -1 : myInfo.NbVolumes();
1660 myCellFactory->GetIterator< TIterator >( new SMDS_MeshElement::TypeFilter( SMDSAbs_Volume ),
1664 SMDS_NodeIteratorPtr SMDS_Mesh::shapeNodesIterator(int shapeID,
1665 size_t nbElemsToReturn,
1666 const SMDS_MeshNode* sm1stNode) const
1668 return myNodeFactory->GetShapeIterator< SMDS_NodeIterator >( shapeID, nbElemsToReturn, sm1stNode );
1671 SMDS_ElemIteratorPtr SMDS_Mesh::shapeElementsIterator(int shapeID,
1672 size_t nbElemsToReturn,
1673 const SMDS_MeshElement* sm1stElem) const
1675 return myCellFactory->GetShapeIterator< SMDS_ElemIterator >( shapeID, nbElemsToReturn, sm1stElem );
1678 ///////////////////////////////////////////////////////////////////////////////
1679 /// Do intersection of sets (more than 2)
1680 ///////////////////////////////////////////////////////////////////////////////
1681 static std::set<const SMDS_MeshElement*> *
1682 intersectionOfSets( std::set<const SMDS_MeshElement*> vs[], int numberOfSets )
1684 std::set<const SMDS_MeshElement*>* rsetA = new std::set<const SMDS_MeshElement*>(vs[0]);
1685 std::set<const SMDS_MeshElement*>* rsetB;
1687 for(int i=0; i<numberOfSets-1; i++)
1689 rsetB = new std::set<const SMDS_MeshElement*>();
1690 set_intersection(rsetA->begin(), rsetA->end(),
1691 vs[i+1].begin(), vs[i+1].end(),
1692 inserter(*rsetB, rsetB->begin()));
1698 ///////////////////////////////////////////////////////////////////////////////
1699 /// Return the list of finite elements owning the given element: elements
1700 /// containing all the nodes of the given element, for instance faces and
1701 /// volumes containing a given edge.
1702 ///////////////////////////////////////////////////////////////////////////////
1703 static std::set<const SMDS_MeshElement*> * getFinitElements(const SMDS_MeshElement * element)
1705 int numberOfSets=element->NbNodes();
1706 std::set<const SMDS_MeshElement*> *initSet = new std::set<const SMDS_MeshElement*>[numberOfSets];
1708 SMDS_NodeIteratorPtr itNodes = element->nodeIterator();
1711 while ( itNodes->more() )
1713 const SMDS_MeshNode * n = itNodes->next();
1714 for ( SMDS_ElemIteratorPtr itFe = n->GetInverseElementIterator(); itFe->more(); )
1715 initSet[i].insert( itFe->next() );
1718 std::set<const SMDS_MeshElement*> *retSet = intersectionOfSets( initSet, numberOfSets );
1723 ///////////////////////////////////////////////////////////////////////////////
1724 /// Return the std::list of nodes used only by the given elements
1725 ///////////////////////////////////////////////////////////////////////////////
1727 std::set<const SMDS_MeshElement*> *getExclusiveNodes(std::set<const SMDS_MeshElement*>& elements)
1729 std::set<const SMDS_MeshElement*> * toReturn = new std::set<const SMDS_MeshElement*>();
1730 std::set<const SMDS_MeshElement*>::iterator itElements = elements.begin();
1732 while( itElements != elements.end() )
1734 SMDS_NodeIteratorPtr itNodes = (*itElements)->nodeIterator();
1737 while( itNodes->more() )
1739 const SMDS_MeshNode * n = itNodes->next();
1740 SMDS_ElemIteratorPtr itFe = n->GetInverseElementIterator();
1741 std::set<const SMDS_MeshElement*> s;
1742 while ( itFe->more() )
1743 s.insert( itFe->next() );
1744 if ( s == elements ) toReturn->insert(n);
1750 ///////////////////////////////////////////////////////////////////////////////
1751 ///Find the children of an element that are made of given nodes
1752 ///@param setOfChildren The set in which matching children will be inserted
1753 ///@param element The element were to search matching children
1754 ///@param nodes The nodes that the children must have to be selected
1755 ///////////////////////////////////////////////////////////////////////////////
1756 void SMDS_Mesh::addChildrenWithNodes(std::set<const SMDS_MeshElement*>& setOfChildren,
1757 const SMDS_MeshElement * element,
1758 std::set<const SMDS_MeshElement*>& nodes)
1760 switch(element->GetType())
1763 throw SALOME_Exception("Internal Error: This should not happen");
1765 case SMDSAbs_0DElement:
1772 SMDS_ElemIteratorPtr itn=element->nodesIterator();
1775 const SMDS_MeshElement * e=itn->next();
1776 if(nodes.find(e)!=nodes.end())
1778 setOfChildren.insert(element);
1785 SMDS_ElemIteratorPtr itn=element->nodesIterator();
1788 const SMDS_MeshElement * e=itn->next();
1789 if(nodes.find(e)!=nodes.end())
1791 setOfChildren.insert(element);
1796 case SMDSAbs_Volume:
1797 case SMDSAbs_NbElementTypes:
1798 case SMDSAbs_All: break;
1802 ///////////////////////////////////////////////////////////////////////////////
1803 ///@param elem The element to delete
1804 ///@param removenodes if true remaining nodes will be removed
1805 ///////////////////////////////////////////////////////////////////////////////
1806 void SMDS_Mesh::RemoveElement(const SMDS_MeshElement * elem,
1807 const bool removenodes)
1809 std::vector<const SMDS_MeshElement *> removedElems;
1810 std::vector<const SMDS_MeshElement *> removedNodes;
1811 RemoveElement( elem, removedElems, removedNodes, removenodes );
1814 ///////////////////////////////////////////////////////////////////////////////
1815 ///@param elem The element to delete
1816 ///@param removedElems to be filled with all removed elements
1817 ///@param removedNodes to be filled with all removed nodes
1818 ///@param removenodes if true remaining nodes will be removed
1819 ///////////////////////////////////////////////////////////////////////////////
1820 void SMDS_Mesh::RemoveElement(const SMDS_MeshElement * elem,
1821 std::vector<const SMDS_MeshElement *>& removedElems,
1822 std::vector<const SMDS_MeshElement *>& removedNodes,
1825 // get finite elements built on elem
1826 std::set<const SMDS_MeshElement*> * s1;
1827 if ( (elem->GetType() == SMDSAbs_0DElement)
1828 || (elem->GetType() == SMDSAbs_Ball)
1829 || (elem->GetType() == SMDSAbs_Edge)
1830 || (elem->GetType() == SMDSAbs_Face)
1831 || (elem->GetType() == SMDSAbs_Volume) )
1833 s1 = new std::set<const SMDS_MeshElement*> ();
1837 s1 = getFinitElements(elem);
1839 // get exclusive nodes (which would become free afterwards)
1840 std::set<const SMDS_MeshElement*> * s2;
1841 if (elem->GetType() == SMDSAbs_Node) // a node is removed
1843 // do not remove nodes except elem
1844 s2 = new std::set<const SMDS_MeshElement*> ();
1849 s2 = getExclusiveNodes(*s1);
1851 // form the set of finite and construction elements to remove
1852 std::set<const SMDS_MeshElement*> s3;
1853 std::set<const SMDS_MeshElement*>::iterator it = s1->begin();
1854 while (it != s1->end())
1856 addChildrenWithNodes(s3, *it, *s2);
1860 if (elem->GetType() != SMDSAbs_Node)
1863 // remove finite and construction elements
1864 for( it = s3.begin();it != s3.end(); ++it )
1866 // Remove element from <InverseElements> of its nodes
1867 SMDS_NodeIteratorPtr itn = (*it)->nodeIterator();
1870 SMDS_MeshNode * n = const_cast<SMDS_MeshNode *> (itn->next());
1871 n->RemoveInverseElement((*it));
1874 int vtkid = (*it)->GetVtkID();
1876 switch ((*it)->GetType()) {
1878 throw SALOME_Exception(LOCALIZED("Internal Error: This should not happen"));
1880 case SMDSAbs_Edge: myInfo.RemoveEdge(*it); break;
1881 case SMDSAbs_Face: myInfo.RemoveFace(*it); break;
1882 case SMDSAbs_Volume: myInfo.RemoveVolume(*it); break;
1883 case SMDSAbs_Ball: myInfo.myNbBalls--; break;
1884 case SMDSAbs_0DElement: myInfo.myNb0DElements--; break;
1885 case SMDSAbs_All: // avoid compilation warning
1886 case SMDSAbs_NbElementTypes: break;
1888 removedElems.push_back( *it);
1890 myCellFactory->Free( static_cast< const SMDS_MeshCell*>( *it ));
1894 this->myGrid->GetCellTypesArray()->SetValue(vtkid, VTK_EMPTY_CELL);
1898 // remove exclusive (free) nodes
1901 for ( it = s2->begin(); it != s2->end(); ++it )
1904 myNodeFactory->Free( (*it) );
1905 removedNodes.push_back((*it));
1914 ///////////////////////////////////////////////////////////////////////////////
1915 ///@param elem The element to delete
1916 ///////////////////////////////////////////////////////////////////////////////
1917 void SMDS_Mesh::RemoveFreeElement(const SMDS_MeshElement * elem)
1919 const int vtkId = elem->GetVtkID();
1920 SMDSAbs_ElementType aType = elem->GetType();
1921 if ( aType == SMDSAbs_Node )
1923 // only free node can be removed by this method
1924 const SMDS_MeshNode* n = static_cast<const SMDS_MeshNode*>( elem );
1925 if ( n->NbInverseElements() == 0 ) { // free node
1927 myNodeFactory->Free( n );
1931 throw SALOME_Exception( LOCALIZED( "RemoveFreeElement: not a free node" ));
1936 // Remove element from <InverseElements> of its nodes
1937 SMDS_NodeIteratorPtr itn = elem->nodeIterator();
1938 while (itn->more()) {
1939 SMDS_MeshNode * n = const_cast<SMDS_MeshNode *>(itn->next());
1940 n->RemoveInverseElement(elem);
1943 // in meshes without descendants elements are always free
1945 case SMDSAbs_0DElement: myInfo.remove(elem); break;
1946 case SMDSAbs_Edge: myInfo.RemoveEdge(elem); break;
1947 case SMDSAbs_Face: myInfo.RemoveFace(elem); break;
1948 case SMDSAbs_Volume: myInfo.RemoveVolume(elem); break;
1949 case SMDSAbs_Ball: myInfo.remove(elem); break;
1952 myCellFactory->Free( elem );
1954 this->myGrid->GetCellTypesArray()->SetValue(vtkId, VTK_EMPTY_CELL);
1958 //=======================================================================
1960 * Checks if the element is present in mesh.
1962 //=======================================================================
1964 bool SMDS_Mesh::Contains (const SMDS_MeshElement* elem) const
1966 if ( !elem || elem->IsNull() )
1969 if ( elem->GetType() == SMDSAbs_Node )
1970 return ( elem == myNodeFactory->FindElement( elem->GetID() ));
1972 return ( elem == myCellFactory->FindElement( elem->GetID() ));
1975 //=======================================================================
1976 //function : MaxNodeID
1978 //=======================================================================
1980 int SMDS_Mesh::MaxNodeID() const
1982 return myNodeFactory->GetMaxID();
1985 //=======================================================================
1986 //function : MinNodeID
1988 //=======================================================================
1990 int SMDS_Mesh::MinNodeID() const
1992 return myNodeFactory->GetMinID();
1995 //=======================================================================
1996 //function : MaxElementID
1998 //=======================================================================
2000 int SMDS_Mesh::MaxElementID() const
2002 return myCellFactory->GetMaxID();
2005 //=======================================================================
2006 //function : MinElementID
2008 //=======================================================================
2010 int SMDS_Mesh::MinElementID() const
2012 return myCellFactory->GetMinID();
2015 //=======================================================================
2016 //function : Renumber
2017 //purpose : Renumber all nodes or elements.
2018 //=======================================================================
2020 // void SMDS_Mesh::Renumber (const bool isNodes, const int startID, const int deltaID)
2022 // if ( deltaID == 0 )
2027 //=======================================================================
2028 //function : GetElementType
2029 //purpose : Return type of element or node with id
2030 //=======================================================================
2032 SMDSAbs_ElementType SMDS_Mesh::GetElementType( const int id, const bool iselem ) const
2034 const SMDS_MeshElement* elem = 0;
2036 elem = myCellFactory->FindElement( id );
2038 elem = myNodeFactory->FindElement( id );
2040 return elem ? elem->GetType() : SMDSAbs_All;
2045 //********************************************************************
2046 //********************************************************************
2047 //******** *********
2048 //***** Methods for addition of quadratic elements ******
2049 //******** *********
2050 //********************************************************************
2051 //********************************************************************
2053 //=======================================================================
2054 //function : AddEdgeWithID
2056 //=======================================================================
2057 SMDS_MeshEdge* SMDS_Mesh::AddEdgeWithID(int n1, int n2, int n12, int ID)
2059 return SMDS_Mesh::AddEdgeWithID (myNodeFactory->FindNode(n1),
2060 myNodeFactory->FindNode(n2),
2061 myNodeFactory->FindNode(n12),
2065 //=======================================================================
2066 //function : AddEdge
2068 //=======================================================================
2069 SMDS_MeshEdge* SMDS_Mesh::AddEdge(const SMDS_MeshNode* n1,
2070 const SMDS_MeshNode* n2,
2071 const SMDS_MeshNode* n12)
2073 return SMDS_Mesh::AddEdgeWithID(n1, n2, n12, myCellFactory->GetFreeID());
2076 //=======================================================================
2077 //function : AddEdgeWithID
2079 //=======================================================================
2080 SMDS_MeshEdge* SMDS_Mesh::AddEdgeWithID(const SMDS_MeshNode * n1,
2081 const SMDS_MeshNode * n2,
2082 const SMDS_MeshNode * n12,
2085 if ( !n1 || !n2 || !n12 ) return 0;
2087 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2089 cell->init( SMDSEntity_Quad_Edge, /*nbNodes=*/3, n1, n2, n12 );
2090 myInfo.myNbQuadEdges++;
2091 return static_cast<SMDS_MeshEdge*>( cell );
2097 //=======================================================================
2098 //function : AddFace
2100 //=======================================================================
2101 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshNode * n1,
2102 const SMDS_MeshNode * n2,
2103 const SMDS_MeshNode * n3,
2104 const SMDS_MeshNode * n12,
2105 const SMDS_MeshNode * n23,
2106 const SMDS_MeshNode * n31)
2108 return SMDS_Mesh::AddFaceWithID(n1,n2,n3,n12,n23,n31,
2109 myCellFactory->GetFreeID());
2112 //=======================================================================
2113 //function : AddFaceWithID
2115 //=======================================================================
2116 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(int n1, int n2, int n3,
2117 int n12,int n23,int n31, int ID)
2119 return SMDS_Mesh::AddFaceWithID (myNodeFactory->FindNode(n1) ,
2120 myNodeFactory->FindNode(n2) ,
2121 myNodeFactory->FindNode(n3) ,
2122 myNodeFactory->FindNode(n12),
2123 myNodeFactory->FindNode(n23),
2124 myNodeFactory->FindNode(n31),
2128 //=======================================================================
2129 //function : AddFaceWithID
2131 //=======================================================================
2132 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshNode * n1,
2133 const SMDS_MeshNode * n2,
2134 const SMDS_MeshNode * n3,
2135 const SMDS_MeshNode * n12,
2136 const SMDS_MeshNode * n23,
2137 const SMDS_MeshNode * n31,
2140 if ( !n1 || !n2 || !n3 || !n12 || !n23 || !n31 ) return 0;
2141 if ( NbFaces() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2143 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2145 cell->init( SMDSEntity_Quad_Triangle, /*nbNodes=*/6, n1, n2, n3, n12, n23, n31 );
2146 myInfo.myNbQuadTriangles++;
2147 return static_cast<SMDS_MeshFace*>( cell );
2153 //=======================================================================
2154 //function : AddFace
2156 //=======================================================================
2157 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshNode * n1,
2158 const SMDS_MeshNode * n2,
2159 const SMDS_MeshNode * n3,
2160 const SMDS_MeshNode * n12,
2161 const SMDS_MeshNode * n23,
2162 const SMDS_MeshNode * n31,
2163 const SMDS_MeshNode * nCenter)
2165 return SMDS_Mesh::AddFaceWithID(n1,n2,n3,n12,n23,n31,nCenter,
2166 myCellFactory->GetFreeID());
2169 //=======================================================================
2170 //function : AddFaceWithID
2172 //=======================================================================
2173 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(int n1, int n2, int n3,
2174 int n12,int n23,int n31, int nCenter, int ID)
2176 return SMDS_Mesh::AddFaceWithID (myNodeFactory->FindNode(n1) ,
2177 myNodeFactory->FindNode(n2) ,
2178 myNodeFactory->FindNode(n3) ,
2179 myNodeFactory->FindNode(n12),
2180 myNodeFactory->FindNode(n23),
2181 myNodeFactory->FindNode(n31),
2182 myNodeFactory->FindNode(nCenter),
2186 //=======================================================================
2187 //function : AddFaceWithID
2189 //=======================================================================
2190 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshNode * n1,
2191 const SMDS_MeshNode * n2,
2192 const SMDS_MeshNode * n3,
2193 const SMDS_MeshNode * n12,
2194 const SMDS_MeshNode * n23,
2195 const SMDS_MeshNode * n31,
2196 const SMDS_MeshNode * nCenter,
2199 if ( !n1 || !n2 || !n3 || !n12 || !n23 || !n31 || !nCenter) return 0;
2200 if ( NbFaces() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2202 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2204 cell->init( SMDSEntity_BiQuad_Triangle, /*nbNodes=*/7, n1, n2, n3, n12, n23, n31, nCenter );
2205 myInfo.myNbBiQuadTriangles++;
2206 return static_cast<SMDS_MeshFace*>( cell );
2212 //=======================================================================
2213 //function : AddFace
2215 //=======================================================================
2216 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshNode * n1,
2217 const SMDS_MeshNode * n2,
2218 const SMDS_MeshNode * n3,
2219 const SMDS_MeshNode * n4,
2220 const SMDS_MeshNode * n12,
2221 const SMDS_MeshNode * n23,
2222 const SMDS_MeshNode * n34,
2223 const SMDS_MeshNode * n41)
2225 return SMDS_Mesh::AddFaceWithID(n1,n2,n3,n4,n12,n23,n34,n41,
2226 myCellFactory->GetFreeID());
2229 //=======================================================================
2230 //function : AddFaceWithID
2232 //=======================================================================
2233 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(int n1, int n2, int n3, int n4,
2234 int n12,int n23,int n34,int n41, int ID)
2236 return SMDS_Mesh::AddFaceWithID (myNodeFactory->FindNode(n1) ,
2237 myNodeFactory->FindNode(n2) ,
2238 myNodeFactory->FindNode(n3) ,
2239 myNodeFactory->FindNode(n4) ,
2240 myNodeFactory->FindNode(n12),
2241 myNodeFactory->FindNode(n23),
2242 myNodeFactory->FindNode(n34),
2243 myNodeFactory->FindNode(n41),
2247 //=======================================================================
2248 //function : AddFaceWithID
2250 //=======================================================================
2251 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshNode * n1,
2252 const SMDS_MeshNode * n2,
2253 const SMDS_MeshNode * n3,
2254 const SMDS_MeshNode * n4,
2255 const SMDS_MeshNode * n12,
2256 const SMDS_MeshNode * n23,
2257 const SMDS_MeshNode * n34,
2258 const SMDS_MeshNode * n41,
2261 if ( !n1 || !n2 || !n3 || !n4 || !n12 || !n23 || !n34 || !n41) return 0;
2262 if ( NbFaces() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2264 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2266 cell->init( SMDSEntity_Quad_Quadrangle, /*nbNodes=*/8, n1, n2, n3, n4, n12, n23, n34, n41 );
2267 myInfo.myNbQuadQuadrangles++;
2268 return static_cast<SMDS_MeshFace*>( cell );
2273 //=======================================================================
2274 //function : AddFace
2276 //=======================================================================
2277 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshNode * n1,
2278 const SMDS_MeshNode * n2,
2279 const SMDS_MeshNode * n3,
2280 const SMDS_MeshNode * n4,
2281 const SMDS_MeshNode * n12,
2282 const SMDS_MeshNode * n23,
2283 const SMDS_MeshNode * n34,
2284 const SMDS_MeshNode * n41,
2285 const SMDS_MeshNode * nCenter)
2287 return SMDS_Mesh::AddFaceWithID(n1,n2,n3,n4,n12,n23,n34,n41,nCenter,
2288 myCellFactory->GetFreeID());
2291 //=======================================================================
2292 //function : AddFaceWithID
2294 //=======================================================================
2295 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(int n1, int n2, int n3, int n4,
2296 int n12,int n23,int n34,int n41, int nCenter, int ID)
2298 return SMDS_Mesh::AddFaceWithID (myNodeFactory->FindNode(n1) ,
2299 myNodeFactory->FindNode(n2) ,
2300 myNodeFactory->FindNode(n3) ,
2301 myNodeFactory->FindNode(n4) ,
2302 myNodeFactory->FindNode(n12),
2303 myNodeFactory->FindNode(n23),
2304 myNodeFactory->FindNode(n34),
2305 myNodeFactory->FindNode(n41),
2306 myNodeFactory->FindNode(nCenter),
2310 //=======================================================================
2311 //function : AddFaceWithID
2313 //=======================================================================
2314 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshNode * n1,
2315 const SMDS_MeshNode * n2,
2316 const SMDS_MeshNode * n3,
2317 const SMDS_MeshNode * n4,
2318 const SMDS_MeshNode * n12,
2319 const SMDS_MeshNode * n23,
2320 const SMDS_MeshNode * n34,
2321 const SMDS_MeshNode * n41,
2322 const SMDS_MeshNode * nCenter,
2325 if ( !n1 || !n2 || !n3 || !n4 || !n12 || !n23 || !n34 || !n41 || !nCenter) return 0;
2326 if ( NbFaces() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2328 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2330 cell->init( SMDSEntity_BiQuad_Quadrangle,
2331 /*nbNodes=*/9, n1, n2, n3, n4, n12, n23, n34, n41, nCenter );
2332 myInfo.myNbBiQuadQuadrangles++;
2333 return static_cast<SMDS_MeshFace*>( cell );
2339 //=======================================================================
2340 //function : AddVolume
2342 //=======================================================================
2343 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
2344 const SMDS_MeshNode * n2,
2345 const SMDS_MeshNode * n3,
2346 const SMDS_MeshNode * n4,
2347 const SMDS_MeshNode * n12,
2348 const SMDS_MeshNode * n23,
2349 const SMDS_MeshNode * n31,
2350 const SMDS_MeshNode * n14,
2351 const SMDS_MeshNode * n24,
2352 const SMDS_MeshNode * n34)
2354 return SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n12, n23,
2355 n31, n14, n24, n34, myCellFactory->GetFreeID());
2358 //=======================================================================
2359 //function : AddVolumeWithID
2361 //=======================================================================
2362 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3, int n4,
2363 int n12,int n23,int n31,
2364 int n14,int n24,int n34, int ID)
2366 return SMDS_Mesh::AddVolumeWithID (myNodeFactory->FindNode(n1) ,
2367 myNodeFactory->FindNode(n2) ,
2368 myNodeFactory->FindNode(n3) ,
2369 myNodeFactory->FindNode(n4) ,
2370 myNodeFactory->FindNode(n12),
2371 myNodeFactory->FindNode(n23),
2372 myNodeFactory->FindNode(n31),
2373 myNodeFactory->FindNode(n14),
2374 myNodeFactory->FindNode(n24),
2375 myNodeFactory->FindNode(n34),
2379 //=======================================================================
2380 //function : AddVolumeWithID
2381 //purpose : 2d order tetrahedron of 10 nodes
2382 //=======================================================================
2383 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
2384 const SMDS_MeshNode * n2,
2385 const SMDS_MeshNode * n3,
2386 const SMDS_MeshNode * n4,
2387 const SMDS_MeshNode * n12,
2388 const SMDS_MeshNode * n23,
2389 const SMDS_MeshNode * n31,
2390 const SMDS_MeshNode * n14,
2391 const SMDS_MeshNode * n24,
2392 const SMDS_MeshNode * n34,
2395 if ( !n1 || !n2 || !n3 || !n4 || !n12 || !n23 || !n31 || !n14 || !n24 || !n34)
2397 if ( NbVolumes() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2399 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2401 cell->init( SMDSEntity_Quad_Tetra,
2402 /*nbNodes=*/10, n1, n2, n3, n4, n12, n23, n31, n14, n24, n34 );
2403 myInfo.myNbQuadTetras++;
2404 return static_cast<SMDS_MeshVolume*>( cell );
2410 //=======================================================================
2411 //function : AddVolume
2413 //=======================================================================
2414 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
2415 const SMDS_MeshNode * n2,
2416 const SMDS_MeshNode * n3,
2417 const SMDS_MeshNode * n4,
2418 const SMDS_MeshNode * n5,
2419 const SMDS_MeshNode * n12,
2420 const SMDS_MeshNode * n23,
2421 const SMDS_MeshNode * n34,
2422 const SMDS_MeshNode * n41,
2423 const SMDS_MeshNode * n15,
2424 const SMDS_MeshNode * n25,
2425 const SMDS_MeshNode * n35,
2426 const SMDS_MeshNode * n45)
2428 return SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n12, n23, n34, n41,
2429 n15, n25, n35, n45, myCellFactory->GetFreeID());
2432 //=======================================================================
2433 //function : AddVolumeWithID
2435 //=======================================================================
2436 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3, int n4, int n5,
2437 int n12,int n23,int n34,int n41,
2438 int n15,int n25,int n35,int n45, int ID)
2440 return SMDS_Mesh::AddVolumeWithID (myNodeFactory->FindNode(n1) ,
2441 myNodeFactory->FindNode(n2) ,
2442 myNodeFactory->FindNode(n3) ,
2443 myNodeFactory->FindNode(n4) ,
2444 myNodeFactory->FindNode(n5) ,
2445 myNodeFactory->FindNode(n12),
2446 myNodeFactory->FindNode(n23),
2447 myNodeFactory->FindNode(n34),
2448 myNodeFactory->FindNode(n41),
2449 myNodeFactory->FindNode(n15),
2450 myNodeFactory->FindNode(n25),
2451 myNodeFactory->FindNode(n35),
2452 myNodeFactory->FindNode(n45),
2456 //=======================================================================
2457 //function : AddVolumeWithID
2458 //purpose : 2d order pyramid of 13 nodes
2459 //=======================================================================
2460 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
2461 const SMDS_MeshNode * n2,
2462 const SMDS_MeshNode * n3,
2463 const SMDS_MeshNode * n4,
2464 const SMDS_MeshNode * n5,
2465 const SMDS_MeshNode * n12,
2466 const SMDS_MeshNode * n23,
2467 const SMDS_MeshNode * n34,
2468 const SMDS_MeshNode * n41,
2469 const SMDS_MeshNode * n15,
2470 const SMDS_MeshNode * n25,
2471 const SMDS_MeshNode * n35,
2472 const SMDS_MeshNode * n45,
2475 if (!n1 || !n2 || !n3 || !n4 || !n5 || !n12 || !n23 ||
2476 !n34 || !n41 || !n15 || !n25 || !n35 || !n45)
2478 if ( NbVolumes() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2480 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2482 cell->init( SMDSEntity_Quad_Pyramid,
2483 /*nbNodes=*/13, n1, n2, n3, n4, n5, n12, n23, n34, n41, n15, n25, n35, n45);
2484 myInfo.myNbQuadPyramids++;
2485 return static_cast<SMDS_MeshVolume*>( cell );
2491 //=======================================================================
2492 //function : AddVolume
2493 //purpose : 2d order Pentahedron (prism) with 15 nodes
2494 //=======================================================================
2495 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
2496 const SMDS_MeshNode * n2,
2497 const SMDS_MeshNode * n3,
2498 const SMDS_MeshNode * n4,
2499 const SMDS_MeshNode * n5,
2500 const SMDS_MeshNode * n6,
2501 const SMDS_MeshNode * n12,
2502 const SMDS_MeshNode * n23,
2503 const SMDS_MeshNode * n31,
2504 const SMDS_MeshNode * n45,
2505 const SMDS_MeshNode * n56,
2506 const SMDS_MeshNode * n64,
2507 const SMDS_MeshNode * n14,
2508 const SMDS_MeshNode * n25,
2509 const SMDS_MeshNode * n36)
2511 return SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n6, n12, n23, n31,
2512 n45, n56, n64, n14, n25, n36, myCellFactory->GetFreeID());
2515 //=======================================================================
2516 //function : AddVolumeWithID
2517 //purpose : 2d order Pentahedron (prism) with 15 nodes
2518 //=======================================================================
2519 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3,
2520 int n4, int n5, int n6,
2521 int n12,int n23,int n31,
2522 int n45,int n56,int n64,
2523 int n14,int n25,int n36, int ID)
2525 return SMDS_Mesh::AddVolumeWithID (myNodeFactory->FindNode(n1) ,
2526 myNodeFactory->FindNode(n2) ,
2527 myNodeFactory->FindNode(n3) ,
2528 myNodeFactory->FindNode(n4) ,
2529 myNodeFactory->FindNode(n5) ,
2530 myNodeFactory->FindNode(n6) ,
2531 myNodeFactory->FindNode(n12),
2532 myNodeFactory->FindNode(n23),
2533 myNodeFactory->FindNode(n31),
2534 myNodeFactory->FindNode(n45),
2535 myNodeFactory->FindNode(n56),
2536 myNodeFactory->FindNode(n64),
2537 myNodeFactory->FindNode(n14),
2538 myNodeFactory->FindNode(n25),
2539 myNodeFactory->FindNode(n36),
2543 //=======================================================================
2544 //function : AddVolumeWithID
2545 //purpose : 2d order Pentahedron (prism) with 15 nodes
2546 //=======================================================================
2547 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
2548 const SMDS_MeshNode * n2,
2549 const SMDS_MeshNode * n3,
2550 const SMDS_MeshNode * n4,
2551 const SMDS_MeshNode * n5,
2552 const SMDS_MeshNode * n6,
2553 const SMDS_MeshNode * n12,
2554 const SMDS_MeshNode * n23,
2555 const SMDS_MeshNode * n31,
2556 const SMDS_MeshNode * n45,
2557 const SMDS_MeshNode * n56,
2558 const SMDS_MeshNode * n64,
2559 const SMDS_MeshNode * n14,
2560 const SMDS_MeshNode * n25,
2561 const SMDS_MeshNode * n36,
2564 if (!n1 || !n2 || !n3 || !n4 || !n5 || !n6 || !n12 || !n23 ||
2565 !n31 || !n45 || !n56 || !n64 || !n14 || !n25 || !n36)
2567 if ( NbVolumes() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2569 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2571 cell->init( SMDSEntity_Quad_Penta, /*nbNodes=*/15,
2572 n1, n2, n3, n4, n5, n6, n12, n23, n31, n45, n56, n64, n14, n25, n36 );
2573 myInfo.myNbQuadPrisms++;
2574 return static_cast<SMDS_MeshVolume*>( cell );
2579 //=======================================================================
2580 //function : AddVolume
2581 //purpose : 2d order Pentahedron (prism) with 18 nodes
2582 //=======================================================================
2583 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
2584 const SMDS_MeshNode * n2,
2585 const SMDS_MeshNode * n3,
2586 const SMDS_MeshNode * n4,
2587 const SMDS_MeshNode * n5,
2588 const SMDS_MeshNode * n6,
2589 const SMDS_MeshNode * n12,
2590 const SMDS_MeshNode * n23,
2591 const SMDS_MeshNode * n31,
2592 const SMDS_MeshNode * n45,
2593 const SMDS_MeshNode * n56,
2594 const SMDS_MeshNode * n64,
2595 const SMDS_MeshNode * n14,
2596 const SMDS_MeshNode * n25,
2597 const SMDS_MeshNode * n36,
2598 const SMDS_MeshNode * n1245,
2599 const SMDS_MeshNode * n2356,
2600 const SMDS_MeshNode * n1346)
2602 return SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n6, n12, n23, n31,
2603 n45, n56, n64, n14, n25, n36, n1245, n2356, n1346,
2604 myCellFactory->GetFreeID());
2607 //=======================================================================
2608 //function : AddVolumeWithID
2609 //purpose : 2d order Pentahedron (prism) with 18 nodes
2610 //=======================================================================
2611 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3,
2612 int n4, int n5, int n6,
2613 int n12,int n23,int n31,
2614 int n45,int n56,int n64,
2615 int n14,int n25,int n36,
2616 int n1245, int n2356, int n1346, int ID)
2618 return SMDS_Mesh::AddVolumeWithID (myNodeFactory->FindNode(n1) ,
2619 myNodeFactory->FindNode(n2) ,
2620 myNodeFactory->FindNode(n3) ,
2621 myNodeFactory->FindNode(n4) ,
2622 myNodeFactory->FindNode(n5) ,
2623 myNodeFactory->FindNode(n6) ,
2624 myNodeFactory->FindNode(n12),
2625 myNodeFactory->FindNode(n23),
2626 myNodeFactory->FindNode(n31),
2627 myNodeFactory->FindNode(n45),
2628 myNodeFactory->FindNode(n56),
2629 myNodeFactory->FindNode(n64),
2630 myNodeFactory->FindNode(n14),
2631 myNodeFactory->FindNode(n25),
2632 myNodeFactory->FindNode(n36),
2633 myNodeFactory->FindNode(n1245),
2634 myNodeFactory->FindNode(n2356),
2635 myNodeFactory->FindNode(n1346),
2639 //=======================================================================
2640 //function : AddVolumeWithID
2641 //purpose : 2d order Pentahedron (prism) with 18 nodes
2642 //=======================================================================
2643 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
2644 const SMDS_MeshNode * n2,
2645 const SMDS_MeshNode * n3,
2646 const SMDS_MeshNode * n4,
2647 const SMDS_MeshNode * n5,
2648 const SMDS_MeshNode * n6,
2649 const SMDS_MeshNode * n12,
2650 const SMDS_MeshNode * n23,
2651 const SMDS_MeshNode * n31,
2652 const SMDS_MeshNode * n45,
2653 const SMDS_MeshNode * n56,
2654 const SMDS_MeshNode * n64,
2655 const SMDS_MeshNode * n14,
2656 const SMDS_MeshNode * n25,
2657 const SMDS_MeshNode * n36,
2658 const SMDS_MeshNode * n1245,
2659 const SMDS_MeshNode * n2356,
2660 const SMDS_MeshNode * n1346,
2663 //MESSAGE("AddVolumeWithID penta18 "<< ID);
2664 if (!n1 || !n2 || !n3 || !n4 || !n5 || !n6 || !n12 || !n23 ||
2665 !n31 || !n45 || !n56 || !n64 || !n14 || !n25 || !n36 || !n1245 || !n2356 || !n1346)
2667 if ( NbVolumes() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2669 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2671 cell->init( SMDSEntity_BiQuad_Penta, /*nbNodes=*/18, n1, n2, n3, n4, n5, n6,
2672 n12, n23, n31, n45, n56, n64, n14, n25, n36, n1245, n2356, n1346 );
2673 myInfo.myNbBiQuadPrisms++;
2674 return static_cast<SMDS_MeshVolume*>( cell );
2680 //=======================================================================
2681 //function : AddVolume
2683 //=======================================================================
2684 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
2685 const SMDS_MeshNode * n2,
2686 const SMDS_MeshNode * n3,
2687 const SMDS_MeshNode * n4,
2688 const SMDS_MeshNode * n5,
2689 const SMDS_MeshNode * n6,
2690 const SMDS_MeshNode * n7,
2691 const SMDS_MeshNode * n8,
2692 const SMDS_MeshNode * n12,
2693 const SMDS_MeshNode * n23,
2694 const SMDS_MeshNode * n34,
2695 const SMDS_MeshNode * n41,
2696 const SMDS_MeshNode * n56,
2697 const SMDS_MeshNode * n67,
2698 const SMDS_MeshNode * n78,
2699 const SMDS_MeshNode * n85,
2700 const SMDS_MeshNode * n15,
2701 const SMDS_MeshNode * n26,
2702 const SMDS_MeshNode * n37,
2703 const SMDS_MeshNode * n48)
2705 return SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n6, n7, n8, n12, n23, n34, n41,
2706 n56, n67, n78, n85, n15, n26, n37, n48,
2707 myCellFactory->GetFreeID());
2710 //=======================================================================
2711 //function : AddVolumeWithID
2713 //=======================================================================
2714 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3, int n4,
2715 int n5, int n6, int n7, int n8,
2716 int n12,int n23,int n34,int n41,
2717 int n56,int n67,int n78,int n85,
2718 int n15,int n26,int n37,int n48, int ID)
2720 return SMDS_Mesh::AddVolumeWithID (myNodeFactory->FindNode(n1),
2721 myNodeFactory->FindNode(n2),
2722 myNodeFactory->FindNode(n3),
2723 myNodeFactory->FindNode(n4),
2724 myNodeFactory->FindNode(n5),
2725 myNodeFactory->FindNode(n6),
2726 myNodeFactory->FindNode(n7),
2727 myNodeFactory->FindNode(n8),
2728 myNodeFactory->FindNode(n12),
2729 myNodeFactory->FindNode(n23),
2730 myNodeFactory->FindNode(n34),
2731 myNodeFactory->FindNode(n41),
2732 myNodeFactory->FindNode(n56),
2733 myNodeFactory->FindNode(n67),
2734 myNodeFactory->FindNode(n78),
2735 myNodeFactory->FindNode(n85),
2736 myNodeFactory->FindNode(n15),
2737 myNodeFactory->FindNode(n26),
2738 myNodeFactory->FindNode(n37),
2739 myNodeFactory->FindNode(n48),
2743 //=======================================================================
2744 //function : AddVolumeWithID
2745 //purpose : 2d order Hexahedrons with 20 nodes
2746 //=======================================================================
2747 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
2748 const SMDS_MeshNode * n2,
2749 const SMDS_MeshNode * n3,
2750 const SMDS_MeshNode * n4,
2751 const SMDS_MeshNode * n5,
2752 const SMDS_MeshNode * n6,
2753 const SMDS_MeshNode * n7,
2754 const SMDS_MeshNode * n8,
2755 const SMDS_MeshNode * n12,
2756 const SMDS_MeshNode * n23,
2757 const SMDS_MeshNode * n34,
2758 const SMDS_MeshNode * n41,
2759 const SMDS_MeshNode * n56,
2760 const SMDS_MeshNode * n67,
2761 const SMDS_MeshNode * n78,
2762 const SMDS_MeshNode * n85,
2763 const SMDS_MeshNode * n15,
2764 const SMDS_MeshNode * n26,
2765 const SMDS_MeshNode * n37,
2766 const SMDS_MeshNode * n48,
2769 if (!n1 || !n2 || !n3 || !n4 || !n5 || !n6 || !n7 || !n8 || !n12 || !n23 ||
2770 !n34 || !n41 || !n56 || !n67 || !n78 || !n85 || !n15 || !n26 || !n37 || !n48)
2772 if ( NbVolumes() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2774 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2776 cell->init( SMDSEntity_Quad_Hexa, /*nbNodes=*/20, n1, n2, n3, n4, n5, n6, n7, n8,
2777 n12, n23, n34, n41, n56, n67, n78, n85, n15, n26, n37, n48 );
2778 myInfo.myNbQuadHexas++;
2779 return static_cast<SMDS_MeshVolume*>( cell );
2784 //=======================================================================
2785 //function : AddVolume
2787 //=======================================================================
2788 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
2789 const SMDS_MeshNode * n2,
2790 const SMDS_MeshNode * n3,
2791 const SMDS_MeshNode * n4,
2792 const SMDS_MeshNode * n5,
2793 const SMDS_MeshNode * n6,
2794 const SMDS_MeshNode * n7,
2795 const SMDS_MeshNode * n8,
2796 const SMDS_MeshNode * n12,
2797 const SMDS_MeshNode * n23,
2798 const SMDS_MeshNode * n34,
2799 const SMDS_MeshNode * n41,
2800 const SMDS_MeshNode * n56,
2801 const SMDS_MeshNode * n67,
2802 const SMDS_MeshNode * n78,
2803 const SMDS_MeshNode * n85,
2804 const SMDS_MeshNode * n15,
2805 const SMDS_MeshNode * n26,
2806 const SMDS_MeshNode * n37,
2807 const SMDS_MeshNode * n48,
2808 const SMDS_MeshNode * n1234,
2809 const SMDS_MeshNode * n1256,
2810 const SMDS_MeshNode * n2367,
2811 const SMDS_MeshNode * n3478,
2812 const SMDS_MeshNode * n1458,
2813 const SMDS_MeshNode * n5678,
2814 const SMDS_MeshNode * nCenter)
2816 return SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n6, n7, n8, n12, n23, n34, n41,
2817 n56, n67, n78, n85, n15, n26, n37, n48,
2818 n1234, n1256, n2367, n3478, n1458, n5678, nCenter,
2819 myCellFactory->GetFreeID());
2822 //=======================================================================
2823 //function : AddVolumeWithID
2825 //=======================================================================
2826 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3, int n4,
2827 int n5, int n6, int n7, int n8,
2828 int n12,int n23,int n34,int n41,
2829 int n56,int n67,int n78,int n85,
2830 int n15,int n26,int n37,int n48,
2831 int n1234,int n1256,int n2367,int n3478,
2832 int n1458,int n5678,int nCenter, int ID)
2834 return SMDS_Mesh::AddVolumeWithID (myNodeFactory->FindNode(n1),
2835 myNodeFactory->FindNode(n2),
2836 myNodeFactory->FindNode(n3),
2837 myNodeFactory->FindNode(n4),
2838 myNodeFactory->FindNode(n5),
2839 myNodeFactory->FindNode(n6),
2840 myNodeFactory->FindNode(n7),
2841 myNodeFactory->FindNode(n8),
2842 myNodeFactory->FindNode(n12),
2843 myNodeFactory->FindNode(n23),
2844 myNodeFactory->FindNode(n34),
2845 myNodeFactory->FindNode(n41),
2846 myNodeFactory->FindNode(n56),
2847 myNodeFactory->FindNode(n67),
2848 myNodeFactory->FindNode(n78),
2849 myNodeFactory->FindNode(n85),
2850 myNodeFactory->FindNode(n15),
2851 myNodeFactory->FindNode(n26),
2852 myNodeFactory->FindNode(n37),
2853 myNodeFactory->FindNode(n48),
2854 myNodeFactory->FindNode(n1234),
2855 myNodeFactory->FindNode(n1256),
2856 myNodeFactory->FindNode(n2367),
2857 myNodeFactory->FindNode(n3478),
2858 myNodeFactory->FindNode(n1458),
2859 myNodeFactory->FindNode(n5678),
2860 myNodeFactory->FindNode(nCenter),
2864 //=======================================================================
2865 //function : AddVolumeWithID
2866 //purpose : 2d order Hexahedrons with 27 nodes
2867 //=======================================================================
2868 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
2869 const SMDS_MeshNode * n2,
2870 const SMDS_MeshNode * n3,
2871 const SMDS_MeshNode * n4,
2872 const SMDS_MeshNode * n5,
2873 const SMDS_MeshNode * n6,
2874 const SMDS_MeshNode * n7,
2875 const SMDS_MeshNode * n8,
2876 const SMDS_MeshNode * n12,
2877 const SMDS_MeshNode * n23,
2878 const SMDS_MeshNode * n34,
2879 const SMDS_MeshNode * n41,
2880 const SMDS_MeshNode * n56,
2881 const SMDS_MeshNode * n67,
2882 const SMDS_MeshNode * n78,
2883 const SMDS_MeshNode * n85,
2884 const SMDS_MeshNode * n15,
2885 const SMDS_MeshNode * n26,
2886 const SMDS_MeshNode * n37,
2887 const SMDS_MeshNode * n48,
2888 const SMDS_MeshNode * n1234,
2889 const SMDS_MeshNode * n1256,
2890 const SMDS_MeshNode * n2367,
2891 const SMDS_MeshNode * n3478,
2892 const SMDS_MeshNode * n1458,
2893 const SMDS_MeshNode * n5678,
2894 const SMDS_MeshNode * nCenter,
2897 if (!n1 || !n2 || !n3 || !n4 || !n5 || !n6 || !n7 || !n8 || !n12 || !n23 ||
2898 !n34 || !n41 || !n56 || !n67 || !n78 || !n85 || !n15 || !n26 || !n37 || !n48 ||
2899 !n1234 || !n1256 || !n2367 || !n3478 || !n1458 || !n5678 || !nCenter )
2901 if ( NbVolumes() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2903 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2905 cell->init( SMDSEntity_TriQuad_Hexa, /*nbNodes=*/27, n1, n2, n3, n4, n5, n6, n7, n8,
2906 n12, n23, n34, n41, n56, n67, n78, n85, n15, n26, n37, n48,
2907 n1234, n1256, n2367, n3478, n1458, n5678, nCenter);
2908 myInfo.myNbTriQuadHexas++;
2909 return static_cast<SMDS_MeshVolume*>( cell );
2914 void SMDS_Mesh::dumpGrid(std::string ficdump)
2916 // vtkUnstructuredGridWriter* aWriter = vtkUnstructuredGridWriter::New();
2917 // aWriter->SetFileName(ficdump.c_str());
2918 // aWriter->SetInput(myGrid);
2919 // if(myGrid->GetNumberOfCells())
2921 // aWriter->Write();
2923 // aWriter->Delete();
2924 ficdump = ficdump + "_connectivity";
2925 std::ofstream ficcon(ficdump.c_str(), ios::out);
2926 int nbPoints = myGrid->GetNumberOfPoints();
2927 ficcon << "-------------------------------- points " << nbPoints << endl;
2928 for (int i=0; i<nbPoints; i++)
2930 ficcon << i << " " << *(myGrid->GetPoint(i)) << " " << *(myGrid->GetPoint(i)+1) << " " << " " << *(myGrid->GetPoint(i)+2) << endl;
2932 int nbCells = myGrid->GetNumberOfCells();
2933 ficcon << "-------------------------------- cells " << nbCells << endl;
2934 for (int i=0; i<nbCells; i++)
2936 ficcon << i << " - " << myGrid->GetCell(i)->GetCellType() << " -";
2937 int nbptcell = myGrid->GetCell(i)->GetNumberOfPoints();
2938 vtkIdList *listid = myGrid->GetCell(i)->GetPointIds();
2939 for (int j=0; j<nbptcell; j++)
2941 ficcon << " " << listid->GetId(j);
2945 ficcon << "-------------------------------- connectivity " << nbPoints << endl;
2946 vtkCellLinks *links = myGrid->GetLinks();
2947 for (int i=0; i<nbPoints; i++)
2949 int ncells = links->GetNcells(i);
2950 vtkIdType *cells = links->GetCells(i);
2951 ficcon << i << " - " << ncells << " -";
2952 for (int j=0; j<ncells; j++)
2954 ficcon << " " << cells[j];
2962 void SMDS_Mesh::CompactMesh()
2964 this->myCompactTime = this->myModifTime;
2966 bool idsChange = HasNumerationHoles();
2969 std::set< SMDS_ElementHolder* >::iterator holder = myElemHolders.begin();
2970 for ( ; holder != myElemHolders.end(); ++holder )
2971 (*holder)->beforeCompacting();
2973 int oldCellSize = myCellFactory->GetMaxID();
2975 // remove "holes" in SMDS numeration
2976 std::vector<int> idNodesOldToNew, idCellsNewToOld, idCellsOldToNew;
2977 myNodeFactory->Compact( idNodesOldToNew );
2978 myCellFactory->Compact( idCellsNewToOld );
2980 // make VTK IDs correspond to SMDS IDs
2981 int newNodeSize = myNodeFactory->NbUsedElements();
2982 int newCellSize = myCellFactory->NbUsedElements();
2983 myGrid->compactGrid( idNodesOldToNew, newNodeSize, idCellsNewToOld, newCellSize );
2985 if ( idsChange && !myElemHolders.empty() )
2987 // idCellsNewToOld -> idCellsOldToNew
2988 idCellsOldToNew.resize( oldCellSize, oldCellSize );
2989 for ( size_t iNew = 0; iNew < idCellsNewToOld.size(); ++iNew )
2991 if ( idCellsNewToOld[ iNew ] >= (int) idCellsOldToNew.size() )
2992 idCellsOldToNew.resize( ( 1 + idCellsNewToOld[ iNew ]) * 1.5, oldCellSize );
2993 idCellsOldToNew[ idCellsNewToOld[ iNew ]] = iNew;
2997 std::set< SMDS_ElementHolder* >::iterator holder = myElemHolders.begin();
2998 for ( ; holder != myElemHolders.end(); ++holder )
3000 (*holder)->restoreElements( idNodesOldToNew, idCellsOldToNew );
3002 (*holder)->compact();
3007 int SMDS_Mesh::FromVtkToSmds( int vtkid ) const
3009 return myCellFactory->FromVtkToSmds( vtkid );
3012 double SMDS_Mesh::getMaxDim()
3014 double dmax = 1.e-3;
3015 if ((xmax - xmin) > dmax) dmax = xmax -xmin;
3016 if ((ymax - ymin) > dmax) dmax = ymax -ymin;
3017 if ((zmax - zmin) > dmax) dmax = zmax -zmin;
3021 //! modification that needs compact structure and redraw
3022 void SMDS_Mesh::Modified()
3024 if (this->myModified)
3027 this->myModifTime++;
3032 //! get last modification timeStamp
3033 vtkMTimeType SMDS_Mesh::GetMTime() const
3035 return this->myModifTime;
3038 bool SMDS_Mesh::IsCompacted()
3040 return ( this->myCompactTime == this->myModifTime );
3043 //! are there holes in elements or nodes numeration
3044 bool SMDS_Mesh::HasNumerationHoles()
3046 return ( myNodeFactory->CompactChangePointers() ||
3047 myCellFactory->CompactChangePointers() );
3050 void SMDS_Mesh::setNbShapes( size_t nbShapes )
3052 myNodeFactory->SetNbShapes( nbShapes );