Salome HOME
PR: double nodes and flat elements for ASTER calculations in progress
[modules/smesh.git] / src / SMDS / SMDS_UnstructuredGrid.hxx
1 /* 
2  * File:   SMDS_UnstructuredGrid.hxx
3  * Author: prascle
4  *
5  * Created on September 16, 2009, 10:28 PM
6  */
7
8 #ifndef _SMDS_UNSTRUCTUREDGRID_HXX
9 #define _SMDS_UNSTRUCTUREDGRID_HXX
10
11 #include <vtkUnstructuredGrid.h>
12 #include <vtkCellLinks.h>
13 #include "chrono.hxx"
14
15 #include <vector>
16 #include <set>
17 #include <map>
18
19 //#define VTK_HAVE_POLYHEDRON
20 //#ifdef VTK_HAVE_POLYHEDRON
21   #define VTK_MAXTYPE VTK_POLYHEDRON
22 //#else
23 //  #define VTK_MAXTYPE VTK_QUADRATIC_PYRAMID
24 //#endif
25
26 #define NBMAXNEIGHBORS 100
27
28 // allow very huge polyhedrons in tests
29 #define NBMAXNODESINCELL 5000
30
31 class SMDS_Downward;
32 class SMDS_Mesh;
33
34 class SMDS_CellLinks: public vtkCellLinks
35 {
36 public:
37   vtkCellLinks::Link* ResizeL(vtkIdType sz);
38   vtkIdType GetLinksSize();
39   static SMDS_CellLinks* New();
40 protected:
41   SMDS_CellLinks();
42   ~SMDS_CellLinks();
43 };
44
45 class SMDS_UnstructuredGrid: public vtkUnstructuredGrid
46 {
47 public:
48   void setSMDS_mesh(SMDS_Mesh *mesh);
49   void compactGrid(std::vector<int>& idNodesOldToNew, int newNodeSize, std::vector<int>& idCellsOldToNew,
50                    int newCellSize);
51
52   virtual unsigned long GetMTime();
53   virtual void Update();
54   virtual void UpdateInformation();
55   virtual vtkPoints *GetPoints();
56
57 //#ifdef VTK_HAVE_POLYHEDRON
58   int InsertNextLinkedCell(int type, int npts, vtkIdType *pts);
59 //#endif
60
61   int CellIdToDownId(int vtkCellId);
62   void setCellIdToDownId(int vtkCellId, int downId);
63   void BuildDownwardConnectivity(bool withEdges);
64   int GetNeighbors(int* neighborsVtkIds, int* downIds, unsigned char* downTypes, int vtkId);
65   int GetParentVolumes(int* volVtkIds, int vtkId);
66   int GetParentVolumes(int* volVtkIds, int downId, unsigned char downType);
67   void GetNodeIds(std::set<int>& nodeSet, int downId, unsigned char downType);
68   void ModifyCellNodes(int vtkVolId, std::map<int, int> localClonedNodeIds);
69   int getOrderedNodesOfFace(int vtkVolId, std::vector<vtkIdType>& orderedNodes);
70   void BuildLinks();
71   bool extrudeVolumeFromFace(int vtkVolId,
72                              int domain1,
73                              int domain2,
74                              std::set<int>& originalNodes,
75                              std::map<int,std::map<int,int> >& nodeDomains,
76                              std::map<int,std::map<long,int> >& nodeQuadDomains);
77   vtkCellLinks* GetLinks()
78   {
79     return Links;
80   }
81   SMDS_Downward* getDownArray(unsigned char vtkType)
82   {
83     return _downArray[vtkType];
84   }
85   static SMDS_UnstructuredGrid* New();
86   SMDS_Mesh *_mesh;
87 protected:
88   SMDS_UnstructuredGrid();
89   ~SMDS_UnstructuredGrid();
90   void copyNodes(vtkPoints *newPoints, std::vector<int>& idNodesOldToNew, int& alreadyCopied, int start, int end);
91   void copyBloc(vtkUnsignedCharArray *newTypes, std::vector<int>& idCellsOldToNew, std::vector<int>& idNodesOldToNew,
92                 vtkCellArray* newConnectivity, vtkIdTypeArray* newLocations, vtkIdType* pointsCell, int& alreadyCopied,
93                 int start, int end);
94
95   std::vector<int> _cellIdToDownId; //!< convert vtk Id to downward[vtkType] id, initialized with -1
96   std::vector<unsigned char> _downTypes;
97   std::vector<SMDS_Downward*> _downArray;
98 };
99
100 #endif  /* _SMDS_UNSTRUCTUREDGRID_HXX */
101