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) throw (std::bad_alloc)
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();
834 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
836 cell->init( SMDSEntity_Polygon, nodes );
837 myInfo.myNbPolygons++;
838 return static_cast<SMDS_MeshFace*>( cell );
843 ///////////////////////////////////////////////////////////////////////////////
844 /// Add a polygon defined by its nodes.
845 /// An ID is automatically affected to the created face.
846 ///////////////////////////////////////////////////////////////////////////////
848 SMDS_MeshFace* SMDS_Mesh::AddPolygonalFace (const std::vector<const SMDS_MeshNode*> & nodes)
850 return SMDS_Mesh::AddPolygonalFaceWithID(nodes, myCellFactory->GetFreeID());
853 ///////////////////////////////////////////////////////////////////////////////
854 /// Add a quadratic polygon defined by its nodes IDs
855 ///////////////////////////////////////////////////////////////////////////////
857 SMDS_MeshFace* SMDS_Mesh::AddQuadPolygonalFaceWithID (const std::vector<int> & nodes_ids,
860 std::vector<const SMDS_MeshNode*> nodes( nodes_ids.size() );
861 for ( size_t i = 0; i < nodes.size(); i++) {
862 nodes[i] = myNodeFactory->FindNode(nodes_ids[i]);
863 if (!nodes[i]) return NULL;
865 return SMDS_Mesh::AddQuadPolygonalFaceWithID(nodes, ID);
868 ///////////////////////////////////////////////////////////////////////////////
869 /// Add a quadratic polygon defined by its nodes
870 ///////////////////////////////////////////////////////////////////////////////
873 SMDS_Mesh::AddQuadPolygonalFaceWithID (const std::vector<const SMDS_MeshNode*> & nodes,
876 if ( NbFaces() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
877 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
879 cell->init( SMDSEntity_Quad_Polygon, nodes );
880 myInfo.myNbQuadPolygons++;
881 return static_cast<SMDS_MeshFace*>( cell );
886 ///////////////////////////////////////////////////////////////////////////////
887 /// Add a quadratic polygon defined by its nodes.
888 /// An ID is automatically affected to the created face.
889 ///////////////////////////////////////////////////////////////////////////////
891 SMDS_MeshFace* SMDS_Mesh::AddQuadPolygonalFace (const std::vector<const SMDS_MeshNode*> & nodes)
893 return SMDS_Mesh::AddQuadPolygonalFaceWithID(nodes, myCellFactory->GetFreeID());
896 ///////////////////////////////////////////////////////////////////////////////
897 /// Create a new polyhedral volume and add it to the mesh.
898 /// @param ID The ID of the new volume
899 /// @return The created volume or NULL if an element with this ID already exists
900 /// or if input nodes are not found.
901 ///////////////////////////////////////////////////////////////////////////////
903 SMDS_MeshVolume * SMDS_Mesh::AddPolyhedralVolumeWithID (const std::vector<int> & nodes_ids,
904 const std::vector<int> & quantities,
907 int nbNodes = nodes_ids.size();
908 std::vector<const SMDS_MeshNode*> nodes (nbNodes);
909 for (int i = 0; i < nbNodes; i++) {
910 nodes[i] = myNodeFactory->FindNode(nodes_ids[i]);
911 if (!nodes[i]) return NULL;
913 return SMDS_Mesh::AddPolyhedralVolumeWithID(nodes, quantities, ID);
916 ///////////////////////////////////////////////////////////////////////////////
917 /// Create a new polyhedral volume and add it to the mesh.
918 /// @param ID The ID of the new volume
919 /// @return The created volume
920 ///////////////////////////////////////////////////////////////////////////////
923 SMDS_Mesh::AddPolyhedralVolumeWithID (const std::vector<const SMDS_MeshNode*>& nodes,
924 const std::vector<int> & quantities,
927 if ( nodes.empty() || quantities.empty() )
929 if ( NbVolumes() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
931 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
933 SMDS_MeshVolume* volume = static_cast<SMDS_MeshVolume*>( cell );
934 volume->init( nodes, quantities );
935 myInfo.myNbPolyhedrons++;
941 ///////////////////////////////////////////////////////////////////////////////
942 /// Create a new polyhedral volume and add it to the mesh.
943 /// @return The created volume
944 ///////////////////////////////////////////////////////////////////////////////
946 SMDS_MeshVolume* SMDS_Mesh::AddPolyhedralVolume
947 (const std::vector<const SMDS_MeshNode*> & nodes,
948 const std::vector<int> & quantities)
950 int ID = myCellFactory->GetFreeID();
951 return SMDS_Mesh::AddPolyhedralVolumeWithID(nodes, quantities, ID);
954 SMDS_MeshVolume* SMDS_Mesh::AddVolumeFromVtkIds(const std::vector<vtkIdType>& vtkNodeIds)
956 SMDS_MeshCell* cell = myCellFactory->NewCell( myCellFactory->GetFreeID() );
957 SMDS_MeshVolume * vol = static_cast<SMDS_MeshVolume*>( cell );
958 vol->init( vtkNodeIds );
963 SMDS_MeshFace* SMDS_Mesh::AddFaceFromVtkIds(const std::vector<vtkIdType>& vtkNodeIds)
965 SMDS_MeshCell* cell = myCellFactory->NewCell( myCellFactory->GetFreeID() );
966 SMDS_MeshFace * f = static_cast<SMDS_MeshFace*>( cell );
967 f->init( vtkNodeIds );
972 //=======================================================================
973 //function : MoveNode
975 //=======================================================================
977 void SMDS_Mesh::MoveNode(const SMDS_MeshNode *n, double x, double y, double z)
979 SMDS_MeshNode * node=const_cast<SMDS_MeshNode*>(n);
983 ///////////////////////////////////////////////////////////////////////////////
984 /// Return the node whose SMDS ID is 'ID'.
985 ///////////////////////////////////////////////////////////////////////////////
986 const SMDS_MeshNode * SMDS_Mesh::FindNode(int ID) const
988 return myNodeFactory->FindNode( ID );
991 ///////////////////////////////////////////////////////////////////////////////
992 /// Return the node whose VTK ID is 'vtkId'.
993 ///////////////////////////////////////////////////////////////////////////////
994 const SMDS_MeshNode * SMDS_Mesh::FindNodeVtk(int vtkId) const
996 return myNodeFactory->FindNode( vtkId + 1 );
999 const SMDS_MeshElement * SMDS_Mesh::FindElementVtk(int IDelem) const
1001 return myCellFactory->FindElement( FromVtkToSmds( IDelem ));
1004 ///////////////////////////////////////////////////////////////////////////////
1005 /// Remove a node and all the elements which own this node
1006 ///////////////////////////////////////////////////////////////////////////////
1008 void SMDS_Mesh::RemoveNode(const SMDS_MeshNode * node)
1010 RemoveElement(node, true);
1013 //=======================================================================
1014 //function : RemoveFromParent
1016 //=======================================================================
1018 bool SMDS_Mesh::RemoveFromParent()
1020 if (myParent==NULL) return false;
1021 else return (myParent->RemoveSubMesh(this));
1024 //=======================================================================
1025 //function : RemoveSubMesh
1027 //=======================================================================
1029 bool SMDS_Mesh::RemoveSubMesh(const SMDS_Mesh * aMesh)
1033 std::list<SMDS_Mesh *>::iterator itmsh=myChildren.begin();
1034 for (; itmsh!=myChildren.end() && !found; itmsh++)
1036 SMDS_Mesh * submesh = *itmsh;
1037 if (submesh == aMesh)
1040 myChildren.erase(itmsh);
1047 //=======================================================================
1048 //function : ChangePolyhedronNodes
1050 //=======================================================================
1052 bool SMDS_Mesh::ChangePolyhedronNodes(const SMDS_MeshElement * element,
1053 const std::vector<const SMDS_MeshNode*>& nodes,
1054 const std::vector<int>& quantities)
1056 // keep current nodes of element
1057 std::set<const SMDS_MeshNode*> oldNodes( element->begin_nodes(), element->end_nodes() );
1061 if ( const SMDS_MeshVolume* vol = DownCast<SMDS_MeshVolume>( element ))
1062 Ok = vol->ChangeNodes( nodes, quantities );
1067 updateInverseElements( element, &nodes[0], nodes.size(), oldNodes );
1072 //=======================================================================
1073 //function : ChangeElementNodes
1075 //=======================================================================
1077 bool SMDS_Mesh::ChangeElementNodes(const SMDS_MeshElement * element,
1078 const SMDS_MeshNode * nodes[],
1081 // keep current nodes of element
1082 std::set<const SMDS_MeshNode*> oldNodes( element->begin_nodes(), element->end_nodes() );
1086 if ( SMDS_MeshCell* cell = dynamic_cast<SMDS_MeshCell*>((SMDS_MeshElement*) element))
1087 Ok = cell->ChangeNodes(nodes, nbnodes);
1092 updateInverseElements( element, nodes, nbnodes, oldNodes );
1097 //=======================================================================
1098 //function : updateInverseElements
1099 //purpose : update InverseElements when element changes node
1100 //=======================================================================
1102 void SMDS_Mesh::updateInverseElements( const SMDS_MeshElement * element,
1103 const SMDS_MeshNode* const* nodes,
1105 std::set<const SMDS_MeshNode*>& oldNodes )
1107 if ( GetGrid()->HasLinks() ) // update InverseElements
1109 std::set<const SMDS_MeshNode*>::iterator it;
1111 // AddInverseElement to new nodes
1112 for ( int i = 0; i < nbnodes; i++ )
1114 it = oldNodes.find( nodes[i] );
1115 if ( it == oldNodes.end() )
1117 const_cast<SMDS_MeshNode*>( nodes[i] )->AddInverseElement( element );
1119 // remove from oldNodes a node that remains in elem
1120 oldNodes.erase( it );
1122 // RemoveInverseElement from the nodes removed from elem
1123 for ( it = oldNodes.begin(); it != oldNodes.end(); it++ )
1125 SMDS_MeshNode * n = const_cast<SMDS_MeshNode *>( *it );
1126 n->RemoveInverseElement( element );
1132 const SMDS_Mesh0DElement* SMDS_Mesh::Find0DElement(const SMDS_MeshNode * node)
1134 if (!node) return 0;
1135 const SMDS_Mesh0DElement* toReturn = NULL;
1136 SMDS_ElemIteratorPtr it1 = node->GetInverseElementIterator(SMDSAbs_0DElement);
1137 while (it1->more() && (toReturn == NULL)) {
1138 const SMDS_MeshElement* e = it1->next();
1139 if (e->NbNodes() == 1) {
1140 toReturn = static_cast<const SMDS_Mesh0DElement*>(e);
1146 const SMDS_BallElement* SMDS_Mesh::FindBall(const SMDS_MeshNode * node)
1148 if (!node) return 0;
1149 const SMDS_BallElement* toReturn = NULL;
1150 SMDS_ElemIteratorPtr it1 = node->GetInverseElementIterator(SMDSAbs_Ball);
1151 while (it1->more() && (toReturn == NULL)) {
1152 const SMDS_MeshElement* e = it1->next();
1153 if (e->GetGeomType() == SMDSGeom_BALL)
1154 toReturn = static_cast<const SMDS_BallElement*>(e);
1159 const SMDS_MeshEdge* SMDS_Mesh::FindEdge(const SMDS_MeshNode * node1,
1160 const SMDS_MeshNode * node2)
1162 if ( !node1 ) return 0;
1163 const SMDS_MeshEdge * toReturn=NULL;
1164 SMDS_ElemIteratorPtr it1=node1->GetInverseElementIterator(SMDSAbs_Edge);
1165 while(it1->more()) {
1166 const SMDS_MeshElement * e = it1->next();
1167 if ( e->NbNodes() == 2 && e->GetNodeIndex( node2 ) >= 0 ) {
1168 toReturn = static_cast<const SMDS_MeshEdge*>( e );
1175 const SMDS_MeshEdge* SMDS_Mesh::FindEdge(const SMDS_MeshNode * node1,
1176 const SMDS_MeshNode * node2,
1177 const SMDS_MeshNode * node3)
1179 if ( !node1 ) return 0;
1180 SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Edge);
1181 while(it1->more()) {
1182 const SMDS_MeshElement * e = it1->next();
1183 if ( e->NbNodes() == 3 ) {
1184 SMDS_ElemIteratorPtr it2 = e->nodesIterator();
1185 while(it2->more()) {
1186 const SMDS_MeshElement* n = it2->next();
1196 return static_cast<const SMDS_MeshEdge *> (e);
1202 //=======================================================================
1203 //function : FindFace
1205 //=======================================================================
1207 const SMDS_MeshFace* SMDS_Mesh::FindFace(const SMDS_MeshNode *node1,
1208 const SMDS_MeshNode *node2,
1209 const SMDS_MeshNode *node3)
1211 if ( !node1 ) return 0;
1212 SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Face);
1213 while(it1->more()) {
1214 const SMDS_MeshElement * e = it1->next();
1215 if ( e->NbNodes() == 3 ) {
1216 SMDS_ElemIteratorPtr it2 = e->nodesIterator();
1217 while(it2->more()) {
1218 const SMDS_MeshElement* n = it2->next();
1228 return static_cast<const SMDS_MeshFace *> (e);
1234 //=======================================================================
1235 //function : FindFace
1237 //=======================================================================
1239 const SMDS_MeshFace* SMDS_Mesh::FindFace(const SMDS_MeshNode *node1,
1240 const SMDS_MeshNode *node2,
1241 const SMDS_MeshNode *node3,
1242 const SMDS_MeshNode *node4)
1244 if ( !node1 ) return 0;
1245 SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Face);
1246 while(it1->more()) {
1247 const SMDS_MeshElement * e = it1->next();
1248 if ( e->NbNodes() == 4 ) {
1249 SMDS_ElemIteratorPtr it2 = e->nodesIterator();
1250 while(it2->more()) {
1251 const SMDS_MeshElement* n = it2->next();
1262 return static_cast<const SMDS_MeshFace *> (e);
1268 //=======================================================================
1269 //function : FindFace
1270 //purpose :quadratic triangle
1271 //=======================================================================
1273 const SMDS_MeshFace* SMDS_Mesh::FindFace(const SMDS_MeshNode *node1,
1274 const SMDS_MeshNode *node2,
1275 const SMDS_MeshNode *node3,
1276 const SMDS_MeshNode *node4,
1277 const SMDS_MeshNode *node5,
1278 const SMDS_MeshNode *node6)
1280 if ( !node1 ) return 0;
1281 SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Face);
1282 while(it1->more()) {
1283 const SMDS_MeshElement * e = it1->next();
1284 if ( e->NbNodes() == 6 ) {
1285 SMDS_ElemIteratorPtr it2 = e->nodesIterator();
1286 while(it2->more()) {
1287 const SMDS_MeshElement* n = it2->next();
1300 return static_cast<const SMDS_MeshFace *> (e);
1307 //=======================================================================
1308 //function : FindFace
1309 //purpose : quadratic quadrangle
1310 //=======================================================================
1312 const SMDS_MeshFace* SMDS_Mesh::FindFace(const SMDS_MeshNode *node1,
1313 const SMDS_MeshNode *node2,
1314 const SMDS_MeshNode *node3,
1315 const SMDS_MeshNode *node4,
1316 const SMDS_MeshNode *node5,
1317 const SMDS_MeshNode *node6,
1318 const SMDS_MeshNode *node7,
1319 const SMDS_MeshNode *node8)
1321 if ( !node1 ) return 0;
1322 SMDS_ElemIteratorPtr it1 = node1->GetInverseElementIterator(SMDSAbs_Face);
1323 while(it1->more()) {
1324 const SMDS_MeshElement * e = it1->next();
1325 if ( e->NbNodes() == 8 ) {
1326 SMDS_ElemIteratorPtr it2 = e->nodesIterator();
1327 while(it2->more()) {
1328 const SMDS_MeshElement* n = it2->next();
1343 return static_cast<const SMDS_MeshFace *> (e);
1350 //=======================================================================
1351 //function : FindElement
1353 //=======================================================================
1355 const SMDS_MeshElement* SMDS_Mesh::FindElement(int IDelem) const
1357 return myCellFactory->FindElement( IDelem );
1360 //=======================================================================
1361 //function : FindFace
1362 //purpose : find polygon
1363 //=======================================================================
1366 const SMDS_MeshFace* SMDS_Mesh::FindFace (const std::vector<const SMDS_MeshNode *>& nodes)
1368 return (const SMDS_MeshFace*) FindElement( nodes, SMDSAbs_Face );
1372 //================================================================================
1374 * \brief Return element based on all given nodes
1375 * \param nodes - node of element
1376 * \param type - type of element
1377 * \param noMedium - true if medium nodes of quadratic element are not included in <nodes>
1378 * \retval const SMDS_MeshElement* - found element or NULL
1380 //================================================================================
1382 const SMDS_MeshElement* SMDS_Mesh::FindElement (const std::vector<const SMDS_MeshNode *>& nodes,
1383 const SMDSAbs_ElementType type,
1384 const bool noMedium)
1386 if ( nodes.size() > 0 && nodes[0] )
1388 SMDS_ElemIteratorPtr itF = nodes[0]->GetInverseElementIterator(type);
1391 const SMDS_MeshElement* e = itF->next();
1392 int nbNodesToCheck = noMedium ? e->NbCornerNodes() : e->NbNodes();
1393 if ( nbNodesToCheck == (int)nodes.size() )
1395 for ( size_t i = 1; e && i < nodes.size(); ++i )
1397 int nodeIndex = e->GetNodeIndex( nodes[ i ]);
1398 if ( nodeIndex < 0 || nodeIndex >= nbNodesToCheck )
1409 //================================================================================
1411 * \brief Return elements including all given nodes
1412 * \param [in] nodes - nodes to find elements around
1413 * \param [out] foundElems - the found elements
1414 * \param [in] type - type of elements to find
1415 * \return int - a number of found elements
1417 //================================================================================
1419 int SMDS_Mesh::GetElementsByNodes(const std::vector<const SMDS_MeshNode *>& nodes,
1420 std::vector<const SMDS_MeshElement *>& foundElems,
1421 const SMDSAbs_ElementType type)
1423 // chose a node with minimal number of inverse elements
1424 const SMDS_MeshNode* n0 = nodes[0];
1425 int minNbInverse = n0 ? n0->NbInverseElements( type ) : 1000;
1426 for ( size_t i = 1; i < nodes.size(); ++i )
1427 if ( nodes[i] && nodes[i]->NbInverseElements( type ) < minNbInverse )
1430 minNbInverse = n0->NbInverseElements( type );
1436 foundElems.reserve( minNbInverse );
1437 SMDS_ElemIteratorPtr eIt = n0->GetInverseElementIterator( type );
1438 while ( eIt->more() )
1440 const SMDS_MeshElement* e = eIt->next();
1441 bool includeAll = true;
1442 for ( size_t i = 0; i < nodes.size() && includeAll; ++i )
1443 if ( nodes[i] != n0 && e->GetNodeIndex( nodes[i] ) < 0 )
1446 foundElems.push_back( e );
1449 return foundElems.size();
1452 ///////////////////////////////////////////////////////////////////////////////
1453 /// Return the number of nodes
1454 ///////////////////////////////////////////////////////////////////////////////
1455 int SMDS_Mesh::NbNodes() const
1457 return myInfo.NbNodes();
1460 ///////////////////////////////////////////////////////////////////////////////
1461 /// Return the number of elements
1462 ///////////////////////////////////////////////////////////////////////////////
1463 int SMDS_Mesh::NbElements() const
1465 return myInfo.NbElements();
1467 ///////////////////////////////////////////////////////////////////////////////
1468 /// Return the number of 0D elements
1469 ///////////////////////////////////////////////////////////////////////////////
1470 int SMDS_Mesh::Nb0DElements() const
1472 return myInfo.Nb0DElements();
1475 ///////////////////////////////////////////////////////////////////////////////
1476 /// Return the number of 0D elements
1477 ///////////////////////////////////////////////////////////////////////////////
1478 int SMDS_Mesh::NbBalls() const
1480 return myInfo.NbBalls();
1483 ///////////////////////////////////////////////////////////////////////////////
1484 /// Return the number of edges (including construction edges)
1485 ///////////////////////////////////////////////////////////////////////////////
1486 int SMDS_Mesh::NbEdges() const
1488 return myInfo.NbEdges();
1491 ///////////////////////////////////////////////////////////////////////////////
1492 /// Return the number of faces (including construction faces)
1493 ///////////////////////////////////////////////////////////////////////////////
1494 int SMDS_Mesh::NbFaces() const
1496 return myInfo.NbFaces();
1499 ///////////////////////////////////////////////////////////////////////////////
1500 /// Return the number of volumes
1501 ///////////////////////////////////////////////////////////////////////////////
1502 int SMDS_Mesh::NbVolumes() const
1504 return myInfo.NbVolumes();
1507 ///////////////////////////////////////////////////////////////////////////////
1508 /// Return the number of child mesh of this mesh.
1509 /// Note that the tree structure of SMDS_Mesh is unused in SMESH
1510 ///////////////////////////////////////////////////////////////////////////////
1511 int SMDS_Mesh::NbSubMesh() const
1513 return myChildren.size();
1516 ///////////////////////////////////////////////////////////////////////////////
1517 /// Destroy the mesh and all its elements
1518 /// All pointer on elements owned by this mesh become illegals.
1519 ///////////////////////////////////////////////////////////////////////////////
1520 SMDS_Mesh::~SMDS_Mesh()
1522 std::list<SMDS_Mesh*>::iterator itc=myChildren.begin();
1523 while(itc!=myChildren.end())
1529 delete myNodeFactory;
1530 delete myCellFactory;
1535 //================================================================================
1537 * \brief Clear all data
1539 //================================================================================
1541 void SMDS_Mesh::Clear()
1543 std::set< SMDS_ElementHolder* >::iterator holder = myElemHolders.begin();
1544 for ( ; holder != myElemHolders.end(); ++holder )
1547 myNodeFactory->Clear();
1548 myCellFactory->Clear();
1550 std::list<SMDS_Mesh*>::iterator itc=myChildren.begin();
1551 while(itc!=myChildren.end())
1562 myGrid->Initialize();
1564 vtkPoints* points = vtkPoints::New();
1565 // rnv: to fix bug "21125: EDF 1233 SMESH: Degrardation of precision in a test case for quadratic conversion"
1566 // using double type for storing coordinates of nodes instead float.
1567 points->SetDataType(VTK_DOUBLE);
1568 points->SetNumberOfPoints( 0 );
1569 myGrid->SetPoints( points );
1571 myGrid->DeleteLinks();
1574 ///////////////////////////////////////////////////////////////////////////////
1575 /// Return an iterator on nodes of the current mesh factory
1576 ///////////////////////////////////////////////////////////////////////////////
1578 SMDS_NodeIteratorPtr SMDS_Mesh::nodesIterator() const
1580 return myNodeFactory->GetIterator< SMDS_NodeIterator >( new SMDS_MeshElement::NonNullFilter );
1583 SMDS_ElemIteratorPtr SMDS_Mesh::elementGeomIterator(SMDSAbs_GeometryType type) const
1585 int nbElems = myCellFactory->CompactChangePointers() ? -1 : myInfo.NbElements( type );
1586 return myCellFactory->GetIterator< SMDS_ElemIterator >( new SMDS_MeshElement::GeomFilter( type ),
1590 SMDS_ElemIteratorPtr SMDS_Mesh::elementEntityIterator(SMDSAbs_EntityType type) const
1592 if ( type == SMDSEntity_Node )
1594 return myNodeFactory->GetIterator< SMDS_ElemIterator >( new SMDS_MeshElement::NonNullFilter );
1596 int nbElems = myCellFactory->CompactChangePointers() ? -1 : myInfo.NbElements( type );
1597 return myCellFactory->GetIterator<SMDS_ElemIterator>( new SMDS_MeshElement::EntityFilter( type ),
1601 ///////////////////////////////////////////////////////////////////////////////
1602 /// Return an iterator on elements of the current mesh factory
1603 ///////////////////////////////////////////////////////////////////////////////
1604 SMDS_ElemIteratorPtr SMDS_Mesh::elementsIterator(SMDSAbs_ElementType type) const
1606 typedef SMDS_ElemIterator TIterator;
1610 return myCellFactory->GetIterator< TIterator >( new SMDS_MeshElement::NonNullFilter );
1613 return myNodeFactory->GetIterator< TIterator >( new SMDS_MeshElement::NonNullFilter );
1616 int nbElems = myCellFactory->CompactChangePointers() ? -1 : myInfo.NbElements( type );
1617 return myCellFactory->GetIterator< TIterator >( new SMDS_MeshElement::TypeFilter( type ),
1620 return SMDS_ElemIteratorPtr();
1623 ///////////////////////////////////////////////////////////////////////////////
1624 ///Return an iterator on edges of the current mesh.
1625 ///////////////////////////////////////////////////////////////////////////////
1627 SMDS_EdgeIteratorPtr SMDS_Mesh::edgesIterator() const
1629 typedef SMDS_EdgeIterator TIterator;
1630 int nbElems = myCellFactory->CompactChangePointers() ? -1 : myInfo.NbEdges();
1631 return myCellFactory->GetIterator< TIterator >( new SMDS_MeshElement::TypeFilter( SMDSAbs_Edge ),
1635 ///////////////////////////////////////////////////////////////////////////////
1636 ///Return an iterator on faces of the current mesh.
1637 ///////////////////////////////////////////////////////////////////////////////
1639 SMDS_FaceIteratorPtr SMDS_Mesh::facesIterator() const
1641 typedef SMDS_FaceIterator TIterator;
1642 int nbElems = myCellFactory->CompactChangePointers() ? -1 : myInfo.NbFaces();
1643 return myCellFactory->GetIterator< TIterator >( new SMDS_MeshElement::TypeFilter( SMDSAbs_Face ),
1647 ///////////////////////////////////////////////////////////////////////////////
1648 ///Return an iterator on volumes of the current mesh.
1649 ///////////////////////////////////////////////////////////////////////////////
1651 SMDS_VolumeIteratorPtr SMDS_Mesh::volumesIterator() const
1653 typedef SMDS_VolumeIterator TIterator;
1654 int nbElems = myCellFactory->CompactChangePointers() ? -1 : myInfo.NbVolumes();
1656 myCellFactory->GetIterator< TIterator >( new SMDS_MeshElement::TypeFilter( SMDSAbs_Volume ),
1660 SMDS_NodeIteratorPtr SMDS_Mesh::shapeNodesIterator(int shapeID,
1661 size_t nbElemsToReturn,
1662 const SMDS_MeshNode* sm1stNode) const
1664 return myNodeFactory->GetShapeIterator< SMDS_NodeIterator >( shapeID, nbElemsToReturn, sm1stNode );
1667 SMDS_ElemIteratorPtr SMDS_Mesh::shapeElementsIterator(int shapeID,
1668 size_t nbElemsToReturn,
1669 const SMDS_MeshElement* sm1stElem) const
1671 return myCellFactory->GetShapeIterator< SMDS_ElemIterator >( shapeID, nbElemsToReturn, sm1stElem );
1674 ///////////////////////////////////////////////////////////////////////////////
1675 /// Do intersection of sets (more than 2)
1676 ///////////////////////////////////////////////////////////////////////////////
1677 static std::set<const SMDS_MeshElement*> *
1678 intersectionOfSets( std::set<const SMDS_MeshElement*> vs[], int numberOfSets )
1680 std::set<const SMDS_MeshElement*>* rsetA = new std::set<const SMDS_MeshElement*>(vs[0]);
1681 std::set<const SMDS_MeshElement*>* rsetB;
1683 for(int i=0; i<numberOfSets-1; i++)
1685 rsetB = new std::set<const SMDS_MeshElement*>();
1686 set_intersection(rsetA->begin(), rsetA->end(),
1687 vs[i+1].begin(), vs[i+1].end(),
1688 inserter(*rsetB, rsetB->begin()));
1694 ///////////////////////////////////////////////////////////////////////////////
1695 /// Return the list of finite elements owning the given element: elements
1696 /// containing all the nodes of the given element, for instance faces and
1697 /// volumes containing a given edge.
1698 ///////////////////////////////////////////////////////////////////////////////
1699 static std::set<const SMDS_MeshElement*> * getFinitElements(const SMDS_MeshElement * element)
1701 int numberOfSets=element->NbNodes();
1702 std::set<const SMDS_MeshElement*> *initSet = new std::set<const SMDS_MeshElement*>[numberOfSets];
1704 SMDS_NodeIteratorPtr itNodes = element->nodeIterator();
1707 while ( itNodes->more() )
1709 const SMDS_MeshNode * n = itNodes->next();
1710 for ( SMDS_ElemIteratorPtr itFe = n->GetInverseElementIterator(); itFe->more(); )
1711 initSet[i].insert( itFe->next() );
1714 std::set<const SMDS_MeshElement*> *retSet = intersectionOfSets( initSet, numberOfSets );
1719 ///////////////////////////////////////////////////////////////////////////////
1720 /// Return the std::list of nodes used only by the given elements
1721 ///////////////////////////////////////////////////////////////////////////////
1723 std::set<const SMDS_MeshElement*> *getExclusiveNodes(std::set<const SMDS_MeshElement*>& elements)
1725 std::set<const SMDS_MeshElement*> * toReturn = new std::set<const SMDS_MeshElement*>();
1726 std::set<const SMDS_MeshElement*>::iterator itElements = elements.begin();
1728 while( itElements != elements.end() )
1730 SMDS_NodeIteratorPtr itNodes = (*itElements)->nodeIterator();
1733 while( itNodes->more() )
1735 const SMDS_MeshNode * n = itNodes->next();
1736 SMDS_ElemIteratorPtr itFe = n->GetInverseElementIterator();
1737 std::set<const SMDS_MeshElement*> s;
1738 while ( itFe->more() )
1739 s.insert( itFe->next() );
1740 if ( s == elements ) toReturn->insert(n);
1746 ///////////////////////////////////////////////////////////////////////////////
1747 ///Find the children of an element that are made of given nodes
1748 ///@param setOfChildren The set in which matching children will be inserted
1749 ///@param element The element were to search matching children
1750 ///@param nodes The nodes that the children must have to be selected
1751 ///////////////////////////////////////////////////////////////////////////////
1752 void SMDS_Mesh::addChildrenWithNodes(std::set<const SMDS_MeshElement*>& setOfChildren,
1753 const SMDS_MeshElement * element,
1754 std::set<const SMDS_MeshElement*>& nodes)
1756 switch(element->GetType())
1759 throw SALOME_Exception("Internal Error: This should not happen");
1761 case SMDSAbs_0DElement:
1768 SMDS_ElemIteratorPtr itn=element->nodesIterator();
1771 const SMDS_MeshElement * e=itn->next();
1772 if(nodes.find(e)!=nodes.end())
1774 setOfChildren.insert(element);
1781 SMDS_ElemIteratorPtr itn=element->nodesIterator();
1784 const SMDS_MeshElement * e=itn->next();
1785 if(nodes.find(e)!=nodes.end())
1787 setOfChildren.insert(element);
1792 case SMDSAbs_Volume:
1793 case SMDSAbs_NbElementTypes:
1794 case SMDSAbs_All: break;
1798 ///////////////////////////////////////////////////////////////////////////////
1799 ///@param elem The element to delete
1800 ///@param removenodes if true remaining nodes will be removed
1801 ///////////////////////////////////////////////////////////////////////////////
1802 void SMDS_Mesh::RemoveElement(const SMDS_MeshElement * elem,
1803 const bool removenodes)
1805 std::vector<const SMDS_MeshElement *> removedElems;
1806 std::vector<const SMDS_MeshElement *> removedNodes;
1807 RemoveElement( elem, removedElems, removedNodes, removenodes );
1810 ///////////////////////////////////////////////////////////////////////////////
1811 ///@param elem The element to delete
1812 ///@param removedElems to be filled with all removed elements
1813 ///@param removedNodes to be filled with all removed nodes
1814 ///@param removenodes if true remaining nodes will be removed
1815 ///////////////////////////////////////////////////////////////////////////////
1816 void SMDS_Mesh::RemoveElement(const SMDS_MeshElement * elem,
1817 std::vector<const SMDS_MeshElement *>& removedElems,
1818 std::vector<const SMDS_MeshElement *>& removedNodes,
1821 // get finite elements built on elem
1822 std::set<const SMDS_MeshElement*> * s1;
1823 if ( (elem->GetType() == SMDSAbs_0DElement)
1824 || (elem->GetType() == SMDSAbs_Ball)
1825 || (elem->GetType() == SMDSAbs_Edge)
1826 || (elem->GetType() == SMDSAbs_Face)
1827 || (elem->GetType() == SMDSAbs_Volume) )
1829 s1 = new std::set<const SMDS_MeshElement*> ();
1833 s1 = getFinitElements(elem);
1835 // get exclusive nodes (which would become free afterwards)
1836 std::set<const SMDS_MeshElement*> * s2;
1837 if (elem->GetType() == SMDSAbs_Node) // a node is removed
1839 // do not remove nodes except elem
1840 s2 = new std::set<const SMDS_MeshElement*> ();
1845 s2 = getExclusiveNodes(*s1);
1847 // form the set of finite and construction elements to remove
1848 std::set<const SMDS_MeshElement*> s3;
1849 std::set<const SMDS_MeshElement*>::iterator it = s1->begin();
1850 while (it != s1->end())
1852 addChildrenWithNodes(s3, *it, *s2);
1856 if (elem->GetType() != SMDSAbs_Node)
1859 // remove finite and construction elements
1860 for( it = s3.begin();it != s3.end(); ++it )
1862 // Remove element from <InverseElements> of its nodes
1863 SMDS_NodeIteratorPtr itn = (*it)->nodeIterator();
1866 SMDS_MeshNode * n = const_cast<SMDS_MeshNode *> (itn->next());
1867 n->RemoveInverseElement((*it));
1870 int vtkid = (*it)->GetVtkID();
1872 switch ((*it)->GetType()) {
1874 throw SALOME_Exception(LOCALIZED("Internal Error: This should not happen"));
1876 case SMDSAbs_Edge: myInfo.RemoveEdge(*it); break;
1877 case SMDSAbs_Face: myInfo.RemoveFace(*it); break;
1878 case SMDSAbs_Volume: myInfo.RemoveVolume(*it); break;
1879 case SMDSAbs_Ball: myInfo.myNbBalls--; break;
1880 case SMDSAbs_0DElement: myInfo.myNb0DElements--; break;
1881 case SMDSAbs_All: // avoid compilation warning
1882 case SMDSAbs_NbElementTypes: break;
1884 removedElems.push_back( *it);
1886 myCellFactory->Free( static_cast< const SMDS_MeshCell*>( *it ));
1890 this->myGrid->GetCellTypesArray()->SetValue(vtkid, VTK_EMPTY_CELL);
1894 // remove exclusive (free) nodes
1897 for ( it = s2->begin(); it != s2->end(); ++it )
1900 myNodeFactory->Free( (*it) );
1901 removedNodes.push_back((*it));
1910 ///////////////////////////////////////////////////////////////////////////////
1911 ///@param elem The element to delete
1912 ///////////////////////////////////////////////////////////////////////////////
1913 void SMDS_Mesh::RemoveFreeElement(const SMDS_MeshElement * elem)
1915 const int vtkId = elem->GetVtkID();
1916 SMDSAbs_ElementType aType = elem->GetType();
1917 if ( aType == SMDSAbs_Node )
1919 // only free node can be removed by this method
1920 const SMDS_MeshNode* n = static_cast<const SMDS_MeshNode*>( elem );
1921 if ( n->NbInverseElements() == 0 ) { // free node
1923 myNodeFactory->Free( n );
1927 throw SALOME_Exception( LOCALIZED( "RemoveFreeElement: not a free node" ));
1932 // Remove element from <InverseElements> of its nodes
1933 SMDS_NodeIteratorPtr itn = elem->nodeIterator();
1934 while (itn->more()) {
1935 SMDS_MeshNode * n = const_cast<SMDS_MeshNode *>(itn->next());
1936 n->RemoveInverseElement(elem);
1939 // in meshes without descendants elements are always free
1941 case SMDSAbs_0DElement: myInfo.remove(elem); break;
1942 case SMDSAbs_Edge: myInfo.RemoveEdge(elem); break;
1943 case SMDSAbs_Face: myInfo.RemoveFace(elem); break;
1944 case SMDSAbs_Volume: myInfo.RemoveVolume(elem); break;
1945 case SMDSAbs_Ball: myInfo.remove(elem); break;
1948 myCellFactory->Free( elem );
1950 this->myGrid->GetCellTypesArray()->SetValue(vtkId, VTK_EMPTY_CELL);
1954 //=======================================================================
1956 * Checks if the element is present in mesh.
1958 //=======================================================================
1960 bool SMDS_Mesh::Contains (const SMDS_MeshElement* elem) const
1962 if ( !elem || elem->IsNull() )
1965 if ( elem->GetType() == SMDSAbs_Node )
1966 return ( elem == myNodeFactory->FindElement( elem->GetID() ));
1968 return ( elem == myCellFactory->FindElement( elem->GetID() ));
1971 //=======================================================================
1972 //function : MaxNodeID
1974 //=======================================================================
1976 int SMDS_Mesh::MaxNodeID() const
1978 return myNodeFactory->GetMaxID();
1981 //=======================================================================
1982 //function : MinNodeID
1984 //=======================================================================
1986 int SMDS_Mesh::MinNodeID() const
1988 return myNodeFactory->GetMinID();
1991 //=======================================================================
1992 //function : MaxElementID
1994 //=======================================================================
1996 int SMDS_Mesh::MaxElementID() const
1998 return myCellFactory->GetMaxID();
2001 //=======================================================================
2002 //function : MinElementID
2004 //=======================================================================
2006 int SMDS_Mesh::MinElementID() const
2008 return myCellFactory->GetMinID();
2011 //=======================================================================
2012 //function : Renumber
2013 //purpose : Renumber all nodes or elements.
2014 //=======================================================================
2016 // void SMDS_Mesh::Renumber (const bool isNodes, const int startID, const int deltaID)
2018 // if ( deltaID == 0 )
2023 //=======================================================================
2024 //function : GetElementType
2025 //purpose : Return type of element or node with id
2026 //=======================================================================
2028 SMDSAbs_ElementType SMDS_Mesh::GetElementType( const int id, const bool iselem ) const
2030 const SMDS_MeshElement* elem = 0;
2032 elem = myCellFactory->FindElement( id );
2034 elem = myNodeFactory->FindElement( id );
2036 return elem ? elem->GetType() : SMDSAbs_All;
2041 //********************************************************************
2042 //********************************************************************
2043 //******** *********
2044 //***** Methods for addition of quadratic elements ******
2045 //******** *********
2046 //********************************************************************
2047 //********************************************************************
2049 //=======================================================================
2050 //function : AddEdgeWithID
2052 //=======================================================================
2053 SMDS_MeshEdge* SMDS_Mesh::AddEdgeWithID(int n1, int n2, int n12, int ID)
2055 return SMDS_Mesh::AddEdgeWithID (myNodeFactory->FindNode(n1),
2056 myNodeFactory->FindNode(n2),
2057 myNodeFactory->FindNode(n12),
2061 //=======================================================================
2062 //function : AddEdge
2064 //=======================================================================
2065 SMDS_MeshEdge* SMDS_Mesh::AddEdge(const SMDS_MeshNode* n1,
2066 const SMDS_MeshNode* n2,
2067 const SMDS_MeshNode* n12)
2069 return SMDS_Mesh::AddEdgeWithID(n1, n2, n12, myCellFactory->GetFreeID());
2072 //=======================================================================
2073 //function : AddEdgeWithID
2075 //=======================================================================
2076 SMDS_MeshEdge* SMDS_Mesh::AddEdgeWithID(const SMDS_MeshNode * n1,
2077 const SMDS_MeshNode * n2,
2078 const SMDS_MeshNode * n12,
2081 if ( !n1 || !n2 || !n12 ) return 0;
2083 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2085 cell->init( SMDSEntity_Quad_Edge, /*nbNodes=*/3, n1, n2, n12 );
2086 myInfo.myNbQuadEdges++;
2087 return static_cast<SMDS_MeshEdge*>( cell );
2093 //=======================================================================
2094 //function : AddFace
2096 //=======================================================================
2097 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshNode * n1,
2098 const SMDS_MeshNode * n2,
2099 const SMDS_MeshNode * n3,
2100 const SMDS_MeshNode * n12,
2101 const SMDS_MeshNode * n23,
2102 const SMDS_MeshNode * n31)
2104 return SMDS_Mesh::AddFaceWithID(n1,n2,n3,n12,n23,n31,
2105 myCellFactory->GetFreeID());
2108 //=======================================================================
2109 //function : AddFaceWithID
2111 //=======================================================================
2112 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(int n1, int n2, int n3,
2113 int n12,int n23,int n31, int ID)
2115 return SMDS_Mesh::AddFaceWithID (myNodeFactory->FindNode(n1) ,
2116 myNodeFactory->FindNode(n2) ,
2117 myNodeFactory->FindNode(n3) ,
2118 myNodeFactory->FindNode(n12),
2119 myNodeFactory->FindNode(n23),
2120 myNodeFactory->FindNode(n31),
2124 //=======================================================================
2125 //function : AddFaceWithID
2127 //=======================================================================
2128 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshNode * n1,
2129 const SMDS_MeshNode * n2,
2130 const SMDS_MeshNode * n3,
2131 const SMDS_MeshNode * n12,
2132 const SMDS_MeshNode * n23,
2133 const SMDS_MeshNode * n31,
2136 if ( !n1 || !n2 || !n3 || !n12 || !n23 || !n31 ) return 0;
2137 if ( NbFaces() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2139 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2141 cell->init( SMDSEntity_Quad_Triangle, /*nbNodes=*/6, n1, n2, n3, n12, n23, n31 );
2142 myInfo.myNbQuadTriangles++;
2143 return static_cast<SMDS_MeshFace*>( cell );
2149 //=======================================================================
2150 //function : AddFace
2152 //=======================================================================
2153 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshNode * n1,
2154 const SMDS_MeshNode * n2,
2155 const SMDS_MeshNode * n3,
2156 const SMDS_MeshNode * n12,
2157 const SMDS_MeshNode * n23,
2158 const SMDS_MeshNode * n31,
2159 const SMDS_MeshNode * nCenter)
2161 return SMDS_Mesh::AddFaceWithID(n1,n2,n3,n12,n23,n31,nCenter,
2162 myCellFactory->GetFreeID());
2165 //=======================================================================
2166 //function : AddFaceWithID
2168 //=======================================================================
2169 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(int n1, int n2, int n3,
2170 int n12,int n23,int n31, int nCenter, int ID)
2172 return SMDS_Mesh::AddFaceWithID (myNodeFactory->FindNode(n1) ,
2173 myNodeFactory->FindNode(n2) ,
2174 myNodeFactory->FindNode(n3) ,
2175 myNodeFactory->FindNode(n12),
2176 myNodeFactory->FindNode(n23),
2177 myNodeFactory->FindNode(n31),
2178 myNodeFactory->FindNode(nCenter),
2182 //=======================================================================
2183 //function : AddFaceWithID
2185 //=======================================================================
2186 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshNode * n1,
2187 const SMDS_MeshNode * n2,
2188 const SMDS_MeshNode * n3,
2189 const SMDS_MeshNode * n12,
2190 const SMDS_MeshNode * n23,
2191 const SMDS_MeshNode * n31,
2192 const SMDS_MeshNode * nCenter,
2195 if ( !n1 || !n2 || !n3 || !n12 || !n23 || !n31 || !nCenter) return 0;
2196 if ( NbFaces() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2198 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2200 cell->init( SMDSEntity_BiQuad_Triangle, /*nbNodes=*/7, n1, n2, n3, n12, n23, n31, nCenter );
2201 myInfo.myNbBiQuadTriangles++;
2202 return static_cast<SMDS_MeshFace*>( cell );
2208 //=======================================================================
2209 //function : AddFace
2211 //=======================================================================
2212 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshNode * n1,
2213 const SMDS_MeshNode * n2,
2214 const SMDS_MeshNode * n3,
2215 const SMDS_MeshNode * n4,
2216 const SMDS_MeshNode * n12,
2217 const SMDS_MeshNode * n23,
2218 const SMDS_MeshNode * n34,
2219 const SMDS_MeshNode * n41)
2221 return SMDS_Mesh::AddFaceWithID(n1,n2,n3,n4,n12,n23,n34,n41,
2222 myCellFactory->GetFreeID());
2225 //=======================================================================
2226 //function : AddFaceWithID
2228 //=======================================================================
2229 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(int n1, int n2, int n3, int n4,
2230 int n12,int n23,int n34,int n41, int ID)
2232 return SMDS_Mesh::AddFaceWithID (myNodeFactory->FindNode(n1) ,
2233 myNodeFactory->FindNode(n2) ,
2234 myNodeFactory->FindNode(n3) ,
2235 myNodeFactory->FindNode(n4) ,
2236 myNodeFactory->FindNode(n12),
2237 myNodeFactory->FindNode(n23),
2238 myNodeFactory->FindNode(n34),
2239 myNodeFactory->FindNode(n41),
2243 //=======================================================================
2244 //function : AddFaceWithID
2246 //=======================================================================
2247 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshNode * n1,
2248 const SMDS_MeshNode * n2,
2249 const SMDS_MeshNode * n3,
2250 const SMDS_MeshNode * n4,
2251 const SMDS_MeshNode * n12,
2252 const SMDS_MeshNode * n23,
2253 const SMDS_MeshNode * n34,
2254 const SMDS_MeshNode * n41,
2257 if ( !n1 || !n2 || !n3 || !n4 || !n12 || !n23 || !n34 || !n41) return 0;
2258 if ( NbFaces() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2260 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2262 cell->init( SMDSEntity_Quad_Quadrangle, /*nbNodes=*/8, n1, n2, n3, n4, n12, n23, n34, n41 );
2263 myInfo.myNbQuadQuadrangles++;
2264 return static_cast<SMDS_MeshFace*>( cell );
2269 //=======================================================================
2270 //function : AddFace
2272 //=======================================================================
2273 SMDS_MeshFace* SMDS_Mesh::AddFace(const SMDS_MeshNode * n1,
2274 const SMDS_MeshNode * n2,
2275 const SMDS_MeshNode * n3,
2276 const SMDS_MeshNode * n4,
2277 const SMDS_MeshNode * n12,
2278 const SMDS_MeshNode * n23,
2279 const SMDS_MeshNode * n34,
2280 const SMDS_MeshNode * n41,
2281 const SMDS_MeshNode * nCenter)
2283 return SMDS_Mesh::AddFaceWithID(n1,n2,n3,n4,n12,n23,n34,n41,nCenter,
2284 myCellFactory->GetFreeID());
2287 //=======================================================================
2288 //function : AddFaceWithID
2290 //=======================================================================
2291 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(int n1, int n2, int n3, int n4,
2292 int n12,int n23,int n34,int n41, int nCenter, int ID)
2294 return SMDS_Mesh::AddFaceWithID (myNodeFactory->FindNode(n1) ,
2295 myNodeFactory->FindNode(n2) ,
2296 myNodeFactory->FindNode(n3) ,
2297 myNodeFactory->FindNode(n4) ,
2298 myNodeFactory->FindNode(n12),
2299 myNodeFactory->FindNode(n23),
2300 myNodeFactory->FindNode(n34),
2301 myNodeFactory->FindNode(n41),
2302 myNodeFactory->FindNode(nCenter),
2306 //=======================================================================
2307 //function : AddFaceWithID
2309 //=======================================================================
2310 SMDS_MeshFace* SMDS_Mesh::AddFaceWithID(const SMDS_MeshNode * n1,
2311 const SMDS_MeshNode * n2,
2312 const SMDS_MeshNode * n3,
2313 const SMDS_MeshNode * n4,
2314 const SMDS_MeshNode * n12,
2315 const SMDS_MeshNode * n23,
2316 const SMDS_MeshNode * n34,
2317 const SMDS_MeshNode * n41,
2318 const SMDS_MeshNode * nCenter,
2321 if ( !n1 || !n2 || !n3 || !n4 || !n12 || !n23 || !n34 || !n41 || !nCenter) return 0;
2322 if ( NbFaces() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2324 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2326 cell->init( SMDSEntity_BiQuad_Quadrangle,
2327 /*nbNodes=*/9, n1, n2, n3, n4, n12, n23, n34, n41, nCenter );
2328 myInfo.myNbBiQuadQuadrangles++;
2329 return static_cast<SMDS_MeshFace*>( cell );
2335 //=======================================================================
2336 //function : AddVolume
2338 //=======================================================================
2339 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
2340 const SMDS_MeshNode * n2,
2341 const SMDS_MeshNode * n3,
2342 const SMDS_MeshNode * n4,
2343 const SMDS_MeshNode * n12,
2344 const SMDS_MeshNode * n23,
2345 const SMDS_MeshNode * n31,
2346 const SMDS_MeshNode * n14,
2347 const SMDS_MeshNode * n24,
2348 const SMDS_MeshNode * n34)
2350 return SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n12, n23,
2351 n31, n14, n24, n34, myCellFactory->GetFreeID());
2354 //=======================================================================
2355 //function : AddVolumeWithID
2357 //=======================================================================
2358 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3, int n4,
2359 int n12,int n23,int n31,
2360 int n14,int n24,int n34, int ID)
2362 return SMDS_Mesh::AddVolumeWithID (myNodeFactory->FindNode(n1) ,
2363 myNodeFactory->FindNode(n2) ,
2364 myNodeFactory->FindNode(n3) ,
2365 myNodeFactory->FindNode(n4) ,
2366 myNodeFactory->FindNode(n12),
2367 myNodeFactory->FindNode(n23),
2368 myNodeFactory->FindNode(n31),
2369 myNodeFactory->FindNode(n14),
2370 myNodeFactory->FindNode(n24),
2371 myNodeFactory->FindNode(n34),
2375 //=======================================================================
2376 //function : AddVolumeWithID
2377 //purpose : 2d order tetrahedron of 10 nodes
2378 //=======================================================================
2379 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
2380 const SMDS_MeshNode * n2,
2381 const SMDS_MeshNode * n3,
2382 const SMDS_MeshNode * n4,
2383 const SMDS_MeshNode * n12,
2384 const SMDS_MeshNode * n23,
2385 const SMDS_MeshNode * n31,
2386 const SMDS_MeshNode * n14,
2387 const SMDS_MeshNode * n24,
2388 const SMDS_MeshNode * n34,
2391 if ( !n1 || !n2 || !n3 || !n4 || !n12 || !n23 || !n31 || !n14 || !n24 || !n34)
2393 if ( NbVolumes() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2395 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2397 cell->init( SMDSEntity_Quad_Tetra,
2398 /*nbNodes=*/10, n1, n2, n3, n4, n12, n23, n31, n14, n24, n34 );
2399 myInfo.myNbQuadTetras++;
2400 return static_cast<SMDS_MeshVolume*>( cell );
2406 //=======================================================================
2407 //function : AddVolume
2409 //=======================================================================
2410 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
2411 const SMDS_MeshNode * n2,
2412 const SMDS_MeshNode * n3,
2413 const SMDS_MeshNode * n4,
2414 const SMDS_MeshNode * n5,
2415 const SMDS_MeshNode * n12,
2416 const SMDS_MeshNode * n23,
2417 const SMDS_MeshNode * n34,
2418 const SMDS_MeshNode * n41,
2419 const SMDS_MeshNode * n15,
2420 const SMDS_MeshNode * n25,
2421 const SMDS_MeshNode * n35,
2422 const SMDS_MeshNode * n45)
2424 return SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n12, n23, n34, n41,
2425 n15, n25, n35, n45, myCellFactory->GetFreeID());
2428 //=======================================================================
2429 //function : AddVolumeWithID
2431 //=======================================================================
2432 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3, int n4, int n5,
2433 int n12,int n23,int n34,int n41,
2434 int n15,int n25,int n35,int n45, int ID)
2436 return SMDS_Mesh::AddVolumeWithID (myNodeFactory->FindNode(n1) ,
2437 myNodeFactory->FindNode(n2) ,
2438 myNodeFactory->FindNode(n3) ,
2439 myNodeFactory->FindNode(n4) ,
2440 myNodeFactory->FindNode(n5) ,
2441 myNodeFactory->FindNode(n12),
2442 myNodeFactory->FindNode(n23),
2443 myNodeFactory->FindNode(n34),
2444 myNodeFactory->FindNode(n41),
2445 myNodeFactory->FindNode(n15),
2446 myNodeFactory->FindNode(n25),
2447 myNodeFactory->FindNode(n35),
2448 myNodeFactory->FindNode(n45),
2452 //=======================================================================
2453 //function : AddVolumeWithID
2454 //purpose : 2d order pyramid of 13 nodes
2455 //=======================================================================
2456 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
2457 const SMDS_MeshNode * n2,
2458 const SMDS_MeshNode * n3,
2459 const SMDS_MeshNode * n4,
2460 const SMDS_MeshNode * n5,
2461 const SMDS_MeshNode * n12,
2462 const SMDS_MeshNode * n23,
2463 const SMDS_MeshNode * n34,
2464 const SMDS_MeshNode * n41,
2465 const SMDS_MeshNode * n15,
2466 const SMDS_MeshNode * n25,
2467 const SMDS_MeshNode * n35,
2468 const SMDS_MeshNode * n45,
2471 if (!n1 || !n2 || !n3 || !n4 || !n5 || !n12 || !n23 ||
2472 !n34 || !n41 || !n15 || !n25 || !n35 || !n45)
2474 if ( NbVolumes() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2476 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2478 cell->init( SMDSEntity_Quad_Pyramid,
2479 /*nbNodes=*/13, n1, n2, n3, n4, n5, n12, n23, n34, n41, n15, n25, n35, n45);
2480 myInfo.myNbQuadPyramids++;
2481 return static_cast<SMDS_MeshVolume*>( cell );
2487 //=======================================================================
2488 //function : AddVolume
2489 //purpose : 2d order Pentahedron (prism) with 15 nodes
2490 //=======================================================================
2491 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
2492 const SMDS_MeshNode * n2,
2493 const SMDS_MeshNode * n3,
2494 const SMDS_MeshNode * n4,
2495 const SMDS_MeshNode * n5,
2496 const SMDS_MeshNode * n6,
2497 const SMDS_MeshNode * n12,
2498 const SMDS_MeshNode * n23,
2499 const SMDS_MeshNode * n31,
2500 const SMDS_MeshNode * n45,
2501 const SMDS_MeshNode * n56,
2502 const SMDS_MeshNode * n64,
2503 const SMDS_MeshNode * n14,
2504 const SMDS_MeshNode * n25,
2505 const SMDS_MeshNode * n36)
2507 return SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n6, n12, n23, n31,
2508 n45, n56, n64, n14, n25, n36, myCellFactory->GetFreeID());
2511 //=======================================================================
2512 //function : AddVolumeWithID
2513 //purpose : 2d order Pentahedron (prism) with 15 nodes
2514 //=======================================================================
2515 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3,
2516 int n4, int n5, int n6,
2517 int n12,int n23,int n31,
2518 int n45,int n56,int n64,
2519 int n14,int n25,int n36, int ID)
2521 return SMDS_Mesh::AddVolumeWithID (myNodeFactory->FindNode(n1) ,
2522 myNodeFactory->FindNode(n2) ,
2523 myNodeFactory->FindNode(n3) ,
2524 myNodeFactory->FindNode(n4) ,
2525 myNodeFactory->FindNode(n5) ,
2526 myNodeFactory->FindNode(n6) ,
2527 myNodeFactory->FindNode(n12),
2528 myNodeFactory->FindNode(n23),
2529 myNodeFactory->FindNode(n31),
2530 myNodeFactory->FindNode(n45),
2531 myNodeFactory->FindNode(n56),
2532 myNodeFactory->FindNode(n64),
2533 myNodeFactory->FindNode(n14),
2534 myNodeFactory->FindNode(n25),
2535 myNodeFactory->FindNode(n36),
2539 //=======================================================================
2540 //function : AddVolumeWithID
2541 //purpose : 2d order Pentahedron (prism) with 15 nodes
2542 //=======================================================================
2543 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
2544 const SMDS_MeshNode * n2,
2545 const SMDS_MeshNode * n3,
2546 const SMDS_MeshNode * n4,
2547 const SMDS_MeshNode * n5,
2548 const SMDS_MeshNode * n6,
2549 const SMDS_MeshNode * n12,
2550 const SMDS_MeshNode * n23,
2551 const SMDS_MeshNode * n31,
2552 const SMDS_MeshNode * n45,
2553 const SMDS_MeshNode * n56,
2554 const SMDS_MeshNode * n64,
2555 const SMDS_MeshNode * n14,
2556 const SMDS_MeshNode * n25,
2557 const SMDS_MeshNode * n36,
2560 if (!n1 || !n2 || !n3 || !n4 || !n5 || !n6 || !n12 || !n23 ||
2561 !n31 || !n45 || !n56 || !n64 || !n14 || !n25 || !n36)
2563 if ( NbVolumes() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2565 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2567 cell->init( SMDSEntity_Quad_Penta, /*nbNodes=*/15,
2568 n1, n2, n3, n4, n5, n6, n12, n23, n31, n45, n56, n64, n14, n25, n36 );
2569 myInfo.myNbQuadPrisms++;
2570 return static_cast<SMDS_MeshVolume*>( cell );
2575 //=======================================================================
2576 //function : AddVolume
2577 //purpose : 2d order Pentahedron (prism) with 18 nodes
2578 //=======================================================================
2579 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
2580 const SMDS_MeshNode * n2,
2581 const SMDS_MeshNode * n3,
2582 const SMDS_MeshNode * n4,
2583 const SMDS_MeshNode * n5,
2584 const SMDS_MeshNode * n6,
2585 const SMDS_MeshNode * n12,
2586 const SMDS_MeshNode * n23,
2587 const SMDS_MeshNode * n31,
2588 const SMDS_MeshNode * n45,
2589 const SMDS_MeshNode * n56,
2590 const SMDS_MeshNode * n64,
2591 const SMDS_MeshNode * n14,
2592 const SMDS_MeshNode * n25,
2593 const SMDS_MeshNode * n36,
2594 const SMDS_MeshNode * n1245,
2595 const SMDS_MeshNode * n2356,
2596 const SMDS_MeshNode * n1346)
2598 return SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n6, n12, n23, n31,
2599 n45, n56, n64, n14, n25, n36, n1245, n2356, n1346,
2600 myCellFactory->GetFreeID());
2603 //=======================================================================
2604 //function : AddVolumeWithID
2605 //purpose : 2d order Pentahedron (prism) with 18 nodes
2606 //=======================================================================
2607 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3,
2608 int n4, int n5, int n6,
2609 int n12,int n23,int n31,
2610 int n45,int n56,int n64,
2611 int n14,int n25,int n36,
2612 int n1245, int n2356, int n1346, int ID)
2614 return SMDS_Mesh::AddVolumeWithID (myNodeFactory->FindNode(n1) ,
2615 myNodeFactory->FindNode(n2) ,
2616 myNodeFactory->FindNode(n3) ,
2617 myNodeFactory->FindNode(n4) ,
2618 myNodeFactory->FindNode(n5) ,
2619 myNodeFactory->FindNode(n6) ,
2620 myNodeFactory->FindNode(n12),
2621 myNodeFactory->FindNode(n23),
2622 myNodeFactory->FindNode(n31),
2623 myNodeFactory->FindNode(n45),
2624 myNodeFactory->FindNode(n56),
2625 myNodeFactory->FindNode(n64),
2626 myNodeFactory->FindNode(n14),
2627 myNodeFactory->FindNode(n25),
2628 myNodeFactory->FindNode(n36),
2629 myNodeFactory->FindNode(n1245),
2630 myNodeFactory->FindNode(n2356),
2631 myNodeFactory->FindNode(n1346),
2635 //=======================================================================
2636 //function : AddVolumeWithID
2637 //purpose : 2d order Pentahedron (prism) with 18 nodes
2638 //=======================================================================
2639 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
2640 const SMDS_MeshNode * n2,
2641 const SMDS_MeshNode * n3,
2642 const SMDS_MeshNode * n4,
2643 const SMDS_MeshNode * n5,
2644 const SMDS_MeshNode * n6,
2645 const SMDS_MeshNode * n12,
2646 const SMDS_MeshNode * n23,
2647 const SMDS_MeshNode * n31,
2648 const SMDS_MeshNode * n45,
2649 const SMDS_MeshNode * n56,
2650 const SMDS_MeshNode * n64,
2651 const SMDS_MeshNode * n14,
2652 const SMDS_MeshNode * n25,
2653 const SMDS_MeshNode * n36,
2654 const SMDS_MeshNode * n1245,
2655 const SMDS_MeshNode * n2356,
2656 const SMDS_MeshNode * n1346,
2659 //MESSAGE("AddVolumeWithID penta18 "<< ID);
2660 if (!n1 || !n2 || !n3 || !n4 || !n5 || !n6 || !n12 || !n23 ||
2661 !n31 || !n45 || !n56 || !n64 || !n14 || !n25 || !n36 || !n1245 || !n2356 || !n1346)
2663 if ( NbVolumes() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2665 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2667 cell->init( SMDSEntity_BiQuad_Penta, /*nbNodes=*/18, n1, n2, n3, n4, n5, n6,
2668 n12, n23, n31, n45, n56, n64, n14, n25, n36, n1245, n2356, n1346 );
2669 myInfo.myNbBiQuadPrisms++;
2670 return static_cast<SMDS_MeshVolume*>( cell );
2676 //=======================================================================
2677 //function : AddVolume
2679 //=======================================================================
2680 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
2681 const SMDS_MeshNode * n2,
2682 const SMDS_MeshNode * n3,
2683 const SMDS_MeshNode * n4,
2684 const SMDS_MeshNode * n5,
2685 const SMDS_MeshNode * n6,
2686 const SMDS_MeshNode * n7,
2687 const SMDS_MeshNode * n8,
2688 const SMDS_MeshNode * n12,
2689 const SMDS_MeshNode * n23,
2690 const SMDS_MeshNode * n34,
2691 const SMDS_MeshNode * n41,
2692 const SMDS_MeshNode * n56,
2693 const SMDS_MeshNode * n67,
2694 const SMDS_MeshNode * n78,
2695 const SMDS_MeshNode * n85,
2696 const SMDS_MeshNode * n15,
2697 const SMDS_MeshNode * n26,
2698 const SMDS_MeshNode * n37,
2699 const SMDS_MeshNode * n48)
2701 return SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n6, n7, n8, n12, n23, n34, n41,
2702 n56, n67, n78, n85, n15, n26, n37, n48,
2703 myCellFactory->GetFreeID());
2706 //=======================================================================
2707 //function : AddVolumeWithID
2709 //=======================================================================
2710 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3, int n4,
2711 int n5, int n6, int n7, int n8,
2712 int n12,int n23,int n34,int n41,
2713 int n56,int n67,int n78,int n85,
2714 int n15,int n26,int n37,int n48, int ID)
2716 return SMDS_Mesh::AddVolumeWithID (myNodeFactory->FindNode(n1),
2717 myNodeFactory->FindNode(n2),
2718 myNodeFactory->FindNode(n3),
2719 myNodeFactory->FindNode(n4),
2720 myNodeFactory->FindNode(n5),
2721 myNodeFactory->FindNode(n6),
2722 myNodeFactory->FindNode(n7),
2723 myNodeFactory->FindNode(n8),
2724 myNodeFactory->FindNode(n12),
2725 myNodeFactory->FindNode(n23),
2726 myNodeFactory->FindNode(n34),
2727 myNodeFactory->FindNode(n41),
2728 myNodeFactory->FindNode(n56),
2729 myNodeFactory->FindNode(n67),
2730 myNodeFactory->FindNode(n78),
2731 myNodeFactory->FindNode(n85),
2732 myNodeFactory->FindNode(n15),
2733 myNodeFactory->FindNode(n26),
2734 myNodeFactory->FindNode(n37),
2735 myNodeFactory->FindNode(n48),
2739 //=======================================================================
2740 //function : AddVolumeWithID
2741 //purpose : 2d order Hexahedrons with 20 nodes
2742 //=======================================================================
2743 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
2744 const SMDS_MeshNode * n2,
2745 const SMDS_MeshNode * n3,
2746 const SMDS_MeshNode * n4,
2747 const SMDS_MeshNode * n5,
2748 const SMDS_MeshNode * n6,
2749 const SMDS_MeshNode * n7,
2750 const SMDS_MeshNode * n8,
2751 const SMDS_MeshNode * n12,
2752 const SMDS_MeshNode * n23,
2753 const SMDS_MeshNode * n34,
2754 const SMDS_MeshNode * n41,
2755 const SMDS_MeshNode * n56,
2756 const SMDS_MeshNode * n67,
2757 const SMDS_MeshNode * n78,
2758 const SMDS_MeshNode * n85,
2759 const SMDS_MeshNode * n15,
2760 const SMDS_MeshNode * n26,
2761 const SMDS_MeshNode * n37,
2762 const SMDS_MeshNode * n48,
2765 if (!n1 || !n2 || !n3 || !n4 || !n5 || !n6 || !n7 || !n8 || !n12 || !n23 ||
2766 !n34 || !n41 || !n56 || !n67 || !n78 || !n85 || !n15 || !n26 || !n37 || !n48)
2768 if ( NbVolumes() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2770 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2772 cell->init( SMDSEntity_Quad_Hexa, /*nbNodes=*/20, n1, n2, n3, n4, n5, n6, n7, n8,
2773 n12, n23, n34, n41, n56, n67, n78, n85, n15, n26, n37, n48 );
2774 myInfo.myNbQuadHexas++;
2775 return static_cast<SMDS_MeshVolume*>( cell );
2780 //=======================================================================
2781 //function : AddVolume
2783 //=======================================================================
2784 SMDS_MeshVolume* SMDS_Mesh::AddVolume(const SMDS_MeshNode * n1,
2785 const SMDS_MeshNode * n2,
2786 const SMDS_MeshNode * n3,
2787 const SMDS_MeshNode * n4,
2788 const SMDS_MeshNode * n5,
2789 const SMDS_MeshNode * n6,
2790 const SMDS_MeshNode * n7,
2791 const SMDS_MeshNode * n8,
2792 const SMDS_MeshNode * n12,
2793 const SMDS_MeshNode * n23,
2794 const SMDS_MeshNode * n34,
2795 const SMDS_MeshNode * n41,
2796 const SMDS_MeshNode * n56,
2797 const SMDS_MeshNode * n67,
2798 const SMDS_MeshNode * n78,
2799 const SMDS_MeshNode * n85,
2800 const SMDS_MeshNode * n15,
2801 const SMDS_MeshNode * n26,
2802 const SMDS_MeshNode * n37,
2803 const SMDS_MeshNode * n48,
2804 const SMDS_MeshNode * n1234,
2805 const SMDS_MeshNode * n1256,
2806 const SMDS_MeshNode * n2367,
2807 const SMDS_MeshNode * n3478,
2808 const SMDS_MeshNode * n1458,
2809 const SMDS_MeshNode * n5678,
2810 const SMDS_MeshNode * nCenter)
2812 return SMDS_Mesh::AddVolumeWithID(n1, n2, n3, n4, n5, n6, n7, n8, n12, n23, n34, n41,
2813 n56, n67, n78, n85, n15, n26, n37, n48,
2814 n1234, n1256, n2367, n3478, n1458, n5678, nCenter,
2815 myCellFactory->GetFreeID());
2818 //=======================================================================
2819 //function : AddVolumeWithID
2821 //=======================================================================
2822 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(int n1, int n2, int n3, int n4,
2823 int n5, int n6, int n7, int n8,
2824 int n12,int n23,int n34,int n41,
2825 int n56,int n67,int n78,int n85,
2826 int n15,int n26,int n37,int n48,
2827 int n1234,int n1256,int n2367,int n3478,
2828 int n1458,int n5678,int nCenter, int ID)
2830 return SMDS_Mesh::AddVolumeWithID (myNodeFactory->FindNode(n1),
2831 myNodeFactory->FindNode(n2),
2832 myNodeFactory->FindNode(n3),
2833 myNodeFactory->FindNode(n4),
2834 myNodeFactory->FindNode(n5),
2835 myNodeFactory->FindNode(n6),
2836 myNodeFactory->FindNode(n7),
2837 myNodeFactory->FindNode(n8),
2838 myNodeFactory->FindNode(n12),
2839 myNodeFactory->FindNode(n23),
2840 myNodeFactory->FindNode(n34),
2841 myNodeFactory->FindNode(n41),
2842 myNodeFactory->FindNode(n56),
2843 myNodeFactory->FindNode(n67),
2844 myNodeFactory->FindNode(n78),
2845 myNodeFactory->FindNode(n85),
2846 myNodeFactory->FindNode(n15),
2847 myNodeFactory->FindNode(n26),
2848 myNodeFactory->FindNode(n37),
2849 myNodeFactory->FindNode(n48),
2850 myNodeFactory->FindNode(n1234),
2851 myNodeFactory->FindNode(n1256),
2852 myNodeFactory->FindNode(n2367),
2853 myNodeFactory->FindNode(n3478),
2854 myNodeFactory->FindNode(n1458),
2855 myNodeFactory->FindNode(n5678),
2856 myNodeFactory->FindNode(nCenter),
2860 //=======================================================================
2861 //function : AddVolumeWithID
2862 //purpose : 2d order Hexahedrons with 27 nodes
2863 //=======================================================================
2864 SMDS_MeshVolume* SMDS_Mesh::AddVolumeWithID(const SMDS_MeshNode * n1,
2865 const SMDS_MeshNode * n2,
2866 const SMDS_MeshNode * n3,
2867 const SMDS_MeshNode * n4,
2868 const SMDS_MeshNode * n5,
2869 const SMDS_MeshNode * n6,
2870 const SMDS_MeshNode * n7,
2871 const SMDS_MeshNode * n8,
2872 const SMDS_MeshNode * n12,
2873 const SMDS_MeshNode * n23,
2874 const SMDS_MeshNode * n34,
2875 const SMDS_MeshNode * n41,
2876 const SMDS_MeshNode * n56,
2877 const SMDS_MeshNode * n67,
2878 const SMDS_MeshNode * n78,
2879 const SMDS_MeshNode * n85,
2880 const SMDS_MeshNode * n15,
2881 const SMDS_MeshNode * n26,
2882 const SMDS_MeshNode * n37,
2883 const SMDS_MeshNode * n48,
2884 const SMDS_MeshNode * n1234,
2885 const SMDS_MeshNode * n1256,
2886 const SMDS_MeshNode * n2367,
2887 const SMDS_MeshNode * n3478,
2888 const SMDS_MeshNode * n1458,
2889 const SMDS_MeshNode * n5678,
2890 const SMDS_MeshNode * nCenter,
2893 if (!n1 || !n2 || !n3 || !n4 || !n5 || !n6 || !n7 || !n8 || !n12 || !n23 ||
2894 !n34 || !n41 || !n56 || !n67 || !n78 || !n85 || !n15 || !n26 || !n37 || !n48 ||
2895 !n1234 || !n1256 || !n2367 || !n3478 || !n1458 || !n5678 || !nCenter )
2897 if ( NbVolumes() % CHECKMEMORY_INTERVAL == 0 ) CheckMemory();
2899 if ( SMDS_MeshCell* cell = myCellFactory->NewCell( ID ))
2901 cell->init( SMDSEntity_TriQuad_Hexa, /*nbNodes=*/27, n1, n2, n3, n4, n5, n6, n7, n8,
2902 n12, n23, n34, n41, n56, n67, n78, n85, n15, n26, n37, n48,
2903 n1234, n1256, n2367, n3478, n1458, n5678, nCenter);
2904 myInfo.myNbTriQuadHexas++;
2905 return static_cast<SMDS_MeshVolume*>( cell );
2910 void SMDS_Mesh::dumpGrid(std::string ficdump)
2912 // vtkUnstructuredGridWriter* aWriter = vtkUnstructuredGridWriter::New();
2913 // aWriter->SetFileName(ficdump.c_str());
2914 // aWriter->SetInput(myGrid);
2915 // if(myGrid->GetNumberOfCells())
2917 // aWriter->Write();
2919 // aWriter->Delete();
2920 ficdump = ficdump + "_connectivity";
2921 std::ofstream ficcon(ficdump.c_str(), ios::out);
2922 int nbPoints = myGrid->GetNumberOfPoints();
2923 ficcon << "-------------------------------- points " << nbPoints << endl;
2924 for (int i=0; i<nbPoints; i++)
2926 ficcon << i << " " << *(myGrid->GetPoint(i)) << " " << *(myGrid->GetPoint(i)+1) << " " << " " << *(myGrid->GetPoint(i)+2) << endl;
2928 int nbCells = myGrid->GetNumberOfCells();
2929 ficcon << "-------------------------------- cells " << nbCells << endl;
2930 for (int i=0; i<nbCells; i++)
2932 ficcon << i << " - " << myGrid->GetCell(i)->GetCellType() << " -";
2933 int nbptcell = myGrid->GetCell(i)->GetNumberOfPoints();
2934 vtkIdList *listid = myGrid->GetCell(i)->GetPointIds();
2935 for (int j=0; j<nbptcell; j++)
2937 ficcon << " " << listid->GetId(j);
2941 ficcon << "-------------------------------- connectivity " << nbPoints << endl;
2942 vtkCellLinks *links = myGrid->GetLinks();
2943 for (int i=0; i<nbPoints; i++)
2945 int ncells = links->GetNcells(i);
2946 vtkIdType *cells = links->GetCells(i);
2947 ficcon << i << " - " << ncells << " -";
2948 for (int j=0; j<ncells; j++)
2950 ficcon << " " << cells[j];
2958 void SMDS_Mesh::CompactMesh()
2960 this->myCompactTime = this->myModifTime;
2962 bool idsChange = HasNumerationHoles();
2965 std::set< SMDS_ElementHolder* >::iterator holder = myElemHolders.begin();
2966 for ( ; holder != myElemHolders.end(); ++holder )
2967 (*holder)->beforeCompacting();
2969 int oldCellSize = myCellFactory->GetMaxID();
2971 // remove "holes" in SMDS numeration
2972 std::vector<int> idNodesOldToNew, idCellsNewToOld, idCellsOldToNew;
2973 myNodeFactory->Compact( idNodesOldToNew );
2974 myCellFactory->Compact( idCellsNewToOld );
2976 // make VTK IDs correspond to SMDS IDs
2977 int newNodeSize = myNodeFactory->NbUsedElements();
2978 int newCellSize = myCellFactory->NbUsedElements();
2979 myGrid->compactGrid( idNodesOldToNew, newNodeSize, idCellsNewToOld, newCellSize );
2981 if ( idsChange && !myElemHolders.empty() )
2983 // idCellsNewToOld -> idCellsOldToNew
2984 idCellsOldToNew.resize( oldCellSize, oldCellSize );
2985 for ( size_t iNew = 0; iNew < idCellsNewToOld.size(); ++iNew )
2987 if ( idCellsNewToOld[ iNew ] >= (int) idCellsOldToNew.size() )
2988 idCellsOldToNew.resize( ( 1 + idCellsNewToOld[ iNew ]) * 1.5, oldCellSize );
2989 idCellsOldToNew[ idCellsNewToOld[ iNew ]] = iNew;
2993 std::set< SMDS_ElementHolder* >::iterator holder = myElemHolders.begin();
2994 for ( ; holder != myElemHolders.end(); ++holder )
2996 (*holder)->restoreElements( idNodesOldToNew, idCellsOldToNew );
2998 (*holder)->compact();
3003 int SMDS_Mesh::FromVtkToSmds( int vtkid ) const
3005 return myCellFactory->FromVtkToSmds( vtkid );
3008 double SMDS_Mesh::getMaxDim()
3010 double dmax = 1.e-3;
3011 if ((xmax - xmin) > dmax) dmax = xmax -xmin;
3012 if ((ymax - ymin) > dmax) dmax = ymax -ymin;
3013 if ((zmax - zmin) > dmax) dmax = zmax -zmin;
3017 //! modification that needs compact structure and redraw
3018 void SMDS_Mesh::Modified()
3020 if (this->myModified)
3023 this->myModifTime++;
3028 //! get last modification timeStamp
3029 vtkMTimeType SMDS_Mesh::GetMTime() const
3031 return this->myModifTime;
3034 bool SMDS_Mesh::IsCompacted()
3036 return ( this->myCompactTime == this->myModifTime );
3039 //! are there holes in elements or nodes numeration
3040 bool SMDS_Mesh::HasNumerationHoles()
3042 return ( myNodeFactory->CompactChangePointers() ||
3043 myCellFactory->CompactChangePointers() );
3046 void SMDS_Mesh::setNbShapes( size_t nbShapes )
3048 myNodeFactory->SetNbShapes( nbShapes );