1 // Copyright (C) 2007-2008 CEA/DEN, EDF R&D, OPEN CASCADE
3 // Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
6 // This library is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 2.1 of the License.
11 // This library is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14 // Lesser General Public License for more details.
16 // You should have received a copy of the GNU Lesser General Public
17 // License along with this library; if not, write to the Free Software
18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
22 // SMESH SMDS : implementaion of Salome mesh data structure
23 // File : SMDS_Mesh.hxx
26 #ifndef _SMDS_Mesh_HeaderFile
27 #define _SMDS_Mesh_HeaderFile
29 #include "SMESH_SMDS.hxx"
31 #include "SMDS_MeshNode.hxx"
32 #include "SMDS_MeshCell.hxx"
33 #include "SMDS_Mesh0DElement.hxx"
34 #include "SMDS_MeshEdge.hxx"
35 #include "SMDS_MeshFace.hxx"
36 #include "SMDS_MeshVolume.hxx"
37 #include "SMDS_MeshNodeIDFactory.hxx"
38 #include "SMDS_MeshElementIDFactory.hxx"
39 #include "SMDS_MeshInfo.hxx"
40 #include "SMDS_ElemIterator.hxx"
41 #include "SMDS_VolumeOfNodes.hxx"
42 #include "SMDS_VtkEdge.hxx"
43 #include "SMDS_VtkFace.hxx"
44 #include "SMDS_VtkVolume.hxx"
45 #include "ObjectPool.hxx"
47 #include <boost/shared_ptr.hpp>
51 #include <vtkSystemIncludes.h>
53 #include "Utils_SALOME_Exception.hxx"
54 #define MYASSERT(val) if (!(val)) throw SALOME_Exception(LOCALIZED("assertion not verified"));
56 class vtkUnstructuredGrid;
58 class SMDS_EXPORT SMDS_Mesh:public SMDS_MeshObject{
60 friend class SMDS_MeshElementIDFactory;
61 friend class SMDS_MeshVolumeVtkNodes;
65 //! to retreive this SMDS_Mesh instance from its elements (index stored in SMDS_Elements)
66 static std::vector<SMDS_Mesh*> _meshList;
68 //! actual nodes coordinates, cells definition and reverse connectivity are stored in a vtkUnstructuredGrid
69 inline vtkUnstructuredGrid* getGrid() {return myGrid; };
70 inline int getMeshId() {return myMeshId; };
72 SMDS_NodeIteratorPtr nodesIterator() const;
73 SMDS_0DElementIteratorPtr elements0dIterator() const;
74 SMDS_EdgeIteratorPtr edgesIterator() const;
75 SMDS_FaceIteratorPtr facesIterator() const;
76 SMDS_VolumeIteratorPtr volumesIterator() const;
77 SMDS_ElemIteratorPtr elementsIterator(SMDSAbs_ElementType type=SMDSAbs_All) const;
79 SMDSAbs_ElementType GetElementType( const int id, const bool iselem ) const;
81 SMDS_Mesh *AddSubMesh();
83 virtual SMDS_MeshNode* AddNodeWithID(double x, double y, double z, int ID);
84 virtual SMDS_MeshNode* AddNode(double x, double y, double z);
86 virtual SMDS_Mesh0DElement* Add0DElementWithID(int n, int ID);
87 virtual SMDS_Mesh0DElement* Add0DElementWithID(const SMDS_MeshNode * n, int ID);
88 virtual SMDS_Mesh0DElement* Add0DElement (const SMDS_MeshNode * n);
90 virtual SMDS_MeshEdge* AddEdgeWithID(int n1, int n2, int ID);
91 virtual SMDS_MeshEdge* AddEdgeWithID(const SMDS_MeshNode * n1,
92 const SMDS_MeshNode * n2,
94 virtual SMDS_MeshEdge* AddEdge(const SMDS_MeshNode * n1,
95 const SMDS_MeshNode * n2);
97 // 2d order edge with 3 nodes: n12 - node between n1 and n2
98 virtual SMDS_MeshEdge* AddEdgeWithID(int n1, int n2, int n12, int ID);
99 virtual SMDS_MeshEdge* AddEdgeWithID(const SMDS_MeshNode * n1,
100 const SMDS_MeshNode * n2,
101 const SMDS_MeshNode * n12,
103 virtual SMDS_MeshEdge* AddEdge(const SMDS_MeshNode * n1,
104 const SMDS_MeshNode * n2,
105 const SMDS_MeshNode * n12);
107 virtual SMDS_MeshFace* AddFaceWithID(int n1, int n2, int n3, int ID);
108 virtual SMDS_MeshFace* AddFaceWithID(const SMDS_MeshNode * n1,
109 const SMDS_MeshNode * n2,
110 const SMDS_MeshNode * n3,
112 virtual SMDS_MeshFace* AddFace(const SMDS_MeshNode * n1,
113 const SMDS_MeshNode * n2,
114 const SMDS_MeshNode * n3);
116 virtual SMDS_MeshFace* AddFaceWithID(int n1, int n2, int n3, int n4, int ID);
117 virtual SMDS_MeshFace* AddFaceWithID(const SMDS_MeshNode * n1,
118 const SMDS_MeshNode * n2,
119 const SMDS_MeshNode * n3,
120 const SMDS_MeshNode * n4,
122 virtual SMDS_MeshFace* AddFace(const SMDS_MeshNode * n1,
123 const SMDS_MeshNode * n2,
124 const SMDS_MeshNode * n3,
125 const SMDS_MeshNode * n4);
127 virtual SMDS_MeshFace* AddFaceWithID(const SMDS_MeshEdge * e1,
128 const SMDS_MeshEdge * e2,
129 const SMDS_MeshEdge * e3, int ID);
130 virtual SMDS_MeshFace* AddFace(const SMDS_MeshEdge * e1,
131 const SMDS_MeshEdge * e2,
132 const SMDS_MeshEdge * e3);
134 virtual SMDS_MeshFace* AddFaceWithID(const SMDS_MeshEdge * e1,
135 const SMDS_MeshEdge * e2,
136 const SMDS_MeshEdge * e3,
137 const SMDS_MeshEdge * e4, int ID);
138 virtual SMDS_MeshFace* AddFace(const SMDS_MeshEdge * e1,
139 const SMDS_MeshEdge * e2,
140 const SMDS_MeshEdge * e3,
141 const SMDS_MeshEdge * e4);
143 // 2d order triangle of 6 nodes
144 virtual SMDS_MeshFace* AddFaceWithID(int n1, int n2, int n3,
145 int n12,int n23,int n31, int ID);
146 virtual SMDS_MeshFace* AddFaceWithID(const SMDS_MeshNode * n1,
147 const SMDS_MeshNode * n2,
148 const SMDS_MeshNode * n3,
149 const SMDS_MeshNode * n12,
150 const SMDS_MeshNode * n23,
151 const SMDS_MeshNode * n31,
153 virtual SMDS_MeshFace* AddFace(const SMDS_MeshNode * n1,
154 const SMDS_MeshNode * n2,
155 const SMDS_MeshNode * n3,
156 const SMDS_MeshNode * n12,
157 const SMDS_MeshNode * n23,
158 const SMDS_MeshNode * n31);
160 // 2d order quadrangle
161 virtual SMDS_MeshFace* AddFaceWithID(int n1, int n2, int n3, int n4,
162 int n12,int n23,int n34,int n41, int ID);
163 virtual SMDS_MeshFace* AddFaceWithID(const SMDS_MeshNode * n1,
164 const SMDS_MeshNode * n2,
165 const SMDS_MeshNode * n3,
166 const SMDS_MeshNode * n4,
167 const SMDS_MeshNode * n12,
168 const SMDS_MeshNode * n23,
169 const SMDS_MeshNode * n34,
170 const SMDS_MeshNode * n41,
172 virtual SMDS_MeshFace* AddFace(const SMDS_MeshNode * n1,
173 const SMDS_MeshNode * n2,
174 const SMDS_MeshNode * n3,
175 const SMDS_MeshNode * n4,
176 const SMDS_MeshNode * n12,
177 const SMDS_MeshNode * n23,
178 const SMDS_MeshNode * n34,
179 const SMDS_MeshNode * n41);
181 virtual SMDS_MeshVolume* AddVolumeWithID(int n1, int n2, int n3, int n4, int ID);
182 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshNode * n1,
183 const SMDS_MeshNode * n2,
184 const SMDS_MeshNode * n3,
185 const SMDS_MeshNode * n4,
187 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshNode * n1,
188 const SMDS_MeshNode * n2,
189 const SMDS_MeshNode * n3,
190 const SMDS_MeshNode * n4);
192 virtual SMDS_MeshVolume* AddVolumeWithID(int n1, int n2, int n3, int n4,
194 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshNode * n1,
195 const SMDS_MeshNode * n2,
196 const SMDS_MeshNode * n3,
197 const SMDS_MeshNode * n4,
198 const SMDS_MeshNode * n5,
200 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshNode * n1,
201 const SMDS_MeshNode * n2,
202 const SMDS_MeshNode * n3,
203 const SMDS_MeshNode * n4,
204 const SMDS_MeshNode * n5);
206 virtual SMDS_MeshVolume* AddVolumeWithID(int n1, int n2, int n3, int n4,
207 int n5, int n6, int ID);
208 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshNode * n1,
209 const SMDS_MeshNode * n2,
210 const SMDS_MeshNode * n3,
211 const SMDS_MeshNode * n4,
212 const SMDS_MeshNode * n5,
213 const SMDS_MeshNode * n6,
215 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshNode * n1,
216 const SMDS_MeshNode * n2,
217 const SMDS_MeshNode * n3,
218 const SMDS_MeshNode * n4,
219 const SMDS_MeshNode * n5,
220 const SMDS_MeshNode * n6);
222 virtual SMDS_MeshVolume* AddVolumeWithID(int n1, int n2, int n3, int n4,
223 int n5, int n6, int n7, int n8, int ID);
224 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshNode * n1,
225 const SMDS_MeshNode * n2,
226 const SMDS_MeshNode * n3,
227 const SMDS_MeshNode * n4,
228 const SMDS_MeshNode * n5,
229 const SMDS_MeshNode * n6,
230 const SMDS_MeshNode * n7,
231 const SMDS_MeshNode * n8,
233 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshNode * n1,
234 const SMDS_MeshNode * n2,
235 const SMDS_MeshNode * n3,
236 const SMDS_MeshNode * n4,
237 const SMDS_MeshNode * n5,
238 const SMDS_MeshNode * n6,
239 const SMDS_MeshNode * n7,
240 const SMDS_MeshNode * n8);
242 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshFace * f1,
243 const SMDS_MeshFace * f2,
244 const SMDS_MeshFace * f3,
245 const SMDS_MeshFace * f4, int ID);
246 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshFace * f1,
247 const SMDS_MeshFace * f2,
248 const SMDS_MeshFace * f3,
249 const SMDS_MeshFace * f4);
251 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshFace * f1,
252 const SMDS_MeshFace * f2,
253 const SMDS_MeshFace * f3,
254 const SMDS_MeshFace * f4,
255 const SMDS_MeshFace * f5, int ID);
256 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshFace * f1,
257 const SMDS_MeshFace * f2,
258 const SMDS_MeshFace * f3,
259 const SMDS_MeshFace * f4,
260 const SMDS_MeshFace * f5);
262 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshFace * f1,
263 const SMDS_MeshFace * f2,
264 const SMDS_MeshFace * f3,
265 const SMDS_MeshFace * f4,
266 const SMDS_MeshFace * f5,
267 const SMDS_MeshFace * f6, int ID);
268 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshFace * f1,
269 const SMDS_MeshFace * f2,
270 const SMDS_MeshFace * f3,
271 const SMDS_MeshFace * f4,
272 const SMDS_MeshFace * f5,
273 const SMDS_MeshFace * f6);
275 // 2d order tetrahedron of 10 nodes
276 virtual SMDS_MeshVolume* AddVolumeWithID(int n1, int n2, int n3, int n4,
277 int n12,int n23,int n31,
278 int n14,int n24,int n34, int ID);
279 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshNode * n1,
280 const SMDS_MeshNode * n2,
281 const SMDS_MeshNode * n3,
282 const SMDS_MeshNode * n4,
283 const SMDS_MeshNode * n12,
284 const SMDS_MeshNode * n23,
285 const SMDS_MeshNode * n31,
286 const SMDS_MeshNode * n14,
287 const SMDS_MeshNode * n24,
288 const SMDS_MeshNode * n34,
290 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshNode * n1,
291 const SMDS_MeshNode * n2,
292 const SMDS_MeshNode * n3,
293 const SMDS_MeshNode * n4,
294 const SMDS_MeshNode * n12,
295 const SMDS_MeshNode * n23,
296 const SMDS_MeshNode * n31,
297 const SMDS_MeshNode * n14,
298 const SMDS_MeshNode * n24,
299 const SMDS_MeshNode * n34);
301 // 2d order pyramid of 13 nodes
302 virtual SMDS_MeshVolume* AddVolumeWithID(int n1, int n2, int n3, int n4, int n5,
303 int n12,int n23,int n34,int n41,
304 int n15,int n25,int n35,int n45,
306 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshNode * n1,
307 const SMDS_MeshNode * n2,
308 const SMDS_MeshNode * n3,
309 const SMDS_MeshNode * n4,
310 const SMDS_MeshNode * n5,
311 const SMDS_MeshNode * n12,
312 const SMDS_MeshNode * n23,
313 const SMDS_MeshNode * n34,
314 const SMDS_MeshNode * n41,
315 const SMDS_MeshNode * n15,
316 const SMDS_MeshNode * n25,
317 const SMDS_MeshNode * n35,
318 const SMDS_MeshNode * n45,
320 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshNode * n1,
321 const SMDS_MeshNode * n2,
322 const SMDS_MeshNode * n3,
323 const SMDS_MeshNode * n4,
324 const SMDS_MeshNode * n5,
325 const SMDS_MeshNode * n12,
326 const SMDS_MeshNode * n23,
327 const SMDS_MeshNode * n34,
328 const SMDS_MeshNode * n41,
329 const SMDS_MeshNode * n15,
330 const SMDS_MeshNode * n25,
331 const SMDS_MeshNode * n35,
332 const SMDS_MeshNode * n45);
334 // 2d order Pentahedron with 15 nodes
335 virtual SMDS_MeshVolume* AddVolumeWithID(int n1, int n2, int n3,
336 int n4, int n5, int n6,
337 int n12,int n23,int n31,
338 int n45,int n56,int n64,
339 int n14,int n25,int n36,
341 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshNode * n1,
342 const SMDS_MeshNode * n2,
343 const SMDS_MeshNode * n3,
344 const SMDS_MeshNode * n4,
345 const SMDS_MeshNode * n5,
346 const SMDS_MeshNode * n6,
347 const SMDS_MeshNode * n12,
348 const SMDS_MeshNode * n23,
349 const SMDS_MeshNode * n31,
350 const SMDS_MeshNode * n45,
351 const SMDS_MeshNode * n56,
352 const SMDS_MeshNode * n64,
353 const SMDS_MeshNode * n14,
354 const SMDS_MeshNode * n25,
355 const SMDS_MeshNode * n36,
357 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshNode * n1,
358 const SMDS_MeshNode * n2,
359 const SMDS_MeshNode * n3,
360 const SMDS_MeshNode * n4,
361 const SMDS_MeshNode * n5,
362 const SMDS_MeshNode * n6,
363 const SMDS_MeshNode * n12,
364 const SMDS_MeshNode * n23,
365 const SMDS_MeshNode * n31,
366 const SMDS_MeshNode * n45,
367 const SMDS_MeshNode * n56,
368 const SMDS_MeshNode * n64,
369 const SMDS_MeshNode * n14,
370 const SMDS_MeshNode * n25,
371 const SMDS_MeshNode * n36);
373 // 2d oreder Hexahedrons with 20 nodes
374 virtual SMDS_MeshVolume* AddVolumeWithID(int n1, int n2, int n3, int n4,
375 int n5, int n6, int n7, int n8,
376 int n12,int n23,int n34,int n41,
377 int n56,int n67,int n78,int n85,
378 int n15,int n26,int n37,int n48,
380 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshNode * n1,
381 const SMDS_MeshNode * n2,
382 const SMDS_MeshNode * n3,
383 const SMDS_MeshNode * n4,
384 const SMDS_MeshNode * n5,
385 const SMDS_MeshNode * n6,
386 const SMDS_MeshNode * n7,
387 const SMDS_MeshNode * n8,
388 const SMDS_MeshNode * n12,
389 const SMDS_MeshNode * n23,
390 const SMDS_MeshNode * n34,
391 const SMDS_MeshNode * n41,
392 const SMDS_MeshNode * n56,
393 const SMDS_MeshNode * n67,
394 const SMDS_MeshNode * n78,
395 const SMDS_MeshNode * n85,
396 const SMDS_MeshNode * n15,
397 const SMDS_MeshNode * n26,
398 const SMDS_MeshNode * n37,
399 const SMDS_MeshNode * n48,
401 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshNode * n1,
402 const SMDS_MeshNode * n2,
403 const SMDS_MeshNode * n3,
404 const SMDS_MeshNode * n4,
405 const SMDS_MeshNode * n5,
406 const SMDS_MeshNode * n6,
407 const SMDS_MeshNode * n7,
408 const SMDS_MeshNode * n8,
409 const SMDS_MeshNode * n12,
410 const SMDS_MeshNode * n23,
411 const SMDS_MeshNode * n34,
412 const SMDS_MeshNode * n41,
413 const SMDS_MeshNode * n56,
414 const SMDS_MeshNode * n67,
415 const SMDS_MeshNode * n78,
416 const SMDS_MeshNode * n85,
417 const SMDS_MeshNode * n15,
418 const SMDS_MeshNode * n26,
419 const SMDS_MeshNode * n37,
420 const SMDS_MeshNode * n48);
422 virtual SMDS_MeshFace* AddPolygonalFaceWithID (std::vector<int> nodes_ids,
425 virtual SMDS_MeshFace* AddPolygonalFaceWithID (std::vector<const SMDS_MeshNode*> nodes,
428 virtual SMDS_MeshFace* AddPolygonalFace (std::vector<const SMDS_MeshNode*> nodes);
430 virtual SMDS_MeshVolume* AddPolyhedralVolumeWithID
431 (std::vector<int> nodes_ids,
432 std::vector<int> quantities,
435 virtual SMDS_MeshVolume* AddPolyhedralVolumeWithID
436 (std::vector<const SMDS_MeshNode*> nodes,
437 std::vector<int> quantities,
440 virtual SMDS_MeshVolume* AddPolyhedralVolume
441 (std::vector<const SMDS_MeshNode*> nodes,
442 std::vector<int> quantities);
444 virtual void RemoveElement(const SMDS_MeshElement * elem,
445 std::list<const SMDS_MeshElement *>& removedElems,
446 std::list<const SMDS_MeshElement *>& removedNodes,
447 const bool removenodes = false);
448 virtual void RemoveElement(const SMDS_MeshElement * elem, bool removenodes = false);
449 virtual void RemoveNode(const SMDS_MeshNode * node);
450 virtual void Remove0DElement(const SMDS_Mesh0DElement * elem0d);
451 virtual void RemoveEdge(const SMDS_MeshEdge * edge);
452 virtual void RemoveFace(const SMDS_MeshFace * face);
453 virtual void RemoveVolume(const SMDS_MeshVolume * volume);
455 /*! Remove only the given element and only if it is free.
456 * Method does not work for meshes with descendants.
457 * Implemented for fast cleaning of meshes.
459 virtual void RemoveFreeElement(const SMDS_MeshElement * elem);
461 virtual void Clear();
463 virtual bool RemoveFromParent();
464 virtual bool RemoveSubMesh(const SMDS_Mesh * aMesh);
466 bool ChangeElementNodes(const SMDS_MeshElement * elem,
467 const SMDS_MeshNode * nodes[],
469 bool ChangePolyhedronNodes(const SMDS_MeshElement * elem,
470 const std::vector<const SMDS_MeshNode*>& nodes,
471 const std::vector<int> & quantities);
473 virtual void Renumber (const bool isNodes, const int startID = 1, const int deltaID = 1);
474 // Renumber all nodes or elements.
476 const SMDS_MeshNode *FindNode(int idnode) const;
477 const SMDS_Mesh0DElement* Find0DElement(int idnode) const;
478 const SMDS_MeshEdge *FindEdge(int idnode1, int idnode2) const;
479 const SMDS_MeshEdge *FindEdge(int idnode1, int idnode2, int idnode3) const;
480 const SMDS_MeshFace *FindFace(int idnode1, int idnode2, int idnode3) const;
481 const SMDS_MeshFace *FindFace(int idnode1, int idnode2, int idnode3, int idnode4) const;
482 const SMDS_MeshFace *FindFace(int idnode1, int idnode2, int idnode3,
483 int idnode4, int idnode5, int idnode6) const;
484 const SMDS_MeshFace *FindFace(int idnode1, int idnode2, int idnode3, int idnode4,
485 int idnode5, int idnode6, int idnode7, int idnode8) const;
486 const SMDS_MeshElement *FindElement(int IDelem) const;
487 static const SMDS_Mesh0DElement* Find0DElement(const SMDS_MeshNode * n);
488 static const SMDS_MeshEdge* FindEdge(const SMDS_MeshNode * n1,
489 const SMDS_MeshNode * n2);
490 static const SMDS_MeshEdge* FindEdge(const SMDS_MeshNode * n1,
491 const SMDS_MeshNode * n2,
492 const SMDS_MeshNode * n3);
493 static const SMDS_MeshFace* FindFace(const SMDS_MeshNode *n1,
494 const SMDS_MeshNode *n2,
495 const SMDS_MeshNode *n3);
496 static const SMDS_MeshFace* FindFace(const SMDS_MeshNode *n1,
497 const SMDS_MeshNode *n2,
498 const SMDS_MeshNode *n3,
499 const SMDS_MeshNode *n4);
500 static const SMDS_MeshFace* FindFace(const SMDS_MeshNode *n1,
501 const SMDS_MeshNode *n2,
502 const SMDS_MeshNode *n3,
503 const SMDS_MeshNode *n4,
504 const SMDS_MeshNode *n5,
505 const SMDS_MeshNode *n6);
506 static const SMDS_MeshFace* FindFace(const SMDS_MeshNode *n1,
507 const SMDS_MeshNode *n2,
508 const SMDS_MeshNode *n3,
509 const SMDS_MeshNode *n4,
510 const SMDS_MeshNode *n5,
511 const SMDS_MeshNode *n6,
512 const SMDS_MeshNode *n7,
513 const SMDS_MeshNode *n8);
515 const SMDS_MeshFace *FindFace(std::vector<int> nodes_ids) const;
516 static const SMDS_MeshFace* FindFace(std::vector<const SMDS_MeshNode *> nodes);
519 * \brief Raise an exception if free memory (ram+swap) too low
520 * \param doNotRaise - if true, suppres exception, just return free memory size
521 * \retval int - amount of available memory in MB or negative number in failure case
523 static int CheckMemory(const bool doNotRaise=false) throw (std::bad_alloc);
525 int MaxNodeID() const;
526 int MinNodeID() const;
527 int MaxElementID() const;
528 int MinElementID() const;
530 const SMDS_MeshInfo& GetMeshInfo() const { return myInfo; }
533 int Nb0DElements() const;
536 int NbVolumes() const;
537 int NbSubMesh() const;
538 void DumpNodes() const;
539 void Dump0DElements() const;
540 void DumpEdges() const;
541 void DumpFaces() const;
542 void DumpVolumes() const;
543 void DebugStats() const;
544 SMDS_Mesh *boundaryFaces();
545 SMDS_Mesh *boundaryEdges();
546 virtual ~SMDS_Mesh();
547 bool hasConstructionEdges();
548 bool hasConstructionFaces();
549 bool hasInverseElements();
550 void setConstructionEdges(bool);
551 void setConstructionFaces(bool);
552 void setInverseElements(bool);
555 * Checks if the element is present in mesh.
556 * Useful to determine dead pointers.
557 * Use this function for debug purpose only! Do not check in the code
558 * using it even in _DEBUG_ mode
560 bool Contains (const SMDS_MeshElement* elem) const;
562 typedef std::vector<SMDS_MeshNode *> SetOfNodes;
563 typedef std::vector<SMDS_MeshCell *> SetOfCells;
565 void updateNodeMinMax();
566 inline int fromVtkToSmds(int vtkid) { MYASSERT(vtkid>=0); return myVtkIndex[vtkid]; };
567 inline int fromSmdsToVtk(int smdsid) { MYASSERT(smdsid>=0); return myIDElements[smdsid]; };
569 void incrementNodesCapacity(int nbNodes);
570 void incrementCellsCapacity(int nbCells);
574 static int chunkSize;
577 SMDS_Mesh(SMDS_Mesh * parent);
579 SMDS_MeshFace * createTriangle(const SMDS_MeshNode * node1,
580 const SMDS_MeshNode * node2,
581 const SMDS_MeshNode * node3,
583 SMDS_MeshFace * createQuadrangle(const SMDS_MeshNode * node1,
584 const SMDS_MeshNode * node2,
585 const SMDS_MeshNode * node3,
586 const SMDS_MeshNode * node4,
588 // SMDS_Mesh0DElement* Find0DElementOrCreate(const SMDS_MeshNode * n);
589 SMDS_MeshEdge* FindEdgeOrCreate(const SMDS_MeshNode * n1,
590 const SMDS_MeshNode * n2);
591 SMDS_MeshFace* FindFaceOrCreate(const SMDS_MeshNode *n1,
592 const SMDS_MeshNode *n2,
593 const SMDS_MeshNode *n3);
594 SMDS_MeshFace* FindFaceOrCreate(const SMDS_MeshNode *n1,
595 const SMDS_MeshNode *n2,
596 const SMDS_MeshNode *n3,
597 const SMDS_MeshNode *n4);
599 bool registerElement(int ID, SMDS_MeshElement * element);
601 void addChildrenWithNodes(std::set<const SMDS_MeshElement*>& setOfChildren,
602 const SMDS_MeshElement * element,
603 std::set<const SMDS_MeshElement*>& nodes);
605 inline void adjustmyCellsCapacity(int ID)
608 if (ID >= myCells.size())
609 myCells.resize(ID+SMDS_Mesh::chunkSize,0);
614 //! index of this SMDS_mesh in the static vector<SMDS_Mesh*> _meshList
617 //! actual nodes coordinates, cells definition and reverse connectivity are stored in a vtkUnstructuredGrid
618 vtkUnstructuredGrid* myGrid;
620 //! Small objects like SMDS_MeshNode are allocated by chunks to limit memory costs of new
621 ObjectPool<SMDS_MeshNode>* myNodePool;
623 //! Small objects like SMDS_VtkVolume are allocated by chunks to limit memory costs of new
624 ObjectPool<SMDS_VtkVolume>* myVolumePool;
625 ObjectPool<SMDS_VtkFace>* myFacePool;
626 ObjectPool<SMDS_VtkEdge>* myEdgePool;
628 //! SMDS_MeshNodes refer to vtk nodes (vtk id = index in myNodes),store reference to this mesh, and subshape
631 //! SMDS_MeshCells refer to vtk cells (vtk id != index in myCells),store reference to this mesh, and subshape
634 //! for cells only: index = ID for SMDS users, value = ID in vtkUnstructuredGrid
635 std::vector<int> myIDElements;
637 //! for cells only: index = ID in vtkUnstructuredGrid, value = ID for SMDS users
638 std::vector<int> myVtkIndex;
640 SMDS_Mesh * myParent;
641 std::list<SMDS_Mesh *> myChildren;
642 SMDS_MeshNodeIDFactory *myNodeIDFactory;
643 SMDS_MeshElementIDFactory *myElementIDFactory;
644 SMDS_MeshInfo myInfo;
649 bool myHasConstructionEdges;
650 bool myHasConstructionFaces;
651 bool myHasInverseElements;