This version works for the explicit FEA code provided by Professor Paolino.
[charm.git] / src / libs / ck-libs / ParFUM-Tops / ParFUM_TOPS.C
1 /**
2
3    @file
4    @brief Implementation of ParFUM-TOPS layer, except for Iterators
5    
6    @author Isaac Dooley
7
8    @todo add code to generate ghost layers
9    @todo Support multiple models
10
11    @note FEM_DATA+0 holds the elemAttr or nodeAtt data
12    @note FEM_DATA+1 holds the id 
13    @note FEM_CONN holds the element connectivity
14    @note FEM_COORD holds the nodal coordinates
15
16 */
17
18 #include "ParFUM_TOPS.h"
19 #include "ParFUM.decl.h"
20 #include "ParFUM_internals.h"
21
22 int elem_attr_size, node_attr_size;
23
24 TopModel* topModel_Create_Init(int elem_attr_sz, int node_attr_sz){
25   CkAssert(elem_attr_sz > 0);
26   CkAssert(node_attr_sz > 0);
27   elem_attr_size = elem_attr_sz;
28   node_attr_size = node_attr_sz;
29
30   // This only uses a single mesh, so better not create multiple ones of these
31   int which_mesh=FEM_Mesh_default_write();
32   FEM_Mesh *mesh = FEM_Mesh_lookup(which_mesh,"TopModel::TopModel()");
33
34   char* temp_array = new char[16]; // just some junk array to use below
35   
36   // Allocate element connectivity
37   FEM_Mesh_data(which_mesh,FEM_ELEM+0,FEM_CONN,temp_array, 0, 0, FEM_INDEX_0, 4);
38   // Allocate element attributes
39   FEM_Mesh_data(which_mesh,FEM_ELEM+0,FEM_DATA+0,temp_array, 0, 0, FEM_BYTE, elem_attr_size);
40   // Allocate element Id array
41   FEM_Mesh_data(which_mesh,FEM_ELEM+0,FEM_DATA+1,temp_array, 0, 0, FEM_INT, 1);
42
43   // Allocate node coords
44   FEM_Mesh_data(which_mesh,FEM_NODE,FEM_COORD,temp_array, 0, 0, FEM_DOUBLE, 3);
45   // Allocate node attributes
46   FEM_Mesh_data(which_mesh,FEM_NODE+0,FEM_DATA+0,temp_array, 0, 0, FEM_BYTE, node_attr_size);
47   // Allocate node Id array
48   FEM_Mesh_data(which_mesh,FEM_NODE+0,FEM_DATA+1,temp_array, 0, 0, FEM_INT, 1);
49
50   delete[] temp_array;
51
52   // Don't Allocate the Global Number attribute for the elements and nodes  
53   // It will be automatically created upon calls to void FEM_Entity::setGlobalno(int r,int g) {
54
55   FEM_Mesh_allocate_valid_attr(which_mesh, FEM_NODE);
56   FEM_Mesh_allocate_valid_attr(which_mesh, FEM_ELEM+0);
57   
58
59   // Setup ghost layers, assume Tet4
60   // const int tet4vertices[4] = {0,1,2,3};
61   //  FEM_Add_ghost_layer(1,1);
62   // FEM_Add_ghost_elem(0,4,tet4vertices);
63
64   return mesh;
65 }
66
67 TopModel* topModel_Create_Driver(int elem_attr_sz, int node_attr_sz){
68   // This only uses a single mesh, so don't create multiple TopModels of these
69   CkAssert(elem_attr_sz > 0);
70   CkAssert(node_attr_sz > 0);
71   elem_attr_size = elem_attr_sz;
72   node_attr_size = node_attr_sz;
73   int which_mesh=FEM_Mesh_default_read();
74   FEM_Mesh *mesh = FEM_Mesh_lookup(which_mesh,"TopModel::TopModel()");
75   return mesh;
76 }
77
78 void topModel_Destroy(TopModel* m){
79 }
80
81
82 TopNode topModel_InsertNode(TopModel* m, double x, double y, double z){
83   int newNode = FEM_add_node_local(m,false,false,false);
84   m->node.set_coord(newNode,x,y,z);
85   return newNode;
86 }
87
88
89 /** Set id of a node 
90 @todo Make this work with ghosts
91 */
92 void topNode_SetId(TopModel* m, TopNode n, TopID id){
93   CkAssert(n>=0);
94   FEM_DataAttribute * at = (FEM_DataAttribute*) m->node.lookup(FEM_DATA+1,"topNode_SetId");
95   at->getInt()(n,0)=id;
96 }
97
98 /** Insert an element */
99 TopElement topModel_InsertElem(TopModel*m, TopElemType type, TopNode* nodes){
100   CkAssert(type ==  TOP_ELEMENT_TET4);
101   int conn[4];
102   conn[0] = nodes[0];
103   conn[1] = nodes[1];
104   conn[2] = nodes[2];
105   conn[3] = nodes[3];
106   int newEl = FEM_add_element_local(m, conn, 4, 0, 0, 0);
107   return newEl;
108 }
109
110 /** Set id of an element 
111 @todo Make this work with ghosts
112 */
113 void topElement_SetId(TopModel* m, TopElement e, TopID id){
114   CkAssert(e>=0);
115   FEM_DataAttribute * at = (FEM_DataAttribute*) m->elem[0].lookup(FEM_DATA+1,"topElement_SetId");
116   at->getInt()(e,0)=id;
117 }
118
119
120
121 /** 
122         @brief Set attribute of a node 
123         
124         The attribute passed in must be a contiguous data structure with size equal to the value node_attr_sz passed into topModel_Create_Driver() and topModel_Create_Init() 
125
126         The supplied attribute will be copied into the ParFUM attribute array "FEM_DATA+0". Then ParFUM will own this data. The function topNode_GetAttrib() will return a pointer to the copy owned by ParFUM. If a single material parameter attribute is used for multiple nodes, each node will get a separate copy of the array. Any subsequent modifications to the data will only be reflected at a single node. 
127
128         The user is responsible for deallocating parameter d passed into this function.
129
130 */
131 void topNode_SetAttrib(TopModel* m, TopNode n, void* d){
132   CkAssert(n>=0);
133   FEM_DataAttribute * at = (FEM_DataAttribute*) m->node.lookup(FEM_DATA+0,"topNode_SetAttrib");
134   AllocTable2d<unsigned char> &dataTable  = at->getChar();
135   unsigned char *data = dataTable.getData();
136   memcpy(data + n*node_attr_size, d, node_attr_size);
137 }
138
139 /** @brief Set attribute of an element 
140 See topNode_SetAttrib() for description
141 */
142 void topElement_SetAttrib(TopModel* m, TopElement e, void* d){
143   CkAssert(e>=0);
144   FEM_DataAttribute * at = (FEM_DataAttribute*) m->elem[0].lookup(FEM_DATA+0,"topElem_SetAttrib");
145   AllocTable2d<unsigned char> &dataTable  = at->getChar();
146   unsigned char *data = dataTable.getData();
147   memcpy(data + e*elem_attr_size, d, elem_attr_size);
148 }
149
150
151 /** @brief Get elem attribute 
152 See topNode_SetAttrib() for description
153 */
154 void* topElement_GetAttrib(TopModel* m, TopElement e){
155   if(! m->elem[0].is_valid_any_idx(e))
156         return NULL;
157
158   FEM_DataAttribute * at = (FEM_DataAttribute*) m->elem[0].lookup(FEM_DATA+0,"topElem_GetAttrib");
159   AllocTable2d<unsigned char> &dataTable  = at->getChar();
160   unsigned char *data = dataTable.getData();
161   return (data + e*elem_attr_size);
162 }
163
164 /** @brief Get nodal attribute 
165 See topNode_SetAttrib() for description
166 */
167 void* topNode_GetAttrib(TopModel* m, TopNode n){
168   if(! m->node.is_valid_any_idx(n))
169         return NULL;
170
171   FEM_DataAttribute * at = (FEM_DataAttribute*) m->node.lookup(FEM_DATA+0,"topNode_GetAttrib");
172   AllocTable2d<unsigned char> &dataTable  = at->getChar();
173   unsigned char *data = dataTable.getData();
174   return (data + n*node_attr_size);
175 }
176
177
178
179 /** 
180         Get node via id 
181         @todo Re-implement this function with some hash to make it fast. 
182         @note In the explicit FEA example, this is just used during initialization, so speed is not too important.
183         @todo Does not work with ghosts yet.
184 */
185 TopNode topModel_GetNodeAtId(TopModel* m, TopID id){
186   // lookup node via global ID
187   FEM_DataAttribute * at = (FEM_DataAttribute*) m->node.lookup(FEM_DATA+1,"topModel_GetNodeAtId");
188   for(int i=0;i<at->getInt().size();++i){
189         if(at->getInt()(i,0)==id){
190           return i;
191         }
192   }
193   return -1;
194 }
195
196 /** 
197         Get elem via id
198         @todo Implement this function
199         @note Is this function even supposed to exist?
200  */
201 TopElement topModel_GetElemAtId(TopModel*m,TopID id){
202   CkPrintf("Error: topModel_GetElemAtId() not yet implemented");
203   CkExit();
204 }
205
206
207
208 TopNode topElement_GetNode(TopModel* m,TopElement e,int idx){
209   CkAssert(e>=0);
210   const AllocTable2d<int> &conn = m->elem[0].getConn();
211   CkAssert(idx>=0 && idx<conn.width());
212   CkAssert(idx<conn.size());
213
214   int node = conn(e,idx);
215
216   return conn(e,idx);
217 }
218
219 int topNode_GetId(TopModel* m, TopNode n){
220   CkAssert(n>=0);
221   FEM_DataAttribute * at = (FEM_DataAttribute*) m->node.lookup(FEM_DATA+1,"topNode_SetId");
222   return at->getInt()(n,0);
223 }
224
225
226 /** @todo handle ghost nodes as appropriate */
227 int topModel_GetNNodes(TopModel *model){
228   return model->node.count_valid();
229 }
230
231 /** @todo Fix to return the width of the conn array */
232 int topElement_GetNNodes(TopModel* model, TopElement elem){
233   return 4;
234 }
235
236 /** @todo make sure we are in a getting mesh */
237 void topNode_GetPosition(TopModel*model, TopNode node,double*x,double*y,double*z){
238   CkAssert(node>=0);
239   double coord[3]={7.01,7.01,7.01};
240
241   // lookup node via global ID
242   FEM_DataAttribute * at = (FEM_DataAttribute*) model->node.lookup(FEM_COORD,"topModel_GetNodeAtId");
243   CkAssert(at->getDouble().width()==3);
244   //  CkPrintf("node=%d, size=%d\n",node,at->getDouble().size() );
245   CkAssert(node<at->getDouble().size());
246   *x = at->getDouble()(node,0);
247   *y = at->getDouble()(node,1);
248   *z = at->getDouble()(node,2);
249 }
250
251 void topModel_Sync(TopModel*m){
252   MPI_Barrier(MPI_COMM_WORLD);
253
254  
255   //  CkPrintf("%d: %d local, %d ghost elements\n", FEM_My_partition(), m->elem[0].size(),m->elem[0].ghost->size() );
256   //  CkPrintf("%d: %d local, %d ghost valid elements\n", FEM_My_partition(), m->elem[0].count_valid(),m->elem[0].ghost->count_valid() );
257
258 }
259
260 void topModel_TestIterators(TopModel*m){
261   CkAssert(m->elem[0].ghost!=NULL);
262   CkAssert(m->node.ghost!=NULL);
263
264   int expected_elem_count = m->elem[0].count_valid() + m->elem[0].ghost->count_valid(); 
265   int iterated_elem_count = 0;
266
267   int expected_node_count = m->node.count_valid() + m->node.ghost->count_valid(); 
268   int iterated_node_count = 0;
269
270   int myId = FEM_My_partition();
271
272
273   TopNodeItr* itr = topModel_CreateNodeItr(m);
274   for(topNodeItr_Begin(itr);topNodeItr_IsValid(itr);topNodeItr_Next(itr)){
275         iterated_node_count++;
276         TopNode node = topNodeItr_GetCurr(itr);
277         void* na = topNode_GetAttrib(m,node);
278         CkAssert(na != NULL);
279   }
280     
281   TopElemItr* e_itr = topModel_CreateElemItr(m);
282   for(topElemItr_Begin(e_itr);topElemItr_IsValid(e_itr);topElemItr_Next(e_itr)){
283         iterated_elem_count++;
284         TopElement elem = topElemItr_GetCurr(e_itr);
285         void* ea = topElement_GetAttrib(m,elem);
286         CkAssert(ea != NULL);
287   }
288   
289   CkAssert(iterated_node_count == expected_node_count);
290   CkAssert(iterated_elem_count==expected_elem_count);
291   
292 }
293
294
295 #include "ParFUM_TOPS.def.h"