677d9223e6e7f3ee130be406740d29f041477f66
[charm.git] / src / libs / ck-libs / ParFUM / adapt_adj.h
1 /* Adaptivity Adjacencies: element-to-element adjacencies for use by
2    adaptivity codes only.  Adaptivity codes should keep these
3    up-to-date for each mesh modification primitive.
4
5    Created 11 Sept 2006 - Terry L. Wilmarth
6    
7    Format of adjacency information: (partitionID, localID on partition, type)
8    
9    2D Adjacencies: Adjacencies are across edges
10    3D Adjacencies: Adjacencies can be across edges and faces
11
12    Numberings:
13    
14    TRIANGLES: Given nodes 0, 1, 2, the edges 0, 1, and 2 of a triangle are:
15    (0, 1), (1, 2) and (2, 0), in that order.
16
17               0
18              / \
19             /   \
20           2/     \0
21           /       \
22          /         \
23         2-----------1
24               1
25    
26    QUADS: Given nodes 0, 1, 2, 3, the edges 0, 1, 2 and 3 of a quad are:
27    (0, 1), (1, 2), (2, 3) and (3, 0), in that order.
28    
29    3D Adjacencies: 
30    
31    TETS: Given nodes 0, 1, 2, 3, the faces 0, 1, 2 and 3 of a tetrahedra are:
32    (0, 1, 2), (0, 3, 1), (0, 2, 3), and (1,3,2), in that order
33
34               0
35              /|\
36             / | \
37           1/  0  \2       Back Face: 2
38           / 0 | 1 \       Bottom Face: 3
39          /....|....\       
40         2-_   |   _-3     Back Edge: 5
41            -_ | _-  
42           3  -1-  4
43
44    
45    HEXES: Given nodes 0, 1, 2, 3, 4, 5, 6, 7, the faces 0, 1, 2, 3, 4, 5
46    of a hex are: (0, 1, 2, 3), (1, 5, 6, 2), (2, 6, 7, 3), (3, 7, 4, 0),
47    (0, 4, 5, 1), (5, 4, 6, 7) in that order
48    
49 */
50
51 // NOTE: review for mixed and cohesive element handling
52 #include <set>
53 #include <algorithm>
54
55 #define MAX_ADJELEMS 6
56 #define MAX_NODESET_SIZE 4
57
58 // Each instance of adaptAdj represents an element to 
59 // element adjacency
60 class adaptAdj{
61  public:
62   int partID;   // partition ID
63   int localID;  // local entity ID on partition partID
64   int elemType; // element type (tri, quad, tet, hex, etc.)
65   adaptAdj():partID(-1),localID(-1),elemType(-1){};
66   adaptAdj(int _partID,int _localID,int _elemType) : partID(_partID), localID(_localID), elemType(_elemType){};
67   inline adaptAdj &operator=(const adaptAdj &rhs){
68     partID = rhs.partID;
69     localID = rhs.localID;
70     elemType = rhs.elemType;
71     return *this;
72   }
73   inline bool operator==(const adaptAdj &rhs) const{
74     return (partID==rhs.partID && localID==rhs.localID && elemType==rhs.elemType);
75   }
76   virtual void pup(PUP::er &p){
77     p | partID;
78     p | localID;
79     p | elemType;
80   }
81 };
82
83 // Each adjElem describes an adjacency by enumerating
84 // the nodes that form the "edge" shared by two 
85 // adjacent elements
86 class adjElem { // list entry for an element incident on a node
87  public:
88   int elemID; // local element id
89   int nodeSetID; // which nodeSet in nodeSetMap does this nodeSet refer to
90   CkVec<int> nodeSet; //local node ids
91   adjElem *next;
92   adjElem(int nodeSetSize) : nodeSet(nodeSetSize){};
93 };
94
95 class adjNode { // struct to store each node's adjacency info
96  public:        
97   int *sharedWithPartition; // array of partition IDs on which there is a corresponding
98                             // shared node; this is NULL if this is not a shared node
99   int *sharedWithLocalIdx;  // local Idx in idxl list with the corresponding chunk in sharedWithPartition
100   int numSharedPartitions;
101   int adjElemCount;         // number of entries in adjElemList (below)
102   // max length of adjElemList is 2*nodal degree
103   adjElem *adjElemList;     // list of elems incident on this node
104   adjNode() { 
105     sharedWithPartition = NULL;
106     adjElemList = new adjElem(0); // Create a dummy head node in the list
107     adjElemList->elemID = -1;
108     adjElemList->next = NULL;
109     adjElemCount = 0; 
110     numSharedPartitions=0;
111   }
112   ~adjNode(){ delete [] sharedWithPartition; delete [] sharedWithLocalIdx;}
113 };
114
115 class adjRequest{
116  public:
117   int elemID,chunkID,elemType,nodeSetID;
118   int translatedNodeSet[MAX_NODESET_SIZE];
119   adjRequest(): elemID(-1),chunkID(-1),elemType(-1){};
120   adjRequest(int _elemID,int _chunkID,int _nodeSetID,int _elemType ): elemID(_elemID),chunkID(_chunkID),nodeSetID(_nodeSetID), elemType(_elemType){};
121   adjRequest(const adjRequest &rhs){
122     *this = rhs;
123   }
124   inline adjRequest& operator=(const adjRequest &rhs){
125     elemID = rhs.elemID;
126     chunkID = rhs.chunkID;
127     elemType = rhs.elemType;
128     nodeSetID = rhs.nodeSetID;
129     memcpy(&translatedNodeSet[0],&(rhs.translatedNodeSet[0]),MAX_NODESET_SIZE*sizeof(int));
130     return *this;
131   }
132   virtual void pup(PUP::er &p){
133     p | elemID;
134     p | chunkID;
135     p | elemType;
136     p | nodeSetID;
137     p(translatedNodeSet,MAX_NODESET_SIZE);
138   }
139 };
140
141 class adjReply {
142  public:
143   int requestingElemID,requestingNodeSetID;
144   adaptAdj replyingElem;
145   adjReply(): requestingElemID(-1),requestingNodeSetID(-1), replyingElem(){};
146   adjReply(const adjReply &rhs){
147     *this = rhs;
148   }
149   
150   inline adjReply& operator=(const adjReply &rhs){
151     requestingElemID = rhs.requestingElemID;
152     requestingNodeSetID = rhs.requestingNodeSetID;
153     replyingElem = rhs.replyingElem;
154                 return *this;
155   }
156   virtual void pup(PUP::er &p){
157     p | requestingElemID;
158     p | requestingNodeSetID;
159     replyingElem.pup(p);
160   }
161 };
162
163
164 typedef ElemList<adjRequest> AdjRequestList;
165 typedef MSA1D<AdjRequestList, DefaultListEntry<AdjRequestList,true>,MSA_DEFAULT_ENTRIES_PER_PAGE> MSA1DREQLIST;
166
167 typedef ElemList<adjReply> AdjReplyList;
168 typedef MSA1D<AdjReplyList, DefaultListEntry<AdjReplyList,true>, MSA_DEFAULT_ENTRIES_PER_PAGE> MSA1DREPLYLIST;
169
170 /** Create Adaptivity Adjacencies for elemType; dimension inferred. */
171 void CreateAdaptAdjacencies(int meshid, int elemType);
172
173 // Access functions
174
175 // 2D gets
176 /** Look up elemID in elemType array, access edgeFaceID-th adaptAdj. */
177 adaptAdj *GetAdaptAdj(int meshid, adaptAdj elem, int edgeFaceID);
178 adaptAdj *GetAdaptAdj(FEM_Mesh *meshPtr, adaptAdj elem, int edgeFaceID);
179 /** Look up elemID in elemType array, calculate edgeFaceID from
180     vertexList (with GetEdgeFace below), and access edgeFaceID-th
181     adaptAdj with GetAdaptAdj above. */
182 adaptAdj *GetAdaptAdj(int meshid, adaptAdj elem, int *vertexList);
183
184 // 3D gets
185 adaptAdj *GetAdaptAdjOnEdge(int meshid, adaptAdj elem, int edgeID, int *size);
186 adaptAdj *GetAdaptAdjOnEdge(FEM_Mesh *meshPtr, adaptAdj elem, int edgeID, int *size);
187 adaptAdj *GetAdaptAdjOnFace(int meshid, adaptAdj elem, int faceID);
188 adaptAdj *GetAdaptAdjOnFace(FEM_Mesh *meshPtr, adaptAdj elem, int faceID);
189
190 /** Look up elemID in elemType array and determine the set of vertices
191     associated with the edge or face represented by edgeFaceID. */
192 void GetVertices(int meshid, adaptAdj elem, int edgeFaceID, int *vertexList);
193 /** Look up elemID in elemType array and determine the edge or face ID
194     specified by the set of vertices in vertexList. */
195 int GetEdgeFace(int meshid, adaptAdj elem, int *vertexList);
196
197 // Update functions
198 /** 2D or 3D (faces): Look up elemID in elemType array and set the adjacency on
199     edgeFaceID to nbr. */
200 void SetAdaptAdj(int meshid, adaptAdj elem, int edgeFaceID, adaptAdj nbr);
201
202 /** 3D: Look up elemID in elemType array and add nbr to the adjacency on
203     edgeID. */
204 void AddToAdaptAdj(int meshid, adaptAdj elem, int edgeID, adaptAdj nbr);
205 /** 3D: Look up elemID in elemType array and remove nbr from the adjacency on
206     edgeID. */
207 void RemoveFromAdaptAdj(int meshid, adaptAdj elem, int edgeID, adaptAdj nbr);
208
209 /** Substitute an old neighbor with a new neighbor, assumes 2D or 3D-face neighbor */
210 void ReplaceAdaptAdj(int meshID, adaptAdj elem, adaptAdj originalNbr, 
211                      adaptAdj newNbr);
212 void ReplaceAdaptAdj(FEM_Mesh *meshPtr, adaptAdj elem, adaptAdj originalNbr,
213                      adaptAdj newNbr);
214 /** 3D edge neighbors: Substitution operation needs to know edgeID to reduce search space. */
215 void ReplaceAdaptAdjOnEdge(int meshID, adaptAdj elem, adaptAdj originalNbr, 
216                      adaptAdj newNbr);
217 void ReplaceAdaptAdjOnEdge(FEM_Mesh *meshPtr, adaptAdj elem, adaptAdj originalNbr,
218                      adaptAdj newNbr);
219
220 /**given the dimensions and nodes per element guess whether the element 
221  is a triangle, quad, tet or hex. At the moment these are the 4 shapes
222  that are handled */
223 void guessElementShape(int dim,int nodesPerElem,int *numAdjElems,
224                        int *nodeSetSize,
225                        int nodeSetMap[MAX_ADJELEMS][MAX_NODESET_SIZE]);
226 void dumpAdaptAdjacencies(adaptAdj *adaptAdjacencies,int numElems,
227                           int numAdjElems,int myRank);
228 void getAndDumpAdaptAdjacencies(int meshid, int numElems, int elemType, int myRank);
229 void fillLocalAdaptAdjacencies(int numNodes, FEM_Node *node, adjNode *adaptAdjTable,
230                                adaptAdj *adaptAdjacencies, int nodeSetSize,
231                                int numAdjElems, int myRank, int elemType);
232 void makeAdjacencyRequests(int numNodes, FEM_Node *node, adjNode *adaptAdjTable,
233                            MSA1DREQLIST *requestTable, int nodeSetSize, int myRank,
234                            int elemType);
235 void replyAdjacencyRequests(MSA1DREQLIST *requestTable, MSA1DREPLYLIST *replyTable,
236                             FEM_Node *node, adjNode *adaptAdjTable, 
237                             adaptAdj *adaptAdjacencies, int nodeSetSize,
238                             int numAdjElems, int myRank, int elemType);