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_Mesh0DElement.hxx"
33 #include "SMDS_MeshEdge.hxx"
34 #include "SMDS_MeshFace.hxx"
35 #include "SMDS_MeshVolume.hxx"
36 #include "SMDS_MeshNodeIDFactory.hxx"
37 #include "SMDS_MeshElementIDFactory.hxx"
38 #include "SMDS_MeshInfo.hxx"
39 #include "SMDS_ElemIterator.hxx"
41 #include <boost/shared_ptr.hpp>
45 #include <vtkSystemIncludes.h>
47 class vtkUnstructuredGrid;
49 class SMDS_EXPORT SMDS_Mesh:public SMDS_MeshObject{
52 static std::vector<SMDS_Mesh*> _meshList; // --- to find the SMDS_mesh from its elements
56 inline vtkUnstructuredGrid* getGrid() {return myGrid; };
58 SMDS_NodeIteratorPtr nodesIterator() const;
59 SMDS_0DElementIteratorPtr elements0dIterator() const;
60 SMDS_EdgeIteratorPtr edgesIterator() const;
61 SMDS_FaceIteratorPtr facesIterator() const;
62 SMDS_VolumeIteratorPtr volumesIterator() const;
63 SMDS_ElemIteratorPtr elementsIterator(SMDSAbs_ElementType type=SMDSAbs_All) const;
65 SMDSAbs_ElementType GetElementType( const int id, const bool iselem ) const;
67 SMDS_Mesh *AddSubMesh();
69 virtual SMDS_MeshNode* AddNodeWithID(double x, double y, double z, int ID);
70 virtual SMDS_MeshNode* AddNode(double x, double y, double z);
72 virtual SMDS_Mesh0DElement* Add0DElementWithID(int n, int ID);
73 virtual SMDS_Mesh0DElement* Add0DElementWithID(const SMDS_MeshNode * n, int ID);
74 virtual SMDS_Mesh0DElement* Add0DElement (const SMDS_MeshNode * n);
76 virtual SMDS_MeshEdge* AddEdgeWithID(int n1, int n2, int ID);
77 virtual SMDS_MeshEdge* AddEdgeWithID(const SMDS_MeshNode * n1,
78 const SMDS_MeshNode * n2,
80 virtual SMDS_MeshEdge* AddEdge(const SMDS_MeshNode * n1,
81 const SMDS_MeshNode * n2);
83 // 2d order edge with 3 nodes: n12 - node between n1 and n2
84 virtual SMDS_MeshEdge* AddEdgeWithID(int n1, int n2, int n12, int ID);
85 virtual SMDS_MeshEdge* AddEdgeWithID(const SMDS_MeshNode * n1,
86 const SMDS_MeshNode * n2,
87 const SMDS_MeshNode * n12,
89 virtual SMDS_MeshEdge* AddEdge(const SMDS_MeshNode * n1,
90 const SMDS_MeshNode * n2,
91 const SMDS_MeshNode * n12);
93 virtual SMDS_MeshFace* AddFaceWithID(int n1, int n2, int n3, int ID);
94 virtual SMDS_MeshFace* AddFaceWithID(const SMDS_MeshNode * n1,
95 const SMDS_MeshNode * n2,
96 const SMDS_MeshNode * n3,
98 virtual SMDS_MeshFace* AddFace(const SMDS_MeshNode * n1,
99 const SMDS_MeshNode * n2,
100 const SMDS_MeshNode * n3);
102 virtual SMDS_MeshFace* AddFaceWithID(int n1, int n2, int n3, int n4, int ID);
103 virtual SMDS_MeshFace* AddFaceWithID(const SMDS_MeshNode * n1,
104 const SMDS_MeshNode * n2,
105 const SMDS_MeshNode * n3,
106 const SMDS_MeshNode * n4,
108 virtual SMDS_MeshFace* AddFace(const SMDS_MeshNode * n1,
109 const SMDS_MeshNode * n2,
110 const SMDS_MeshNode * n3,
111 const SMDS_MeshNode * n4);
113 virtual SMDS_MeshFace* AddFaceWithID(const SMDS_MeshEdge * e1,
114 const SMDS_MeshEdge * e2,
115 const SMDS_MeshEdge * e3, int ID);
116 virtual SMDS_MeshFace* AddFace(const SMDS_MeshEdge * e1,
117 const SMDS_MeshEdge * e2,
118 const SMDS_MeshEdge * e3);
120 virtual SMDS_MeshFace* AddFaceWithID(const SMDS_MeshEdge * e1,
121 const SMDS_MeshEdge * e2,
122 const SMDS_MeshEdge * e3,
123 const SMDS_MeshEdge * e4, int ID);
124 virtual SMDS_MeshFace* AddFace(const SMDS_MeshEdge * e1,
125 const SMDS_MeshEdge * e2,
126 const SMDS_MeshEdge * e3,
127 const SMDS_MeshEdge * e4);
129 // 2d order triangle of 6 nodes
130 virtual SMDS_MeshFace* AddFaceWithID(int n1, int n2, int n3,
131 int n12,int n23,int n31, int ID);
132 virtual SMDS_MeshFace* AddFaceWithID(const SMDS_MeshNode * n1,
133 const SMDS_MeshNode * n2,
134 const SMDS_MeshNode * n3,
135 const SMDS_MeshNode * n12,
136 const SMDS_MeshNode * n23,
137 const SMDS_MeshNode * n31,
139 virtual SMDS_MeshFace* AddFace(const SMDS_MeshNode * n1,
140 const SMDS_MeshNode * n2,
141 const SMDS_MeshNode * n3,
142 const SMDS_MeshNode * n12,
143 const SMDS_MeshNode * n23,
144 const SMDS_MeshNode * n31);
146 // 2d order quadrangle
147 virtual SMDS_MeshFace* AddFaceWithID(int n1, int n2, int n3, int n4,
148 int n12,int n23,int n34,int n41, int ID);
149 virtual SMDS_MeshFace* AddFaceWithID(const SMDS_MeshNode * n1,
150 const SMDS_MeshNode * n2,
151 const SMDS_MeshNode * n3,
152 const SMDS_MeshNode * n4,
153 const SMDS_MeshNode * n12,
154 const SMDS_MeshNode * n23,
155 const SMDS_MeshNode * n34,
156 const SMDS_MeshNode * n41,
158 virtual SMDS_MeshFace* AddFace(const SMDS_MeshNode * n1,
159 const SMDS_MeshNode * n2,
160 const SMDS_MeshNode * n3,
161 const SMDS_MeshNode * n4,
162 const SMDS_MeshNode * n12,
163 const SMDS_MeshNode * n23,
164 const SMDS_MeshNode * n34,
165 const SMDS_MeshNode * n41);
167 virtual SMDS_MeshVolume* AddVolumeWithID(int n1, int n2, int n3, int n4, int ID);
168 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshNode * n1,
169 const SMDS_MeshNode * n2,
170 const SMDS_MeshNode * n3,
171 const SMDS_MeshNode * n4,
173 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshNode * n1,
174 const SMDS_MeshNode * n2,
175 const SMDS_MeshNode * n3,
176 const SMDS_MeshNode * n4);
178 virtual SMDS_MeshVolume* AddVolumeWithID(int n1, int n2, int n3, int n4,
180 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshNode * n1,
181 const SMDS_MeshNode * n2,
182 const SMDS_MeshNode * n3,
183 const SMDS_MeshNode * n4,
184 const SMDS_MeshNode * n5,
186 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshNode * n1,
187 const SMDS_MeshNode * n2,
188 const SMDS_MeshNode * n3,
189 const SMDS_MeshNode * n4,
190 const SMDS_MeshNode * n5);
192 virtual SMDS_MeshVolume* AddVolumeWithID(int n1, int n2, int n3, int n4,
193 int n5, int n6, int ID);
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,
199 const SMDS_MeshNode * n6,
201 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshNode * n1,
202 const SMDS_MeshNode * n2,
203 const SMDS_MeshNode * n3,
204 const SMDS_MeshNode * n4,
205 const SMDS_MeshNode * n5,
206 const SMDS_MeshNode * n6);
208 virtual SMDS_MeshVolume* AddVolumeWithID(int n1, int n2, int n3, int n4,
209 int n5, int n6, int n7, int n8, int ID);
210 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshNode * n1,
211 const SMDS_MeshNode * n2,
212 const SMDS_MeshNode * n3,
213 const SMDS_MeshNode * n4,
214 const SMDS_MeshNode * n5,
215 const SMDS_MeshNode * n6,
216 const SMDS_MeshNode * n7,
217 const SMDS_MeshNode * n8,
219 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshNode * n1,
220 const SMDS_MeshNode * n2,
221 const SMDS_MeshNode * n3,
222 const SMDS_MeshNode * n4,
223 const SMDS_MeshNode * n5,
224 const SMDS_MeshNode * n6,
225 const SMDS_MeshNode * n7,
226 const SMDS_MeshNode * n8);
228 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshFace * f1,
229 const SMDS_MeshFace * f2,
230 const SMDS_MeshFace * f3,
231 const SMDS_MeshFace * f4, int ID);
232 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshFace * f1,
233 const SMDS_MeshFace * f2,
234 const SMDS_MeshFace * f3,
235 const SMDS_MeshFace * f4);
237 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshFace * f1,
238 const SMDS_MeshFace * f2,
239 const SMDS_MeshFace * f3,
240 const SMDS_MeshFace * f4,
241 const SMDS_MeshFace * f5, int ID);
242 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshFace * f1,
243 const SMDS_MeshFace * f2,
244 const SMDS_MeshFace * f3,
245 const SMDS_MeshFace * f4,
246 const SMDS_MeshFace * f5);
248 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshFace * f1,
249 const SMDS_MeshFace * f2,
250 const SMDS_MeshFace * f3,
251 const SMDS_MeshFace * f4,
252 const SMDS_MeshFace * f5,
253 const SMDS_MeshFace * f6, int ID);
254 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshFace * f1,
255 const SMDS_MeshFace * f2,
256 const SMDS_MeshFace * f3,
257 const SMDS_MeshFace * f4,
258 const SMDS_MeshFace * f5,
259 const SMDS_MeshFace * f6);
261 // 2d order tetrahedron of 10 nodes
262 virtual SMDS_MeshVolume* AddVolumeWithID(int n1, int n2, int n3, int n4,
263 int n12,int n23,int n31,
264 int n14,int n24,int n34, int ID);
265 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshNode * n1,
266 const SMDS_MeshNode * n2,
267 const SMDS_MeshNode * n3,
268 const SMDS_MeshNode * n4,
269 const SMDS_MeshNode * n12,
270 const SMDS_MeshNode * n23,
271 const SMDS_MeshNode * n31,
272 const SMDS_MeshNode * n14,
273 const SMDS_MeshNode * n24,
274 const SMDS_MeshNode * n34,
276 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshNode * n1,
277 const SMDS_MeshNode * n2,
278 const SMDS_MeshNode * n3,
279 const SMDS_MeshNode * n4,
280 const SMDS_MeshNode * n12,
281 const SMDS_MeshNode * n23,
282 const SMDS_MeshNode * n31,
283 const SMDS_MeshNode * n14,
284 const SMDS_MeshNode * n24,
285 const SMDS_MeshNode * n34);
287 // 2d order pyramid of 13 nodes
288 virtual SMDS_MeshVolume* AddVolumeWithID(int n1, int n2, int n3, int n4, int n5,
289 int n12,int n23,int n34,int n41,
290 int n15,int n25,int n35,int n45,
292 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshNode * n1,
293 const SMDS_MeshNode * n2,
294 const SMDS_MeshNode * n3,
295 const SMDS_MeshNode * n4,
296 const SMDS_MeshNode * n5,
297 const SMDS_MeshNode * n12,
298 const SMDS_MeshNode * n23,
299 const SMDS_MeshNode * n34,
300 const SMDS_MeshNode * n41,
301 const SMDS_MeshNode * n15,
302 const SMDS_MeshNode * n25,
303 const SMDS_MeshNode * n35,
304 const SMDS_MeshNode * n45,
306 virtual SMDS_MeshVolume* AddVolume(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 // 2d order Pentahedron with 15 nodes
321 virtual SMDS_MeshVolume* AddVolumeWithID(int n1, int n2, int n3,
322 int n4, int n5, int n6,
323 int n12,int n23,int n31,
324 int n45,int n56,int n64,
325 int n14,int n25,int n36,
327 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshNode * n1,
328 const SMDS_MeshNode * n2,
329 const SMDS_MeshNode * n3,
330 const SMDS_MeshNode * n4,
331 const SMDS_MeshNode * n5,
332 const SMDS_MeshNode * n6,
333 const SMDS_MeshNode * n12,
334 const SMDS_MeshNode * n23,
335 const SMDS_MeshNode * n31,
336 const SMDS_MeshNode * n45,
337 const SMDS_MeshNode * n56,
338 const SMDS_MeshNode * n64,
339 const SMDS_MeshNode * n14,
340 const SMDS_MeshNode * n25,
341 const SMDS_MeshNode * n36,
343 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshNode * n1,
344 const SMDS_MeshNode * n2,
345 const SMDS_MeshNode * n3,
346 const SMDS_MeshNode * n4,
347 const SMDS_MeshNode * n5,
348 const SMDS_MeshNode * n6,
349 const SMDS_MeshNode * n12,
350 const SMDS_MeshNode * n23,
351 const SMDS_MeshNode * n31,
352 const SMDS_MeshNode * n45,
353 const SMDS_MeshNode * n56,
354 const SMDS_MeshNode * n64,
355 const SMDS_MeshNode * n14,
356 const SMDS_MeshNode * n25,
357 const SMDS_MeshNode * n36);
359 // 2d oreder Hexahedrons with 20 nodes
360 virtual SMDS_MeshVolume* AddVolumeWithID(int n1, int n2, int n3, int n4,
361 int n5, int n6, int n7, int n8,
362 int n12,int n23,int n34,int n41,
363 int n56,int n67,int n78,int n85,
364 int n15,int n26,int n37,int n48,
366 virtual SMDS_MeshVolume* AddVolumeWithID(const SMDS_MeshNode * n1,
367 const SMDS_MeshNode * n2,
368 const SMDS_MeshNode * n3,
369 const SMDS_MeshNode * n4,
370 const SMDS_MeshNode * n5,
371 const SMDS_MeshNode * n6,
372 const SMDS_MeshNode * n7,
373 const SMDS_MeshNode * n8,
374 const SMDS_MeshNode * n12,
375 const SMDS_MeshNode * n23,
376 const SMDS_MeshNode * n34,
377 const SMDS_MeshNode * n41,
378 const SMDS_MeshNode * n56,
379 const SMDS_MeshNode * n67,
380 const SMDS_MeshNode * n78,
381 const SMDS_MeshNode * n85,
382 const SMDS_MeshNode * n15,
383 const SMDS_MeshNode * n26,
384 const SMDS_MeshNode * n37,
385 const SMDS_MeshNode * n48,
387 virtual SMDS_MeshVolume* AddVolume(const SMDS_MeshNode * n1,
388 const SMDS_MeshNode * n2,
389 const SMDS_MeshNode * n3,
390 const SMDS_MeshNode * n4,
391 const SMDS_MeshNode * n5,
392 const SMDS_MeshNode * n6,
393 const SMDS_MeshNode * n7,
394 const SMDS_MeshNode * n8,
395 const SMDS_MeshNode * n12,
396 const SMDS_MeshNode * n23,
397 const SMDS_MeshNode * n34,
398 const SMDS_MeshNode * n41,
399 const SMDS_MeshNode * n56,
400 const SMDS_MeshNode * n67,
401 const SMDS_MeshNode * n78,
402 const SMDS_MeshNode * n85,
403 const SMDS_MeshNode * n15,
404 const SMDS_MeshNode * n26,
405 const SMDS_MeshNode * n37,
406 const SMDS_MeshNode * n48);
408 virtual SMDS_MeshFace* AddPolygonalFaceWithID (std::vector<int> nodes_ids,
411 virtual SMDS_MeshFace* AddPolygonalFaceWithID (std::vector<const SMDS_MeshNode*> nodes,
414 virtual SMDS_MeshFace* AddPolygonalFace (std::vector<const SMDS_MeshNode*> nodes);
416 virtual SMDS_MeshVolume* AddPolyhedralVolumeWithID
417 (std::vector<int> nodes_ids,
418 std::vector<int> quantities,
421 virtual SMDS_MeshVolume* AddPolyhedralVolumeWithID
422 (std::vector<const SMDS_MeshNode*> nodes,
423 std::vector<int> quantities,
426 virtual SMDS_MeshVolume* AddPolyhedralVolume
427 (std::vector<const SMDS_MeshNode*> nodes,
428 std::vector<int> quantities);
430 virtual void RemoveElement(const SMDS_MeshElement * elem,
431 std::list<const SMDS_MeshElement *>& removedElems,
432 std::list<const SMDS_MeshElement *>& removedNodes,
433 const bool removenodes = false);
434 virtual void RemoveElement(const SMDS_MeshElement * elem, bool removenodes = false);
435 virtual void RemoveNode(const SMDS_MeshNode * node);
436 virtual void Remove0DElement(const SMDS_Mesh0DElement * elem0d);
437 virtual void RemoveEdge(const SMDS_MeshEdge * edge);
438 virtual void RemoveFace(const SMDS_MeshFace * face);
439 virtual void RemoveVolume(const SMDS_MeshVolume * volume);
441 /*! Remove only the given element and only if it is free.
442 * Method does not work for meshes with descendants.
443 * Implemented for fast cleaning of meshes.
445 virtual void RemoveFreeElement(const SMDS_MeshElement * elem);
447 virtual void Clear();
449 virtual bool RemoveFromParent();
450 virtual bool RemoveSubMesh(const SMDS_Mesh * aMesh);
452 bool ChangeElementNodes(const SMDS_MeshElement * elem,
453 const SMDS_MeshNode * nodes[],
455 bool ChangePolyhedronNodes(const SMDS_MeshElement * elem,
456 const std::vector<const SMDS_MeshNode*>& nodes,
457 const std::vector<int> & quantities);
459 virtual void Renumber (const bool isNodes, const int startID = 1, const int deltaID = 1);
460 // Renumber all nodes or elements.
462 const SMDS_MeshNode *FindNode(int idnode) const;
463 const SMDS_Mesh0DElement* Find0DElement(int idnode) const;
464 const SMDS_MeshEdge *FindEdge(int idnode1, int idnode2) const;
465 const SMDS_MeshEdge *FindEdge(int idnode1, int idnode2, int idnode3) const;
466 const SMDS_MeshFace *FindFace(int idnode1, int idnode2, int idnode3) const;
467 const SMDS_MeshFace *FindFace(int idnode1, int idnode2, int idnode3, int idnode4) const;
468 const SMDS_MeshFace *FindFace(int idnode1, int idnode2, int idnode3,
469 int idnode4, int idnode5, int idnode6) const;
470 const SMDS_MeshFace *FindFace(int idnode1, int idnode2, int idnode3, int idnode4,
471 int idnode5, int idnode6, int idnode7, int idnode8) const;
472 const SMDS_MeshElement *FindElement(int IDelem) const;
473 static const SMDS_Mesh0DElement* Find0DElement(const SMDS_MeshNode * n);
474 static const SMDS_MeshEdge* FindEdge(const SMDS_MeshNode * n1,
475 const SMDS_MeshNode * n2);
476 static const SMDS_MeshEdge* FindEdge(const SMDS_MeshNode * n1,
477 const SMDS_MeshNode * n2,
478 const SMDS_MeshNode * n3);
479 static const SMDS_MeshFace* FindFace(const SMDS_MeshNode *n1,
480 const SMDS_MeshNode *n2,
481 const SMDS_MeshNode *n3);
482 static const SMDS_MeshFace* FindFace(const SMDS_MeshNode *n1,
483 const SMDS_MeshNode *n2,
484 const SMDS_MeshNode *n3,
485 const SMDS_MeshNode *n4);
486 static const SMDS_MeshFace* FindFace(const SMDS_MeshNode *n1,
487 const SMDS_MeshNode *n2,
488 const SMDS_MeshNode *n3,
489 const SMDS_MeshNode *n4,
490 const SMDS_MeshNode *n5,
491 const SMDS_MeshNode *n6);
492 static const SMDS_MeshFace* FindFace(const SMDS_MeshNode *n1,
493 const SMDS_MeshNode *n2,
494 const SMDS_MeshNode *n3,
495 const SMDS_MeshNode *n4,
496 const SMDS_MeshNode *n5,
497 const SMDS_MeshNode *n6,
498 const SMDS_MeshNode *n7,
499 const SMDS_MeshNode *n8);
501 const SMDS_MeshFace *FindFace(std::vector<int> nodes_ids) const;
502 static const SMDS_MeshFace* FindFace(std::vector<const SMDS_MeshNode *> nodes);
505 * \brief Raise an exception if free memory (ram+swap) too low
506 * \param doNotRaise - if true, suppres exception, just return free memory size
507 * \retval int - amount of available memory in MB or negative number in failure case
509 static int CheckMemory(const bool doNotRaise=false) throw (std::bad_alloc);
511 int MaxNodeID() const;
512 int MinNodeID() const;
513 int MaxElementID() const;
514 int MinElementID() const;
516 const SMDS_MeshInfo& GetMeshInfo() const { return myInfo; }
519 int Nb0DElements() const;
522 int NbVolumes() const;
523 int NbSubMesh() const;
524 void DumpNodes() const;
525 void Dump0DElements() const;
526 void DumpEdges() const;
527 void DumpFaces() const;
528 void DumpVolumes() const;
529 void DebugStats() const;
530 SMDS_Mesh *boundaryFaces();
531 SMDS_Mesh *boundaryEdges();
532 virtual ~SMDS_Mesh();
533 bool hasConstructionEdges();
534 bool hasConstructionFaces();
535 bool hasInverseElements();
536 void setConstructionEdges(bool);
537 void setConstructionFaces(bool);
538 void setInverseElements(bool);
541 * Checks if the element is present in mesh.
542 * Useful to determine dead pointers.
543 * Use this function for debug purpose only! Do not check in the code
544 * using it even in _DEBUG_ mode
546 bool Contains (const SMDS_MeshElement* elem) const;
548 typedef std::vector<SMDS_MeshNode *> SetOfNodes;
549 typedef std::vector<SMDS_MeshCell *> SetOfCells;
551 void updateNodeMinMax();
552 int fromVtkToSmds(int vtkid) { return myElementIDFactory->fromVtkToSmds(vtkid); };
556 static int chunkSize;
559 SMDS_Mesh(SMDS_Mesh * parent);
561 SMDS_MeshFace * createTriangle(const SMDS_MeshNode * node1,
562 const SMDS_MeshNode * node2,
563 const SMDS_MeshNode * node3);
564 SMDS_MeshFace * createQuadrangle(const SMDS_MeshNode * node1,
565 const SMDS_MeshNode * node2,
566 const SMDS_MeshNode * node3,
567 const SMDS_MeshNode * node4,
569 // SMDS_Mesh0DElement* Find0DElementOrCreate(const SMDS_MeshNode * n);
570 SMDS_MeshEdge* FindEdgeOrCreate(const SMDS_MeshNode * n1,
571 const SMDS_MeshNode * n2);
572 SMDS_MeshFace* FindFaceOrCreate(const SMDS_MeshNode *n1,
573 const SMDS_MeshNode *n2,
574 const SMDS_MeshNode *n3);
575 SMDS_MeshFace* FindFaceOrCreate(const SMDS_MeshNode *n1,
576 const SMDS_MeshNode *n2,
577 const SMDS_MeshNode *n3,
578 const SMDS_MeshNode *n4);
580 bool registerElement(int ID, SMDS_MeshElement * element);
582 void addChildrenWithNodes(std::set<const SMDS_MeshElement*>& setOfChildren,
583 const SMDS_MeshElement * element,
584 std::set<const SMDS_MeshElement*>& nodes);
586 inline void adjustmyCellsCapacity(int ID)
589 if (ID >= myCells.size())
590 myCells.resize(ID+SMDS_Mesh::chunkSize,0);
595 int myMeshId; // --- index for this mesh in the vector
596 vtkUnstructuredGrid* myGrid;
602 // SetOf0DElements my0DElements;
603 // SetOfEdges myEdges;
604 // SetOfFaces myFaces;
605 // SetOfVolumes myVolumes;
607 SMDS_Mesh * myParent;
608 std::list<SMDS_Mesh *> myChildren;
609 SMDS_MeshNodeIDFactory *myNodeIDFactory;
610 SMDS_MeshElementIDFactory *myElementIDFactory;
611 SMDS_MeshInfo myInfo;
616 bool myHasConstructionEdges;
617 bool myHasConstructionFaces;
618 bool myHasInverseElements;