7b84e82cab7298987dbe99e33741cdc611eb4cb7
[charm.git] / src / libs / ck-libs / ParFUM / mesh_modify.C
1 /**
2  * \addtogroup ParFUM
3 */
4 /*@{*/
5
6
7 /** File: fem_mesh_modify.C
8  * Authors: Nilesh Choudhury, Isaac Dooley
9  * 
10  * This file contains a set of functions, which allow primitive operations upon meshes in parallel.
11  * See the assumptions listed in fem_mesh_modify.h before using these functions.
12  * 
13  * It also contains a shadow array class that is attached to a chunk and
14  * handles all remote communication between chunks (needed for adaptivity)
15  */
16
17 #include "ParFUM.h"
18 #include "ParFUM_internals.h"
19
20 CProxy_femMeshModify meshMod;
21
22
23 /** A wrapper to simplify the lookup to determine whether a node is shared
24  */
25 inline int is_shared(FEM_Mesh *m, int node){
26   return m->getfmMM()->getfmUtil()->isShared(node);
27 }
28
29
30 /* Some helper print functions that can be used for debugging 
31    helper functions that print the various node/element connectivities
32    calls the corresponding functions in the FEM_MUtil class to do the printing */
33 CDECL void FEM_Print_n2n(int mesh, int nodeid){
34   FEM_Mesh *m=FEM_Mesh_lookup(mesh,"FEM_Print_Mesh_Summary");
35   m->getfmMM()->getfmUtil()->FEM_Print_n2n(m, nodeid);
36 }
37 CDECL void FEM_Print_n2e(int mesh, int eid){
38   FEM_Mesh *m=FEM_Mesh_lookup(mesh,"FEM_Print_Mesh_Summary");
39   m->getfmMM()->getfmUtil()->FEM_Print_n2e(m, eid);
40 }
41 CDECL void FEM_Print_e2n(int mesh, int eid){
42   FEM_Mesh *m=FEM_Mesh_lookup(mesh,"FEM_Print_Mesh_Summary");
43   m->getfmMM()->getfmUtil()->FEM_Print_e2n(m, eid);
44 }
45 CDECL void FEM_Print_e2e(int mesh, int eid){
46   FEM_Mesh *m=FEM_Mesh_lookup(mesh,"FEM_Print_Mesh_Summary");
47   m->getfmMM()->getfmUtil()->FEM_Print_e2e(m, eid);
48 }
49
50 /** Prints the mesh summary, i.e. number of valid nodes/elements (local/ghost)
51  */
52 CDECL void FEM_Print_Mesh_Summary(int mesh){
53   CkPrintf("---------------FEM_Print_Mesh_Summary-------------\n");
54   FEM_Mesh *m=FEM_Mesh_lookup(mesh,"FEM_Print_Mesh_Summary");
55   // Print Node information
56   CkPrintf("     Nodes: %d/%d and ghost nodes: %d/%d\n", m->node.count_valid(), m->node.size(),m->node.getGhost()->count_valid(),m->node.getGhost()->size());
57   // Print Element information
58   CkPrintf("     Element Types: %d\n", m->elem.size());
59   for (int t=0;t<m->elem.size();t++) {// for each element type t
60     if (m->elem.has(t)) {
61       int numEl = m->elem[t].size();
62       int numElG = m->elem[t].getGhost()->size();
63       int numValidEl = m->elem[t].count_valid();
64       int numValidElG = m->elem[t].getGhost()->count_valid();
65       CkPrintf("     Element type %d contains %d/%d elements and %d/%d ghosts\n", t, numValidEl, numEl, numValidElG, numElG);
66     }
67   }
68   CkPrintf("\n");
69 }
70
71
72
73 //========================Basic Adaptivity operations=======================
74 /*These are the basic add/remove node/element operations */
75 //The following functions translate the calls from meshId to mesh pointer
76 CDECL int FEM_add_node(int mesh, int* adjacent_nodes, int num_adjacent_nodes, int *chunks, int numChunks, int forceShared){
77   return FEM_add_node(FEM_Mesh_lookup(mesh,"FEM_add_node"), adjacent_nodes, num_adjacent_nodes, chunks, numChunks, forceShared);
78 }
79 CDECL void FEM_remove_node(int mesh,int node){
80   return FEM_remove_node(FEM_Mesh_lookup(mesh,"FEM_remove_node"), node);
81 }
82 CDECL int FEM_add_element(int mesh, int* conn, int conn_size, int elem_type, int chunkNo){
83   return FEM_add_element(FEM_Mesh_lookup(mesh,"FEM_add_element"), conn, conn_size, elem_type, chunkNo);
84 }
85 CDECL int FEM_remove_element(int mesh, int element, int elem_type, int permanent){
86   return FEM_remove_element(FEM_Mesh_lookup(mesh,"FEM_remove_element"), element, elem_type, permanent);
87 }
88 CDECL int FEM_purge_element(int mesh, int element, int elem_type) {
89   return FEM_purge_element(FEM_Mesh_lookup(mesh,"FEM_remove_element"), element, elem_type);
90 }
91
92
93 /* Implementation of the add/remove/purge node/element functions 
94    (along with helpers)
95    Add node and helpers -- add_local/ add_remote/ add_shared, etc*/
96
97 /** Add a node to this mesh 'm', could be ghost/local decided by 'addGhost'
98     basically, this function finds an empty row in the nodelist
99     and clears its adjacencies
100     then, it creates a lock for this node, if this node is new (i.e. unused)
101     otherwise, the existing lock is reset
102 */
103 int FEM_add_node_local(FEM_Mesh *m, bool addGhost, bool doLocking, bool doAdjacencies){
104   int newNode;
105   if(addGhost){
106     // find a place for new node in tables, reusing old invalid nodes if possible
107     newNode = m->node.getGhost()->get_next_invalid(m,true,true); 
108         if(doAdjacencies){
109           m->n2e_removeAll(FEM_To_ghost_index(newNode));    // initialize element adjacencies
110           m->n2n_removeAll(FEM_To_ghost_index(newNode));    // initialize node adjacencies
111         }
112   }
113   else{
114     newNode = m->node.get_next_invalid(m,true,false);
115     if(doAdjacencies){
116           m->n2e_removeAll(newNode);    // initialize element adjacencies
117           m->n2n_removeAll(newNode);    // initialize node adjacencies
118         }
119     //add a lock
120     if(doLocking){
121           if(newNode >= m->getfmMM()->fmLockN.size()) {
122                 m->getfmMM()->fmLockN.push_back(FEM_lockN(newNode,m->getfmMM()));
123           }
124           else {
125                 m->getfmMM()->fmLockN[newNode].reset(newNode,m->getfmMM());
126           }
127         }
128   }
129 #ifdef DEBUG_2
130   CkPrintf("[%d]Adding Node %d\n",m->getfmMM()->getfmUtil()->getIdx(),(addGhost==0)?newNode:FEM_To_ghost_index(newNode));
131 #endif
132   return newNode;  // return a new index
133 }
134
135
136 /** 'adjacentNodes' specifies the set of nodes between which this new node
137     is to be added.
138     'forceshared' is used to convey extra information. Sometimes, we know that
139     the new node will not be shared, i.e. if the two nodes are shared, but the
140     two elements on two sides of this node belong to the same chunk, then this node
141     will not be shared. (in spite of the two nodes being shared by the same set of chunks
142     'numChunks' is the number of entries in 'chunks'.
143     'chunks' contains the chunk index where this node is to be added.
144     This is populated when the operation originates from edge_bisect and the 
145     operation knows which all chunks this node needs to be added to.
146     Returns the index of the new node that was added
147  */
148 int FEM_add_node(FEM_Mesh *m, int* adjacentNodes, int numAdjacentNodes, int*chunks, int numChunks, int forceShared){
149   // add local node
150   //should be used only when all the adjacent nodes are shared 
151   //but you know that the new node should not
152   //be shared, but should belong to only one of the chunks sharing the face/edge.
153   //usually it is safe to use -1, except in some weird cases.
154   int index = m->getfmMM()->getIdx();
155   if(numChunks==1 && chunks[0]!=index) {
156     //translate the indices.. the new node will be a ghost
157     int chunkNo = chunks[0];
158     addNodeMsg *am = new (numAdjacentNodes,numChunks,0) addNodeMsg;
159     am->chk = index;
160     am->nBetween = numAdjacentNodes;
161     for(int i=0; i<numAdjacentNodes; i++) {
162       am->between[i] = m->getfmMM()->getfmUtil()->exists_in_IDXL(m,adjacentNodes[i],chunkNo,0);
163       CkAssert(am->between[i]!=-1);
164     }
165     am->chunks[0] = chunkNo; //there should be only 1 chunk
166     am->numChunks = numChunks;
167     intMsg *imsg = new intMsg(-1);
168     m->getfmMM()->getfmUtil()->idxllock(m,chunkNo,0);
169     imsg = meshMod[chunkNo].addNodeRemote(am);
170     int newghost = FEM_add_node_local(m,1);
171     m->node.ghost->ghostRecv.addNode(newghost,chunkNo);
172     m->getfmMM()->getfmUtil()->idxlunlock(m,chunkNo,0);
173     //this is the ghostsend index on that chunk, translate it from ghostrecv idxl table
174     //return FEM_To_ghost_index(m->getfmMM()->getfmUtil()->lookup_in_IDXL(m,imsg->i,chunkNo,2));
175     //rather this is same as
176     delete imsg;
177     return FEM_To_ghost_index(newghost);
178   }
179   int newNode = FEM_add_node_local(m, 0);
180   int sharedCount = 0;
181   //interpolate the node data on the new node, based on the ones from which it is created
182   FEM_Interpolate *inp = m->getfmMM()->getfmInp();
183   FEM_Interpolate::NodalArgs nm;
184   nm.n = newNode;
185   for(int i=0; i<numAdjacentNodes; i++) {
186     nm.nodes[i] = adjacentNodes[i];
187   }
188   nm.frac = 0.5;
189   nm.addNode = true;
190   inp->FEM_InterpolateNodeOnEdge(nm);
191   //if any of the adjacent nodes are shared...
192   if(numChunks>1) {
193     // for each adjacent node, if the node is shared
194     for(int i=0;i<numAdjacentNodes;i++){ //a newly added node is shared only if all the 
195       //nodes between which it is added are shared
196       if(is_shared(m, adjacentNodes[i]))
197         {
198           sharedCount++;
199           // lookup adjacent_nodes[i] in IDXL, to find all remote chunks which share this node
200           // call_shared_node_remote() on all chunks for which the shared node exists
201           // we must make sure that we only call the remote entry method once for each remote chunk
202         }
203     }
204     //this is the entry in the IDXL.
205     //just a sanity check, not required really
206     if((sharedCount==numAdjacentNodes && numAdjacentNodes!=0) && forceShared!=-1 ) {
207       //m->getfmMM()->getfmUtil()->splitEntityAll(m, newNode, numAdjacentNodes, adjacentNodes, 0);
208       m->getfmMM()->getfmUtil()->splitEntitySharing(m, newNode, numAdjacentNodes, adjacentNodes, numChunks, chunks);
209     }
210   }
211   return newNode;
212 }
213
214 /** The function called by the entry method on the remote chunk
215     Adds a node on this chunk between the nodes specified in 'between'
216  */
217 void FEM_add_shared_node_remote(FEM_Mesh *m, int chk, int nBetween, int *between){
218   // create local node
219   int newnode = FEM_add_node_local(m, 0);
220   // must negotiate the common IDXL number for the new node, 
221   // and store it in appropriate IDXL tables
222   //note that these are the shared indices, local indices need to be calculated
223   m->getfmMM()->getfmUtil()->splitEntityRemote(m, chk, newnode, nBetween, between);
224 }
225
226
227
228 /*Remove node along with helpers -- remove_local/ remove_remote, etc */
229 /** This function invalidates a node entry and if there are any 
230     entries in the IDXL lists (ghost/local) then, these entries are cleared
231     both on the local as well as the remote chunk
232     and then the node is invalidated and any lock it might be holding is reset
233     The pre-condition which is asserted here that the node should not
234     have any adjacency entries, n2e or n2n.
235 */
236 void FEM_remove_node_local(FEM_Mesh *m, int node) {
237   // if node is local:
238   int numAdjNodes, numAdjElts;
239   int *adjNodes, *adjElts;
240   m->n2n_getAll(node, adjNodes, numAdjNodes);
241   m->n2e_getAll(node, adjElts, numAdjElts);
242   // mark node as deleted/invalid
243   if(FEM_Is_ghost_index(node)){
244     //CkAssert((numAdjNodes==0) && (numAdjElts==0));
245     //otherwise this ghost node is connected to some element in another chunk, 
246     //which the chunk that just informed us doesn't know abt
247     //look up the ghostrecv idxl list & clean up all instances
248     const IDXL_Rec *irec = m->node.ghost->ghostRecv.getRec(FEM_To_ghost_index(node));
249     //this ghost node might be sent from a chunk without any elems
250     int size = 0;
251     if(irec) size = irec->getShared();
252     int *chunks1, *inds1;
253     if(size>0) {
254       chunks1 = new int[size];
255       inds1 = new int[size];
256     }
257     for(int i=0; i<size; i++) {
258       chunks1[i] = irec->getChk(i);
259       inds1[i] = irec->getIdx(i);
260     }
261     int index = m->getfmMM()->getIdx();
262     for(int i=0; i<size; i++) {
263       int chk = chunks1[i];
264       int sharedIdx = inds1[i];
265       m->node.ghost->ghostRecv.removeNode(FEM_To_ghost_index(node),chk);
266       meshMod[chk].removeIDXLRemote(index,sharedIdx,1);
267     }
268     if(size>0) {
269       delete [] chunks1;
270       delete [] inds1;
271     }
272     m->node.ghost->set_invalid(FEM_To_ghost_index(node),true);
273   }
274   else {
275     //look it up on the idxl list and delete any instances in it
276     const IDXL_Rec *irec = m->node.ghostSend.getRec(node);
277     int size = 0;
278     if(irec) size = irec->getShared();
279     if(size > 0) {
280       int *chknos = (int*)malloc(size*sizeof(int));
281       int *sharedIndices = (int*)malloc(size*sizeof(int));
282       for(int i=0; i<size; i++) {
283         chknos[i] = irec->getChk(i);
284         sharedIndices[i] = irec->getIdx(i);
285       }
286       for(int chkno=0; chkno<size; chkno++) {
287         int remoteChunk = chknos[chkno];
288         int sharedIdx = sharedIndices[chkno];
289         //go to that remote chunk & delete this idxl list and this ghost node
290         meshMod[remoteChunk].removeGhostNode(m->getfmMM()->getIdx(), sharedIdx);
291         m->node.ghostSend.removeNode(node, remoteChunk);
292       }
293       free(chknos);
294       free(sharedIndices);
295     }
296     /*if(!((numAdjNodes==0) && (numAdjElts==0))) {
297       CkPrintf("Error: Node %d cannot be removed, it is connected to :\n",node);
298       m->getfmMM()->fmUtil->FEM_Print_n2e(m,node);
299       m->getfmMM()->fmUtil->FEM_Print_n2n(m,node);
300       }*/
301     // we shouldn't be removing a node away that is connected to something
302     CkAssert((numAdjNodes==0) && (numAdjElts==0));
303     //should be done for the locked node only, i.e. the node on the smallest chunk no
304     m->node.set_invalid(node,true);
305     m->getfmMM()->fmLockN[node].reset(node,m->getfmMM());
306   }
307   if(numAdjNodes != 0) delete[] adjNodes;
308   if(numAdjElts != 0) delete[] adjElts;
309 #ifdef DEBUG_2
310   CkPrintf("[%d]Removing Node %d\n",m->getfmMM()->getfmUtil()->getIdx(),node);
311 #endif
312   return;
313 }
314
315 /** remove a local or shared node, but NOT a ghost node
316     Should probably be able to handle ghosts someday, but here
317     we interpret it as deleting the entry in the ghost list.
318 */
319 void FEM_remove_node(FEM_Mesh *m, int node){
320   CkAssert(node >= 0);
321   //someone might actually want to remove a ghost node... 
322   //when there is a ghost edge with both end nodes shared
323   //we will have to notice such a situation and call it on the remotechunk
324   if(FEM_Is_ghost_index(node)) {
325     //this is interpreted as a chunk trying to delete this ghost node, 
326     //this is just to get rid of its version of the node entity
327     CkAssert(m->node.ghost->is_valid(FEM_To_ghost_index(node))); // make sure the node is still there
328     FEM_remove_node_local(m,node);
329   }
330   else {
331     CkAssert(m->node.is_valid(node)); // make sure the node is still there
332     // if node is shared:
333     if(is_shared(m, node)){
334       // verify it is not adjacent to any elements locally
335       int numAdjNodes, numAdjElts;
336       int *adjNodes, *adjElts;
337       m->n2n_getAll(node, adjNodes, numAdjNodes);
338       m->n2e_getAll(node, adjElts, numAdjElts);
339       // we shouldn't be removing a node away that is connected to anything
340       CkAssert((numAdjNodes==0) && (numAdjElts==0)); 
341       // verify it is not adjacent to any elements on any of the associated chunks
342       // delete it on remote chunks(shared and ghost), update IDXL tables
343       m->getfmMM()->getfmUtil()->removeNodeAll(m, node);
344       // mark node as deleted/invalid locally
345       //FEM_remove_node_local(m,node);
346       //m->node.set_invalid(node,true);
347       if(numAdjNodes!=0) delete[] adjNodes;
348       if(numAdjElts!=0) delete[] adjElts;
349     }
350     else {
351       FEM_remove_node_local(m,node);
352     }
353   }
354 }
355
356
357
358 /*Add an element, along with helpers -- add_remote, add_local, update_e2e etc */
359 /** A helper function for FEM_add_element_local below
360     Will only work with the same element type as the one given, may crash otherwise
361 */
362 void update_new_element_e2e(FEM_Mesh *m, int newEl, int elemType){
363   // this function most definitely will not yet work with mixed element types.
364   CkAssert(elemType==0); 
365   // Create tuple table
366   FEM_ElemAdj_Layer *g = m->getElemAdjLayer();
367   CkAssert(g->initialized);
368   const int nodesPerTuple = g->nodesPerTuple;
369   const int tuplesPerElem = g->elem[elemType].tuplesPerElem;
370   tupleTable table(nodesPerTuple);
371   FEM_Symmetries_t allSym;
372   // locate all elements adjacent to the nodes adjacent to the new 
373   // element, including ghosts, and the new element itself
374   const int nodesPerElem = m->elem[elemType].getNodesPer();
375   int *adjnodes = new int[nodesPerElem];
376   CkVec<int> elist;
377   m->e2n_getAll(newEl, adjnodes, elemType);
378   for(int i=0;i<nodesPerElem;i++){
379     int sz;
380     int *adjelements=0;
381     m->n2e_getAll(adjnodes[i], adjelements, sz);
382     for(int j=0;j<sz;j++){
383       int found=0;
384       // only insert if it is not already in the list
385       // we use a slow linear scan of the CkVec
386       for(int i=0;i<elist.length();i++) {
387         if(elist[i] == adjelements[j]) found=1;
388       }
389       if(!found){
390         elist.push_back(adjelements[j]);
391       }
392     }
393     if(sz!=0) delete[] adjelements;
394   }
395   delete[] adjnodes;
396   // Add all the potentially adjacent elements to the tuple table
397   for(int i=0;i<elist.length();i++){
398     int nextElem = elist[i];
399     int tuple[tupleTable::MAX_TUPLE];
400     int *conn;
401     if(FEM_Is_ghost_index(nextElem))
402       conn=((FEM_Elem*)m->elem[elemType].getGhost())->connFor(FEM_To_ghost_index(nextElem));
403     else
404       conn=m->elem[elemType].connFor(nextElem);
405     for (int u=0;u<tuplesPerElem;u++) {
406       for (int i=0;i<nodesPerTuple;i++) {
407         int eidx=g->elem[elemType].elem2tuple[i+u*g->nodesPerTuple];
408         if (eidx==-1)  //"not-there" node--
409           tuple[i]=-1; //Don't map via connectivity
410         else           //Ordinary node
411           tuple[i]=conn[eidx]; 
412       }
413       table.addTuple(tuple,new elemList(0,nextElem,elemType,allSym,u)); 
414     }
415   }
416   // extract true adjacencies from table and update all e2e tables for both newEl and the others
417   table.beginLookup();
418   // look through each elemList that is returned by the tuple table
419   elemList *l;
420   FEM_IndexAttribute *elemAdjTypesAttr = (FEM_IndexAttribute *)m->elem[elemType].lookup(FEM_ELEM_ELEM_ADJ_TYPES,"update_new_element_e2e");
421   FEM_IndexAttribute *elemAdjAttr = (FEM_IndexAttribute *)m->elem[elemType].lookup(FEM_ELEM_ELEM_ADJACENCY,"update_new_element_e2e");
422   FEM_IndexAttribute *elemAdjTypesAttrGhost = (FEM_IndexAttribute *)m->elem[elemType].getGhost()->lookup(FEM_ELEM_ELEM_ADJ_TYPES,"update_new_element_e2e");
423   FEM_IndexAttribute *elemAdjAttrGhost = (FEM_IndexAttribute *)m->elem[elemType].getGhost()->lookup(FEM_ELEM_ELEM_ADJACENCY,"update_new_element_e2e");
424   //allocate some data structures
425   AllocTable2d<int> &adjTable = elemAdjAttr->get();
426   int *adjs = adjTable.getData();
427   AllocTable2d<int> &adjTypesTable = elemAdjTypesAttr->get();
428   int *adjTypes = adjTypesTable.getData();
429   AllocTable2d<int> &adjTableGhost = elemAdjAttrGhost->get();
430   int *adjsGhost = adjTableGhost.getData();
431   AllocTable2d<int> &adjTypesTableGhost = elemAdjTypesAttrGhost->get();
432   int *adjTypesGhost = adjTypesTableGhost.getData();
433   while (NULL!=(l=table.lookupNext())) {
434     if (l->next==NULL) {} // One-entry list: must be a symmetry
435     else { /* Several elements in list: normal case */
436       // for each a,b pair of adjacent edges
437       for (const elemList *a=l;a!=NULL;a=a->next){
438         for (const elemList *b=a->next;b!=NULL;b=b->next){
439           // if a and b are different elements
440           if((a->localNo != b->localNo) || (a->type != b->type)){
441             //  CkPrintf("%d:%d:%d adj to %d:%d:%d\n", a->type, a->localNo, a->tupleNo, b->type, b->localNo, b->tupleNo);
442             // Put b in a's adjacency list
443             if(FEM_Is_ghost_index(a->localNo)){
444               const int j = FEM_To_ghost_index(a->localNo)*tuplesPerElem + a->tupleNo;
445               adjsGhost[j] = b->localNo;
446               adjTypesGhost[j] = b->type;
447             }
448             else{
449               const int j= a->localNo*tuplesPerElem + a->tupleNo;
450               adjs[j] = b->localNo;
451               adjTypes[j] = b->type;
452             }
453             // Put a in b's adjacency list
454             if(FEM_Is_ghost_index(b->localNo)){
455               const int j = FEM_To_ghost_index(b->localNo)*tuplesPerElem + b->tupleNo;
456               adjsGhost[j] = a->localNo;
457               adjTypesGhost[j] = a->type;
458             }
459             else{
460               const int j= b->localNo*tuplesPerElem + b->tupleNo;
461               adjs[j] = a->localNo;
462               adjTypes[j] = a->type;
463             }
464           }
465         }
466       }
467     }
468   }
469 }
470
471 /** A helper function that adds the local element, and updates adjacencies
472  */
473 int FEM_add_element_local(FEM_Mesh *m, int *conn, int connSize, int elemType, bool addGhost, bool create_adjacencies){
474   // lengthen element attributes
475   int newEl;
476   if(addGhost){
477     // find a place in the array for the new el
478     newEl = m->elem[elemType].getGhost()->get_next_invalid(m,false,true); 
479     // update element's conn, i.e. e2n table
480     ((FEM_Elem*)m->elem[elemType].getGhost())->connIs(newEl,conn);
481     // get the signed ghost value
482     newEl = FEM_From_ghost_index(newEl);
483   }
484   else{
485     newEl = m->elem[elemType].get_next_invalid(m,false,false);
486     /*#ifdef FEM_ELEMSORDERED
487     //correct the orientation, which might have been lost if it moves across chunks
488     if(m->getfmMM()->fmAdaptAlgs->getSignedArea(conn[0],conn[1],conn[2])>0) {
489       int tmp = conn[1];
490       conn[1] = conn[2];
491       conn[2] = tmp;
492     }
493     #endif*/
494     // update element's conn, i.e. e2n table
495     m->elem[elemType].connIs(newEl,conn); 
496   }
497
498   if(create_adjacencies){
499         // add to n2e and n2n table
500         for(int i=0;i<connSize;i++){
501           m->n2e_add(conn[i],newEl);
502           for(int j=i+1;j<connSize;j++){
503                 if(! m->n2n_exists(conn[i],conn[j]))
504                   m->n2n_add(conn[i],conn[j]);
505                 if(! m->n2n_exists(conn[j],conn[i]))
506                   m->n2n_add(conn[j],conn[i]);
507           }
508         }
509         // update e2e table -- too complicated, so it gets is own function
510         m->e2e_removeAll(newEl);
511         update_new_element_e2e(m,newEl,elemType);
512         int *adjes = new int[connSize];
513         m->e2e_getAll(newEl, adjes, 0);
514         CkAssert(!((adjes[0]==adjes[1] && adjes[0]!=-1) || (adjes[1]==adjes[2] && adjes[1]!=-1) || (adjes[2]==adjes[0] && adjes[2]!=-1)));
515 #ifdef DEBUG_2
516         CkPrintf("[%d]Adding element %d(%d,%d,%d)\n",m->getfmMM()->getfmUtil()->getIdx(),newEl,conn[0],conn[1],conn[2]);
517 #endif
518 #ifdef FEM_ELEMSORDERED
519         if(!FEM_Is_ghost_index(newEl)) {
520           m->getfmMM()->fmAdaptAlgs->ensureQuality(conn[0],conn[1],conn[2]);
521         }
522 #endif
523         delete[] adjes;
524   }
525
526   return newEl;
527 }
528
529 /** This addds an element with the given 'conn' connecitivity and 'connSize' number of nodes in e2n
530     to the mesh 'm'. 'chunkNo' denotes the chunk the new element should be added to.
531     If it is -1, then this function decides which chunk this new element should belong to,
532     based on its node connectivity.
533     It returns the index of the newly created element
534 */
535 int FEM_add_element(FEM_Mesh *m, int* conn, int connSize, int elemType, int chunkNo){
536   //the first step is to figure out the type of each node (local/shared/ghost)
537   int newEl = -1;
538   int index = m->getfmMM()->getIdx();
539   int buildGhosts = 0; //this determines if we need to modify ghosts locally or remotely
540                        //if all nodes are local, there is no need to update ghosts anywhere
541   int sharedcount=0;
542   int ghostcount=0;
543   int localcount=0;
544   CkAssert(conn[0]!=conn[1] &&conn[1]!=conn[2] &&conn[2]!=conn[0]);
545   int *nodetype = (int *)malloc(connSize *sizeof(int)); //0 -- local, 1 -- shared, 2--ghost
546   for(int i=0;i<connSize;i++){
547     //make sure that every connectivity is valid
548     CkAssert(conn[i]!=-1); 
549     if(is_shared(m,conn[i])) {
550       nodetype[i] = 1;
551       sharedcount++;
552     }
553     else if(FEM_Is_ghost_index(conn[i])) {
554       nodetype[i] = 2;
555       ghostcount++;
556     }
557     else {
558       nodetype[i] = 0;
559     }
560   }
561   //if it is not shared or ghost, it should be local
562   localcount = connSize - (sharedcount + ghostcount);
563
564   if(sharedcount==0 && ghostcount==0){
565     // add a local elem with all local nodes
566     newEl = FEM_add_element_local(m,conn,connSize,elemType,0);
567     //no modifications required for ghostsend or ghostrecv of nodes or elements
568   }
569   else if(ghostcount==0 && sharedcount > 0 && localcount > 0){
570     // a local elem with not ALL shared nodes or ALL local nodes
571     newEl = FEM_add_element_local(m,conn,connSize,elemType,0);
572     buildGhosts = 1;
573   }
574   else if(ghostcount==0 && sharedcount > 0 && localcount == 0){
575     // a local elem with ALL shared nodes
576     //it is interesting to note that such a situation can only occur between two chunks,
577     //in any number of dimensions.
578     //So, the solution is to figure out, which chunk it belongs to...
579     if(!(chunkNo==index || chunkNo==-1)) {
580       //if chunkNo points to a remotechunk, then this element should be added to that chunk
581       addElemMsg *am = new (connSize, 0, 0) addElemMsg;
582       int chk = index;
583       am->chk = chk;
584       am->elemtype = elemType;
585       am->connSize = connSize;
586       am->numGhostIndex = 0;
587       for(int i=0; i<connSize; i++) {
588         //translate the node connectivity local to this chunk to shared idxl entries
589         CkAssert(nodetype[i] == 1);
590         am->conn[i] = m->getfmMM()->getfmUtil()->exists_in_IDXL(m,conn[i],chunkNo,0);
591         CkAssert(am->conn[i] >= 0);
592       }
593       intMsg* imsg = meshMod[chunkNo].addElementRemote(am);
594       int shidx = imsg->i;
595       delete imsg;
596       const IDXL_List ilist = m->elem[elemType].ghost->ghostRecv.addList(chunkNo);
597       newEl = ilist[shidx];
598       newEl = FEM_To_ghost_index(newEl);
599       free(nodetype);
600       return newEl;
601     }
602     //otherwise, it will be added on this chunk
603     newEl = FEM_add_element_local(m,conn,connSize,elemType,0);
604     buildGhosts = 1;
605   }
606   else if(ghostcount > 0 /*&& localcount == 0 && sharedcount > 0*/) { 
607     //I will assume that this is a correct call.. the caller never gives junk ghost elements
608     //it is a remote elem with some shared nodes
609     //this is the part when it eats a ghost element
610     if((chunkNo!=-1) && (chunkNo==index)) { 
611       //this is because the chunk doing this operation is supposed to eat into someone else.
612       //change all ghost nodes to shared nodes
613       //this also involves going to allchunks that it was local/shared 
614       //on and make it shared and add this chunk to the list of chunks it is shared to.
615       for(int i=0; i<connSize; i++) {
616         if(nodetype[i]==2) {
617           //build up the list of chunks it is shared/local to
618           int numchunks;
619           IDXL_Share **chunks1;
620           //make sure that the ghostRecv table is completely populated...
621           //if this is a ghost->local transition, then it will be local on only one
622           //chunk, i.e. on the chunk where it was local before this eat operation started.
623           m->getfmMM()->getfmUtil()->getChunkNos(0,conn[i],&numchunks,&chunks1);
624           //add a new node with the same attributes as this ghost node,
625           //do not remove the ghost node yet
626           int newN = m->getfmMM()->getfmUtil()->Replace_node_local(m, conn[i], -1);
627           //need to know if the newly added node will be local or shared
628           //if shared, add index to the shared list of this node on all the chunks
629           //if local, remove the oldIdx on the other chunk, if the oldIdx
630           //will be local only on this chunk after this operation
631           //(happens during a local->ghost transition)
632           for(int j=0; j<numchunks; j++) {
633             int chk = chunks1[j]->chk;
634             if(chk==index) continue;
635             //coarser lock
636             m->getfmMM()->getfmUtil()->idxllock(m,chk,0);
637             //this new node, might be a shared node or a local node... 
638             m->node.shared.addNode(newN,chk);
639             //find out what chk calls this node (from the ghostrecv idxl list)
640             int idx = m->getfmMM()->getfmUtil()->exists_in_IDXL(m, conn[i], chk, 2);
641             m->node.ghost->ghostRecv.removeNode(FEM_To_ghost_index(conn[i]), chk);
642             //when doing this add, make all idxl lists consistent, and return the list 
643             //of n2e connections for this node, which are local to this chunk, 
644             //and also a list of the connectivity
645             meshMod[chk].addToSharedList(index, idx); 
646             //when adding a new node, always verify if this new node is a corner on
647             //the other chunk, if so add it to the list of fixed nodes (corners)
648             int idx1 = m->getfmMM()->getfmUtil()->exists_in_IDXL(m, newN, chk, 0);
649             boolMsg* bmsg1 = meshMod[chk].isFixedNodeRemote(index,idx1);
650             bool isfixed = bmsg1->b;
651             delete bmsg1;
652             if(isfixed) m->getfmMM()->fmfixedNodes.push_back(newN);
653             //coarser lock
654             m->getfmMM()->getfmUtil()->idxlunlock(m, chk, 0);
655           }
656           nodetype[i] = 1;
657           sharedcount++;
658           ghostcount--;
659           //remove the ghost node
660           FEM_remove_node_local(m,conn[i]);
661           conn[i] = newN;
662           //lock the newly formed node, if it needs to be
663           //FEM_Modify_LockUpdate(m,conn[i]);
664           for(int j=0; j<numchunks; j++) {
665             delete chunks1[j];
666           }
667           if(numchunks != 0) free(chunks1);
668         }
669       }
670       newEl = FEM_add_element_local(m,conn,connSize,elemType,0);
671       buildGhosts = 1;
672     }
673     else {
674       int numSharedChunks = 0;
675       int remoteChunk = -1;
676       if(chunkNo==-1) {
677         CkAssert(false); //shouldn't be here
678         //figure out which chunk it is local to
679         //among all chunks who share some of the nodes or from whom this chunk receives ghost nodes
680         //if there is any chunk which owns a ghost node which is not shared, then that is a local node 
681         //to that chunk and that chunk owns that element. However, only that chunk knows abt it.
682         //So, just go to the owner of every ghost node and figure out who all share that node.
683         //Build up this table of nodes owned by which all chunks.
684         //The chunk that is in the table corresponding to all nodes wins the element.
685         CkVec<int> **allShared;
686         CkVec<int> *allChunks = NULL;
687         int **sharedConn; 
688         m->getfmMM()->getfmUtil()->buildChunkToNodeTable(nodetype, sharedcount, ghostcount, localcount, conn, connSize, &allShared, &numSharedChunks, &allChunks, &sharedConn);   
689         //we are looking for a chunk which does not have a ghost node
690         for(int i=0; i<numSharedChunks; i++) {
691           remoteChunk = i;
692           for(int j=0; j<connSize; j++) {
693             if(sharedConn[i][j] == -1 || sharedConn[i][j] == 2) {
694               remoteChunk = -1;
695               break; //this chunk has a ghost node
696             }
697             //if(sharedConn[i][j] == 0) {
698             //  break; //this is a local node, hence it is the remotechunk
699             //}
700           }
701           if(remoteChunk == i) break;
702           else remoteChunk = -1;
703         }
704         if(remoteChunk==-1) {
705           //every chunk has a ghost node.
706           if(chunkNo != -1) {
707             //upgrade the ghost node to a shared node on chunkNo
708             remoteChunk = chunkNo;
709             for(int k=0; k<numSharedChunks; k++) {
710               if(chunkNo == (*allChunks)[k]) {
711                 for(int l=0; l<connSize; l++) {
712                   if(sharedConn[k][l]==-1 || sharedConn[k][l]==2) {
713                     //FIXME: upgrade this ghost node
714                   }
715                 }
716               }
717             }
718           }
719           else {
720             CkPrintf("[%d]ERROR: Can not derive where (%d,%d,%d) belongs to\n",index,conn[0],conn[1],conn[2]);
721             CkAssert(false);
722           }
723         }
724         remoteChunk = (*allChunks)[remoteChunk];
725         //convert all connections to the shared IDXL indices. We should also tell which are ghost indices
726         CkPrintf("[%d]Error: I derived it should go to chunk %d, which is not %d\n",index,remoteChunk,chunkNo);
727         for(int k=0; k<numSharedChunks; k++) {
728           free(sharedConn[k]);
729         }
730         if(numSharedChunks!=0) free(sharedConn);
731         for(int k=0; k<connSize; k++) {
732           delete allShared[k];
733         }
734         free(allShared);
735         allChunks->free();
736         if(allChunks!=NULL) free(allChunks);
737       }
738       else {
739         remoteChunk=chunkNo;
740       }
741       int numGhostNodes = 0;
742       for(int i=0; i<connSize; i++) {
743         if(nodetype[i] == 2) { //a ghost node
744           numGhostNodes++;
745         }
746       }
747       CkAssert(numGhostNodes > 0);
748       addElemMsg *am = new (connSize, numGhostNodes, 0) addElemMsg;
749       int chk = index;
750       am->chk = chk;
751       am->elemtype = elemType;
752       am->connSize = connSize;
753       am->numGhostIndex = numGhostNodes;
754       int j = 0;
755       for(int i=0; i<connSize; i++) {
756         if(nodetype[i] == 1) {
757           am->conn[i] = m->getfmMM()->getfmUtil()->exists_in_IDXL(m,conn[i],remoteChunk,0);
758           CkAssert(am->conn[i] >= 0);
759         }
760         else if(nodetype[i] == 2) {
761           am->conn[i] = m->getfmMM()->getfmUtil()->exists_in_IDXL(m,conn[i],remoteChunk,2);
762           CkAssert(am->conn[i] >= 0);
763           am->ghostIndices[j] = i;
764           j++;
765         }
766       }
767       intMsg* imsg = meshMod[remoteChunk].addElementRemote(am);
768       int shidx = imsg->i;
769       delete imsg;
770       const IDXL_List ilist = m->elem[elemType].ghost->ghostRecv.addList(remoteChunk);
771       newEl = ilist[shidx];
772       newEl = FEM_To_ghost_index(newEl);
773     }
774   }
775   else if(ghostcount > 0 && localcount == 0 && sharedcount == 0) { 
776     // it is a remote elem with no shared nodes
777     //I guess such a situation will never occur
778     //figure out which chunk it is local to -- this would be really difficult to do in such a case.
779     //so, for now, we do not allow a chunk to add an element 
780     //for which it does not share even a single node.
781     // almost the code in the preceeding else case will work for this, but its better not to allow this
782   }
783   else if(ghostcount > 0 && localcount > 0 && sharedcount > 0){
784     // it is a flip operation
785     //actually this can be generalized as an operation which moves the boundary,
786     //to make the ghost node shared
787     //if one uses FEM_elem_acquire, then this condition gets distributed 
788     //across other conditions already done.
789     //   promote ghosts to shared on others, requesting new ghosts
790     //   grow local element and attribute tables if needed
791     //   add to local elem[elemType] table, and update IDXL if needed
792     //   update remote adjacencies
793     //   update local adjacencies
794   }
795   else if(ghostcount > 0 && localcount > 0 && sharedcount == 0) { 
796     //this is an impossible case
797   }
798   //start building ghosts if required
799   if(buildGhosts==1) {
800     //   make this element ghost on all others, updating all IDXL's
801     //   also in same remote entry method, update adjacencies on all others
802     //   grow local element and attribute tables if needed
803     //   add to local elem[elemType] table, and update IDXL if needed
804     //   update local adjacencies
805     //   return the new element id
806     //build a mapping of all shared chunks to all nodes in this element    
807     CkVec<int> **allShared;
808     int numSharedChunks = 0;
809     CkVec<int> *allChunks = NULL;
810     int **sharedConn; 
811     m->getfmMM()->getfmUtil()->buildChunkToNodeTable(nodetype, sharedcount, ghostcount, localcount, conn, connSize, &allShared, &numSharedChunks, &allChunks, &sharedConn);   
812     //add all the local nodes in this element to the ghost list, if they did not exist already
813     for(int i=0; i<numSharedChunks; i++) {
814       int chk = (*allChunks)[i];
815       if(chk == index) continue; //it is this chunk
816       //it is a new element so it could not have existed as a ghost on that chunk. Just add it.
817       m->getfmMM()->getfmUtil()->idxllock(m, chk, 0);
818       m->elem[elemType].ghostSend.addNode(newEl,chk);
819       //ghost nodes should be added only if they were not already present as ghosts on that chunk.
820       int *sharedIndices = (int *)malloc((connSize)*sizeof(int));
821       int *typeOfIndex = (int *)malloc((connSize)*sizeof(int));
822       for(int j=0; j<connSize; j++) {
823         int sharedNode = m->getfmMM()->getfmUtil()->exists_in_IDXL(m,conn[j],chk,0);
824         if(sharedNode == -1) {
825           //node 'j' is a ghost on chunk 'i'
826           int sharedGhost = m->getfmMM()->getfmUtil()->exists_in_IDXL(m,conn[j],chk,1);
827           if( sharedGhost == -1) {
828             //this might be a new ghost, figure out if any of the chunks sharing
829             //this node has created this as a ghost on 'chk'
830             const IDXL_Rec *irec = m->node.shared.getRec(conn[j]);
831             if(irec) {
832               int noShared = irec->getShared();
833               for(int sharedck=0; sharedck<noShared; sharedck++) {
834                 int ckshared = irec->getChk(sharedck);
835                 int idxshared = irec->getIdx(sharedck);
836                 if(ckshared == chk) continue;
837                 CkAssert(chk!=index && chk!=ckshared && ckshared!=index);
838                 intMsg* imsg = meshMod[ckshared].getIdxGhostSend(index,idxshared,chk);
839                 int idxghostsend = imsg->i;
840                 delete imsg;
841                 if(idxghostsend != -1) {
842                   m->node.ghostSend.addNode(conn[j],chk);
843                   meshMod[chk].updateIdxlList(index,idxghostsend,ckshared);
844                   sharedGhost = m->getfmMM()->getfmUtil()->exists_in_IDXL(m,conn[j],chk,1);
845                   CkAssert(sharedGhost != -1);
846                   break; //found a chunk that sends it out, update my tables
847                 }
848                 //Chunk 'ckshared' does not send this to Chunk 'chk' as ghost
849               }
850             }
851             //else it is a new ghost
852           }
853           if(sharedGhost == -1) {
854             sharedIndices[j] = conn[j];
855             typeOfIndex[j] = -1;
856           }
857           else {
858             //it is a shared ghost
859             sharedIndices[j] = sharedGhost;
860             typeOfIndex[j] = 1;
861           }
862         }
863         else {
864           //it is a shared node
865           sharedIndices[j] = sharedNode;
866           typeOfIndex[j] = 0;
867         }
868       }
869       addGhostElemMsg *fm = new (connSize, connSize, 0)addGhostElemMsg;
870       fm->chk = index;
871       fm->elemType = elemType;
872       for(int j=0; j<connSize; j++) {
873         fm->indices[j] = sharedIndices[j];
874         fm->typeOfIndex[j] = typeOfIndex[j];
875         if(typeOfIndex[j]==-1) {
876           m->node.ghostSend.addNode(sharedIndices[j],chk);
877         }
878       }
879       fm->connSize = connSize;
880       meshMod[chk].addGhostElem(fm); //newEl, m->fmMM->idx, elemType;
881       m->getfmMM()->getfmUtil()->idxlunlock(m, chk, 0);
882       free(sharedIndices);
883       free(typeOfIndex);
884     }
885     for(int k=0; k<numSharedChunks; k++) {
886       free(sharedConn[k]);
887     }
888     if(numSharedChunks!=0) free(sharedConn);
889     for(int k=0; k<connSize; k++) {
890       delete allShared[k];
891     }
892     free(allShared);
893     allChunks->free();
894     if(allChunks!=NULL) free(allChunks);
895   }
896   free(nodetype);
897   return newEl;
898 }
899
900
901
902 /*Remove an element with helpers --remove_local, etc */
903 /** Remove a local element from any of the adjacency tables that use it,
904     namely the n2e, e2e and e2n tables of the neighbors. 
905     e2n of the neighbors does not change
906     all tables for this element get deleted anyway when it is marked invalid, which is 
907     now done in the purge step
908 */
909 void FEM_remove_element_local(FEM_Mesh *m, int element, int etype){
910   // replace this element with -1 in adjacent nodes' adjacencies
911   // nodesPerEl be the number of nodes that can be adjacent to this element
912   const int nodesPerEl = m->elem[etype].getConn().width(); 
913   int *adjnodes = new int[nodesPerEl];
914   m->e2n_getAll(element, adjnodes, etype);
915 #ifdef DEBUG_2
916   CkPrintf("[%d]Deleting element %d(%d,%d,%d)\n",m->getfmMM()->getfmUtil()->getIdx(),element,adjnodes[0],adjnodes[1],adjnodes[2]);
917 #endif
918   for(int i=0;i<nodesPerEl;i++) {
919     //if(adjnodes[i] != -1) //if an element is local, then an adjacent node should not be -1
920     m->n2e_remove(adjnodes[i],element);
921   }
922   // replace this element with -1 in adjacent elements' adjacencies
923   const int numAdjElts = nodesPerEl;  
924   int *adjelts = new int[numAdjElts]; 
925   m->e2e_getAll(element, adjelts, etype);
926   for(int i=0;i<numAdjElts;i++){
927     m->e2e_replace(adjelts[i],element,-1);
928   }
929   // We must now remove any n2n adjacencies which existed because of the 
930   // element that is now being removed. This is done by removing all 
931   // n2n adjacencies for the nodes of this element, and recreating those
932   // that still exist by using the neighboring elements.
933   //FIXME: I think, if we just consider on which all faces (for a 2D, a face is an edge),
934   //we have valid neighboring elements, then none of the edges (n2n conn) on that
935   //face should go off, this would be more efficient.
936   for(int i=0;i<nodesPerEl;i++) {
937     for(int j=i+1;j<nodesPerEl;j++){
938       m->n2n_remove(adjnodes[i],adjnodes[j]);
939       m->n2n_remove(adjnodes[j],adjnodes[i]);
940     }
941   }
942   for(int i=0;i<numAdjElts;i++){ // for each neighboring element
943     if(adjelts[i] != -1){
944       int *adjnodes2 = new int[nodesPerEl];
945       m->e2n_getAll(adjelts[i], adjnodes2, etype);
946       // for each j,k pair of nodes adjacent to the neighboring element
947       for(int j=0;j<nodesPerEl;j++){ 
948         for(int k=j+1;k<nodesPerEl;k++){   
949           if(!m->n2n_exists(adjnodes2[j],adjnodes2[k]))
950             m->n2n_add(adjnodes2[j],adjnodes2[k]);
951           if(!m->n2n_exists(adjnodes2[k],adjnodes2[j]))
952             m->n2n_add(adjnodes2[k],adjnodes2[j]);
953         }
954       }
955       delete[] adjnodes2;
956     }
957   }
958   //done in purge now
959   /*if(FEM_Is_ghost_index(element)){
960     m->elem[etype].getGhost()->set_invalid(FEM_To_ghost_index(element),false);
961   }
962   else {
963     m->elem[etype].set_invalid(element,false);
964     }*/
965   //cleanup
966   delete[] adjelts;
967   delete[] adjnodes;
968   return;
969 }
970
971 /** Can be called on local or ghost elements
972     removes the connectivity of the element on this chunk as well as all chunks on which it is a ghost
973     'permanent' is -1 if the element is to be deleted, e.g. in a coarsen operation, one element
974     is deleted permanently.. for a refine however the existing elements are eaten by the new ones
975     'permanent' is the id of the chunk which wants to be the owner (i.e. wants to eat it)
976     this is done just to save deleting some things when it would be immediately added again, 
977     when the element is deleted on one and added on another chunk 
978     The return value is the chunk on which this element was local
979     Removing an element is a two step process, FEM_remove_element and FEM_purge_element
980     This step will get rid of all adjacencies info, but the idxl information
981     will not be deleted, just ghostsend on this junk and receive info on the remote chunk.
982     The idxl info will be deleted in Purge, that way shared 
983     chunks could still refer to it, needed for solution transfer
984 */
985 int FEM_remove_element(FEM_Mesh *m, int elementid, int elemtype, int permanent){
986   //CkAssert(elementid != -1);
987   if(elementid == -1) return -1;
988   int index = m->getfmMM()->getfmUtil()->getIdx();
989   if(FEM_Is_ghost_index(elementid)){
990     //an element can come as a ghost from only one chunk, so just convert args and call it on that chunk
991     int ghostid = FEM_To_ghost_index(elementid);
992     const IDXL_Rec *irec = m->elem[elemtype].ghost->ghostRecv.getRec(ghostid);
993     int size = irec->getShared();
994     CkAssert(size == 1);
995     int remoteChunk = irec->getChk(0);
996     int sharedIdx = irec->getIdx(0);
997     //the message that should go to the remote chunk, where this element lives
998     removeElemMsg *rm = new removeElemMsg;
999     rm->chk = index;
1000     rm->elementid = sharedIdx;
1001     rm->elemtype = elemtype;
1002     rm->permanent = permanent;
1003     meshMod[remoteChunk].removeElementRemote(rm);
1004     // remove local ghost element, now done in purge
1005     return remoteChunk;
1006   }
1007   else {
1008     //this is a local element on this chunk
1009     //if it is a ghost element on some other chunks, then for all ghost nodes in this element,
1010     //we should individually go to each ghost and figure out if it has any connectivity in
1011     //the ghost layer after the removal. If not, then it will be removed from their ghost layers.
1012     //This is something similar to an add_element, where we do a similar analaysis on each 
1013     //involved chunk neighbor and add ghost nodes.
1014     //get a list of chunks for which this element is a ghost
1015     //this has to be done as the idxl lists change in this loop, so we
1016     //if we call getRec multiple times, we will get a new irec every time.
1017     const IDXL_Rec *irec = m->elem[elemtype].ghostSend.getRec(elementid);
1018     if(irec){
1019       int numSharedChunks = irec->getShared();
1020       int connSize = m->elem[elemtype].getConn().width();
1021       if(numSharedChunks>0) {
1022         int *chknos = (int*)malloc(numSharedChunks*sizeof(int));
1023         int *inds = (int*)malloc(numSharedChunks*sizeof(int));
1024         for(int i=0; i<numSharedChunks; i++) {
1025           chknos[i] = irec->getChk(i);
1026           inds[i] = irec->getIdx(i);
1027         }
1028         int *newghost, *ghostidx, *losingThisNode, *nodes, *willBeGhost;
1029         //will someone acquire this element soon?
1030         if(permanent>=0) { //this element will be reassigned or change connectivity
1031           newghost = new int[connSize];
1032           ghostidx = new int[connSize]; //need to create new ghosts
1033           losingThisNode = new int[connSize];
1034           nodes = new int[connSize];
1035           willBeGhost = new int[connSize]; 
1036           //willBeGhost is needed because not necessarily a node that is lost
1037           //will be a ghost, think of what happens when an isolated element is deleted
1038           int *elems;
1039           int numElems;
1040           m->e2n_getAll(elementid,nodes,elemtype);
1041           for(int i=0; i<connSize; i++) {
1042             newghost[i] = -1;
1043             ghostidx[i] = -1;
1044             //Is this chunk losing this node?
1045             losingThisNode[i] = 1;
1046             //if this node is connected to any element which is local and is not being deleted now
1047             //then this chunk will not lose this node now
1048             m->n2e_getAll(nodes[i], elems, numElems);
1049             for(int k=0; k<numElems; k++) {
1050               if((!FEM_Is_ghost_index(elems[k])) && (elems[k]!=elementid)) {
1051                 losingThisNode[i] = 0;
1052                 break;
1053               }
1054             }
1055             if(losingThisNode[i]) {
1056               willBeGhost[i] = 1; 
1057               //current allotment, it might be modified, see below
1058             }
1059             free(elems);
1060           }
1061           //verify if it is losing all nodes, this means that this was an isolated element
1062           //and this element does not have to be a ghost
1063           if(losingThisNode[0]==1 && losingThisNode[1]==1 && losingThisNode[2]==1) {
1064             //only if it is connected to any other node of this chunk, it will be a ghost
1065             for(int i=0; i<connSize; i++) {
1066               int *ndnbrs;
1067               int numndnbrs;
1068               m->n2n_getAll(nodes[i], ndnbrs, numndnbrs);
1069               willBeGhost[i]=0;
1070               //if this node is connected to any node that is local other than the nodes of this
1071               //element, only then will this node be a ghost.
1072               for(int n1=0; n1<numndnbrs; n1++) {
1073                 if(ndnbrs[n1]>=0 && ndnbrs[n1]!=nodes[(i+1)%connSize] && ndnbrs[n1]!=nodes[(i+2)%connSize]) {
1074                   willBeGhost[i]=1;
1075                 }
1076               }
1077             }
1078           }
1079           //add a new ghost for every node that should be lost and a ghost added
1080           for(int i=0; i<connSize; i++) {
1081             if(losingThisNode[i]) {
1082               //lock it on the min chunk on which this node is local
1083               FEM_Modify_LockAll(m,nodes[i],false);
1084             }
1085             if(losingThisNode[i]==1 && willBeGhost[i]==1) {
1086               newghost[i] = FEM_add_node_local(m,1);
1087               ghostidx[i] = FEM_To_ghost_index(newghost[i]); //translate the index
1088             }
1089           }
1090         }
1091         //the above code just figured out what to do with all the nodes of this element
1092         //and corrected the locks as well
1093         /*//a correction to deal with physical corners
1094         //if losing a node which is a corner, then upgrade that ghostnode on the chunk 
1095         //to a shared node before performing this entire operation
1096         //I can only think of this happening for 1 chunk
1097         for(int i=0; i<connSize; i++) {
1098           if(losingThisNode[i]==1 && m->node.shared.getRec(nodes[i]==NULL)) {
1099           }
1100           }*/
1101         int deleteNodeLater = -1; //if a local node becomes ghost, then one
1102         //has to keep track of that and delete it only after the remote chunk
1103         //has copied the attributes, this interger remembers which node is to be
1104         //deleted later
1105         //Note that the chunk to which it loses is always 'permanent'
1106         for(int i=0; i<numSharedChunks; i++) {
1107           bool lockedRemoteIDXL = false;
1108           int chk = chknos[i];
1109           int sharedIdx = inds[i];
1110           int numGhostNodes = 0;
1111           int *ghostIndices = (int*)malloc(connSize*sizeof(int));
1112           int numGhostRN = 0;
1113           CkVec<int> ghostRNIndices;
1114           int numGhostRE = 0;
1115           CkVec<int> ghostREIndices;
1116           int numSharedNodes = 0;
1117           int *sharedIndices = (int*)malloc(connSize*sizeof(int));
1118           //purge will do this now
1119           //m->elem[elemtype].ghostSend.removeNode(elementid, chk);
1120           if(permanent>=0) {
1121             //get the list of n2e for all nodes of this element. 
1122             //If any node has only this element in its list.
1123             //it no longer should be a ghost on chk
1124             //the next step is to figure out all the ghost elements 
1125             //& nodes that it no longer needs to have
1126             CkVec<int> testelems;
1127             const IDXL_List ll = m->elem[elemtype].ghostSend.addList(chk);
1128             int size = ll.size();
1129             const IDXL_List ln = m->node.ghostSend.addList(chk);
1130             int sizeN = ln.size();
1131             const IDXL_List lre = m->elem[elemtype].ghost->ghostRecv.addList(chk);
1132             int lresize = lre.size();
1133             const IDXL_List lrn = m->node.ghost->ghostRecv.addList(chk);
1134             int lrnsize = lrn.size();
1135             //iterate over the nodes connected to this element
1136             for(int j=0; j<connSize; j++) {
1137               int *elems;
1138               int numElems;
1139               //we are tyring to figure out if this local node, which might be a ghost
1140               //on some chunk be deleted on that chunk
1141               m->n2e_getAll(nodes[j], elems, numElems);
1142               //do not delete ghost nodes on the eating chunk, and be sure 
1143               //that the node is indeed a ghost on the chunk being considered
1144               if((chk != permanent) && (m->getfmMM()->fmUtil->exists_in_IDXL(m,nodes[j],chk,1)!=-1)) {
1145                 //'chk' is not the chunk that is eating this element
1146                 //if any of these elems is a ghost on chk then do not delete this ghost node
1147                 bool shouldBeDeleted = true;
1148                 for(int k=0; k<numElems; k++) {
1149                   if(elems[k]==elementid) continue;
1150                   if(elems[k]>0) {
1151                     if(m->getfmMM()->fmUtil->exists_in_IDXL(m,elems[k],chk,3)!=-1) {
1152                       shouldBeDeleted = false;
1153                       break;
1154                     }
1155                   }
1156                 }
1157                 //find out what the other shared chunks think abt it
1158                 //if every shared chunk believes it should be deleted, get rid of it
1159                 if(shouldBeDeleted) {
1160                   const IDXL_Rec *irecsh = m->node.shared.getRec(nodes[j]);
1161                   if(irecsh!=NULL) {
1162                     for(int k=0; k<irecsh->getShared(); k++) {
1163                       if(shouldBeDeleted) {
1164                         boolMsg* bmsg = meshMod[irecsh->getChk(k)].shouldLoseGhost(index,irecsh->getIdx(k),chk);
1165                         shouldBeDeleted = bmsg->b;
1166                         delete bmsg;
1167                       }
1168                     }
1169                   }
1170                 }
1171                 //since we have a remote call in between, to make the program thread-safe
1172                 //we will need to verify if the connectivity of this node changed
1173                 int *elems1, numElems1;
1174                 m->n2e_getAll(nodes[j], elems1, numElems1);
1175                 if(numElems1!=numElems) {
1176                   CkAssert(false);
1177                 }
1178                 else {
1179                   if(numElems1!=0) delete[] elems1;
1180                 }
1181
1182                 //add this to the list of ghost nodes to be deleted on the remote chunk
1183                 if(shouldBeDeleted) {
1184                   //convert this local index to a shared index
1185                   int shidx = m->getfmMM()->fmUtil->exists_in_IDXL(m,nodes[j],chk,1);
1186                   if(shidx!=-1) {
1187                     m->node.ghostSend.removeNode(nodes[j], chk);
1188                     ghostIndices[numGhostNodes] = shidx;
1189                     numGhostNodes++;
1190                   }
1191                 }
1192               }
1193               //try to find out which ghosts this chunk will be losing
1194               //if I am losing this node, then only will I lose some ghosts
1195               if(losingThisNode[j]) {
1196                 //was either a shared node, but now it should be a ghost node
1197                 //or it was a local node(corner), but now should be a ghost
1198                 for(int k=0; k<numElems; k++) {
1199                   int *nds = (int*)malloc(m->elem[elemtype].getConn().width()*sizeof(int));
1200                   if(FEM_Is_ghost_index(elems[k]) && m->getfmMM()->getfmUtil()->exists_in_IDXL(m,elems[k],chk,4,elemtype)!=-1) {
1201                     //this ghost element might be deleted, if it is not connected to any
1202                     //other node that is local or shared to this chunk
1203                     m->e2n_getAll(elems[k], nds, elemtype);
1204                     int geShouldBeDeleted = 1;
1205                     //look at the connectivity of this element
1206                     for(int l=0; l<connSize; l++) {
1207                       //if a node is not a ghost but is about to be lost anyway, 
1208                       //then it should be deleted
1209                       if(!FEM_Is_ghost_index(nds[l]) && (nodes[j]!=nds[l])) {
1210                         if(losingThisNode[(j+1)%connSize]==1 && nds[l]==nodes[(j+1)%connSize]) continue;
1211                         else if(losingThisNode[(j+2)%connSize]==1 && nds[l]==nodes[(j+2)%connSize]) continue;
1212                         geShouldBeDeleted = 0;
1213                       }
1214                     }
1215                     //the element will be deleted only if all the nodes are about to be lost
1216                     //here we delete this ghost element and clean the idxl lists as well
1217                     //and mark the element invalid in the fem_elem list
1218                     if(geShouldBeDeleted) {
1219                       int sge = m->getfmMM()->getfmUtil()->exists_in_IDXL(m,elems[k],chk,4,elemtype);
1220                       m->elem[elemtype].ghost->ghostRecv.removeNode(FEM_To_ghost_index(elems[k]), chk);
1221                       FEM_remove_element_local(m,elems[k],elemtype);
1222                       m->elem[elemtype].ghost->set_invalid(FEM_From_ghost_index(elems[k]),false);
1223                       ghostREIndices.push_back(sge);
1224                       numGhostRE++;
1225                       //find out if any nodes need to be removed
1226                       for(int l=0; l<connSize; l++) {
1227                         int *elts;
1228                         int numElts;
1229                         m->n2e_getAll(nds[l], elts, numElts);
1230                         //if this is no longer connected to this chunk through 
1231                         //any of its neighboring elements, it should no longer be a ghost
1232                         if(nds[l]<-1) { //it is a ghost
1233                           bool removeflag = true;
1234                           for(int lm=0; lm<numElts; lm++) {
1235                             //if it is connected to any element on this chunk, be it 
1236                             //local or ghost, other than the element that is about
1237                             //to be deleted, then it will still be a ghost
1238                             if(elts[lm]!=elems[k]) {
1239                               removeflag = false;
1240                             }
1241                           }
1242                           if(removeflag) {
1243                             //remove this ghost node on this chunk
1244                             int sgn = m->getfmMM()->getfmUtil()->exists_in_IDXL(m,nds[l],chk,2);
1245                             m->node.ghost->ghostRecv.removeNode(FEM_To_ghost_index(nds[l]), chk);
1246                             if(numElts==0) FEM_remove_node_local(m,nds[l]);
1247                             ghostRNIndices.push_back(sgn);
1248                             numGhostRN++;
1249                           }
1250                         }
1251                         if(numElts!=0) delete[] elts;
1252                       }
1253                     }
1254                   }
1255                   free(nds);
1256                 }
1257                 //if it is shared, tell that chunk, it no longer exists on me
1258                 int ssn = m->getfmMM()->getfmUtil()->exists_in_IDXL(m,nodes[j],chk,0);
1259                 //add a new coarse idxl lock, should lock a remote idxl only once
1260                 //even if it loses more than one nodes
1261                 if(!lockedRemoteIDXL) m->getfmMM()->getfmUtil()->idxllock(m,chk,0);
1262                 lockedRemoteIDXL = true;
1263                 bool losinglocal = false;
1264                 if(ssn==-1 && chk==permanent) {
1265                   //seems to me that if a local node becomes a ghost because of an element
1266                   //being eaten by a chunk, then it has to be a corner
1267                   //the other chunk must have this as a ghost
1268                   losinglocal = true;
1269                   ssn = -1000000000;
1270                   ssn += m->getfmMM()->getfmUtil()->exists_in_IDXL(m,nodes[j],chk,1);
1271                   if(willBeGhost[j]==1) {
1272                     m->node.ghost->ghostRecv.addNode(newghost[j],chk);
1273                   }
1274                   else { //happens when it was local but won't even be a ghost now
1275                     ssn -= 500000000; //to tell the other chunk, not to add it in ghostsend
1276                   }
1277                   //m->node.ghostSend.removeNode(nodes[j], chk);
1278                   sharedIndices[numSharedNodes] = ssn;
1279                   numSharedNodes++;
1280                 }
1281                 else if(ssn!=-1) {
1282                   if(willBeGhost[j]==1) {
1283                     m->node.ghost->ghostRecv.addNode(newghost[j],chk);
1284                   }
1285                   else { //happens when it was shared but won't even be a ghost now
1286                     ssn -= 500000000; //to tell the other chunk, not to add it in ghostsend
1287                   }
1288                   m->node.shared.removeNode(nodes[j], chk);
1289                   sharedIndices[numSharedNodes] = ssn;
1290                   numSharedNodes++;
1291                 }
1292                 if(i==numSharedChunks-1) { //this is the last time
1293                   //so, update this element's adjacencies now
1294                   //since this node is about to be lost, so update the neighboring 
1295                   //elements and nodes that this node has a connectivity with,
1296                   //about the new node that will replace it
1297                   //add only to the elems which still exist
1298                   //if this node will not be a ghost, then these neighboring nodes and elems
1299                   //to this node do not really need updating, but I believe it will not
1300                   //hurt, because the connectivity in this case will be -1
1301                   int *elems1, numElems1;
1302                   m->n2e_getAll(nodes[j], elems1, numElems1);
1303                   for(int k=0; k<numElems1; k++) {
1304                     bool flagElem = false;
1305                     if(FEM_Is_ghost_index(elems1[k])) {
1306                       if(m->elem[0].ghost->is_valid(FEM_From_ghost_index(elems1[k]))) {
1307                         flagElem = true;
1308                       }
1309                     }
1310                     else {
1311                       if(m->elem[0].is_valid(elems1[k])) flagElem = true;
1312                     }
1313                     if(flagElem) {
1314                       m->e2n_replace(elems1[k],nodes[j],ghostidx[j],elemtype);
1315                       m->n2e_remove(nodes[j],elems1[k]);
1316                       m->n2e_add(ghostidx[j],elems1[k]);
1317                       testelems.push_back(elems1[k]);
1318                     }
1319                   }
1320                   if(numElems1!=0) delete[] elems1;
1321                   int *n2ns;
1322                   int numn2ns;
1323                   m->n2n_getAll(nodes[j],n2ns,numn2ns);
1324                   for(int k=0; k<numn2ns; k++) {
1325                     m->n2n_replace(n2ns[k],nodes[j],ghostidx[j]);
1326                     m->n2n_remove(nodes[j],n2ns[k]);
1327                     m->n2n_add(ghostidx[j],n2ns[k]);
1328                   }
1329                   if(numn2ns!=0) delete[] n2ns;
1330                   //if it is not shared, then all shared entries have been deleted
1331                   const IDXL_Rec *irec1 = m->node.shared.getRec(nodes[j]);
1332                   if(!irec1) {
1333                     if(!losinglocal) {
1334                       FEM_remove_node_local(m,nodes[j]);
1335                     }
1336                     else {
1337                       deleteNodeLater = nodes[j];
1338                     }
1339                     //if losing a local node, then can remove it only after its attributes 
1340                     //have been copied, since this is the only copy.. 
1341                     //(as no other chunk has this node as local/shared)
1342                     //so we'll need to delete the ghostsend idxl entry and the node later
1343                   }
1344                 }
1345               }
1346               if(numElems!=0) delete[] elems;
1347             }
1348             //now that all ghost nodes to be removed have been decided, we add the elem,
1349             // and call the entry method test if all the elems added are valid
1350             //sanity testing code START
1351             for(int testelemscnt=0; testelemscnt <testelems.size(); testelemscnt++) {
1352               int el = testelems[testelemscnt];
1353               if(FEM_Is_ghost_index(el)) {
1354                 CkAssert(m->elem[0].ghost->is_valid(FEM_From_ghost_index(el))==1);
1355               }
1356               else {
1357                 CkAssert(m->elem[0].is_valid(el)==1);
1358               }
1359             }
1360             //sanity testing code END
1361           }
1362           removeGhostElemMsg *rm = new (numGhostNodes, numGhostRN, numGhostRE, numSharedNodes, 0) removeGhostElemMsg;
1363           rm->chk = index;
1364           rm->elemtype = elemtype;
1365           rm->elementid = sharedIdx;
1366           rm->numGhostIndex = numGhostNodes;
1367           for(int j=0; j<numGhostNodes; j++) {
1368             rm->ghostIndices[j] = ghostIndices[j];
1369           }
1370           rm->numGhostRNIndex = numGhostRN;
1371           for(int j=0; j<numGhostRN; j++) {
1372             rm->ghostRNIndices[j] = ghostRNIndices[j];
1373           }
1374           rm->numGhostREIndex = numGhostRE;
1375           for(int j=0; j<numGhostRE; j++) {
1376             rm->ghostREIndices[j] = ghostREIndices[j];
1377           }
1378           rm->numSharedIndex = numSharedNodes;
1379           for(int j=0; j<numSharedNodes; j++) {
1380             rm->sharedIndices[j] = sharedIndices[j];
1381           }
1382           meshMod[chk].removeGhostElem(rm);  //update the ghosts on all shared chunks
1383           //can delete the ghostSend and node now, if it was not deleted earlier
1384           if((i==numSharedChunks-1) && (deleteNodeLater!=-1)) {
1385             m->node.ghostSend.removeNode(deleteNodeLater, permanent);
1386             FEM_remove_node_local(m,deleteNodeLater);
1387             //finally if it is on the fmfixednodes list, remove it
1388             //an O(n) search, but fmfixedNodes is a small list and called rarely
1389             for(int cnt=0; cnt<m->getfmMM()->fmfixedNodes.size(); cnt++) {
1390               if(m->getfmMM()->fmfixedNodes[cnt]==deleteNodeLater) {
1391                 m->getfmMM()->fmfixedNodes.remove(cnt);
1392                 break;
1393               }
1394             }
1395           }
1396           if(lockedRemoteIDXL) m->getfmMM()->getfmUtil()->idxlunlock(m,chk,0);
1397           free(ghostIndices);
1398           free(sharedIndices);
1399         }
1400         free(chknos);
1401         free(inds);
1402         if(permanent>=0) {
1403           free(newghost);
1404           free(ghostidx);
1405           free(losingThisNode);
1406           free(nodes);
1407           free(willBeGhost);
1408         }
1409       }
1410     }
1411     // remove local element
1412     FEM_remove_element_local(m, elementid, elemtype);
1413   }
1414   return index;
1415 }
1416
1417
1418
1419 /** Purge an element
1420     works for both local and ghost elements. Gets rid of any IDXL connectivity that might
1421     be left behind and cleans it on all chunks. Then it invalidates the element entry
1422     in the element table
1423 */
1424 int FEM_purge_element(FEM_Mesh *m, int elementid, int elemtype) {
1425   if(elementid==-1) return 1;
1426   int index = m->getfmMM()->idx;
1427   if(FEM_Is_ghost_index(elementid)) {
1428     const IDXL_Rec *irec1 = m->elem[elemtype].ghost->ghostRecv.getRec(FEM_To_ghost_index(elementid));
1429     int remotechk = irec1->getChk(0);
1430     int sharedIdx = irec1->getIdx(0);
1431     meshMod[remotechk].purgeElement(index,sharedIdx);
1432   }
1433   else {
1434     //figure out all chunks where it is a ghost on
1435     const IDXL_Rec *irec = m->elem[elemtype].ghostSend.getRec(elementid);
1436     if(irec){
1437       //copy the (chunk,index) tuple
1438       int numSharedChunks = irec->getShared();
1439       int connSize = m->elem[elemtype].getConn().width();
1440       int *chknos, *inds;
1441       if(numSharedChunks>0) {
1442         chknos = (int*)malloc(numSharedChunks*sizeof(int));
1443         inds = (int*)malloc(numSharedChunks*sizeof(int));
1444         for(int i=0; i<numSharedChunks; i++) {
1445           chknos[i] = irec->getChk(i);
1446           inds[i] = irec->getIdx(i);
1447         }
1448       }
1449       //go and clean it up on each of these chunks
1450       for(int i=0; i<numSharedChunks; i++) {
1451         meshMod[chknos[i]].cleanupIDXL(index,inds[i]);
1452         m->elem[elemtype].ghostSend.removeNode(elementid, chknos[i]);
1453       }
1454       //free up allocated memory
1455       if(numSharedChunks>0) {
1456         free(chknos);
1457         free(inds);
1458       }
1459     }
1460   }
1461   // delete element by marking invalid
1462   if(!FEM_Is_ghost_index(elementid)){
1463     m->elem[elemtype].set_invalid(elementid,false);
1464   }
1465   return 1;
1466 }
1467 //========================Basic Adaptivity operations=======================
1468
1469
1470
1471
1472
1473
1474 //========================Locking operations=======================
1475 /* These are the locking operations on nodes or chunks, depending on the locking scheme */
1476 //these are based on chunk locking: DEPRECATED
1477 CDECL int FEM_Modify_Lock(int mesh, int* affectedNodes, int numAffectedNodes, int* affectedElts, int numAffectedElts, int elemtype){
1478   return FEM_Modify_Lock(FEM_Mesh_lookup(mesh,"FEM_Modify_Lock"), affectedNodes, numAffectedNodes, affectedElts, numAffectedElts, elemtype);
1479 }
1480 CDECL int FEM_Modify_Unlock(int mesh){
1481   return FEM_Modify_Unlock(FEM_Mesh_lookup(mesh,"FEM_Modify_Unlock"));
1482 }
1483 //these are based on node locking
1484 CDECL int FEM_Modify_LockN(int mesh, int nodeId, int readLock) {
1485   return FEM_Modify_LockN(FEM_Mesh_lookup(mesh,"FEM_Modify_LockN"),nodeId, readLock);
1486 }
1487 CDECL int FEM_Modify_UnlockN(int mesh, int nodeId, int readLock) {
1488   return FEM_Modify_UnlockN(FEM_Mesh_lookup(mesh,"FEM_Modify_UnlockN"),nodeId, readLock);
1489 }
1490
1491
1492 /** coarse lock (lock chunks)
1493     find out which all chunks these nodes belong to and lock them all
1494 */
1495 int FEM_Modify_Lock(FEM_Mesh *m, int* affectedNodes, int numAffectedNodes, int* affectedElts, int numAffectedElts, int elemtype) {
1496   return m->getfmMM()->getfmLock()->lock(numAffectedNodes, affectedNodes, numAffectedElts, affectedElts, elemtype);
1497 }
1498
1499 /** Unlock all chunks that was being held by a lock originating from here
1500  */
1501 int FEM_Modify_Unlock(FEM_Mesh *m) {
1502   return m->getfmMM()->getfmLock()->unlock();
1503 }
1504
1505 /** Lock a node
1506     get a 'read/write lock on 'nodeId'
1507     we find the chunk with the smallest id that owns/shares this node and lock this node
1508     on this chunk. This decreases the amount of communication and the number of 
1509     contention points for a single lock, as well as any chance of a deadlock
1510     since there is only a single resource to be acquired for a request
1511 */
1512 int FEM_Modify_LockN(FEM_Mesh *m, int nodeId, int readlock) {
1513   int ret = -1;
1514   if(is_shared(m,nodeId)) {
1515     //find the index of the least chunk it is shared on and try to lock it
1516     int index = m->getfmMM()->getIdx();
1517     int numchunks;
1518     IDXL_Share **chunks1;
1519     m->getfmMM()->getfmUtil()->getChunkNos(0,nodeId,&numchunks,&chunks1);
1520     int minChunk = MAX_CHUNK;
1521     for(int j=0; j<numchunks; j++) {
1522       int pchk = chunks1[j]->chk;
1523       if(pchk < minChunk) minChunk = pchk;
1524     }
1525     for(int j=0; j<numchunks; j++) {
1526       delete chunks1[j];
1527     }
1528     if(numchunks!=0) free(chunks1);
1529     CkAssert(minChunk!=MAX_CHUNK);
1530     if(minChunk==index) {
1531       if(readlock) {
1532         ret = m->getfmMM()->fmLockN[nodeId].rlock();
1533       } else {
1534         ret = m->getfmMM()->fmLockN[nodeId].wlock(index);
1535       }
1536       if(ret==1) {
1537         m->getfmMM()->fmLockN[nodeId].verifyLock();
1538       }
1539       return ret;
1540     }
1541     else {
1542       CkAssert(minChunk!=MAX_CHUNK);
1543       int sharedIdx = m->getfmMM()->getfmUtil()->exists_in_IDXL(m,nodeId,minChunk,0);
1544       if(sharedIdx < 0) return -1;
1545       intMsg* imsg;
1546       if(readlock) {
1547         imsg = meshMod[minChunk].lockRemoteNode(sharedIdx, index, 0, 1);
1548       } else {
1549         imsg = meshMod[minChunk].lockRemoteNode(sharedIdx, index, 0, 0);
1550       }
1551       ret = imsg->i;
1552       delete imsg;
1553       if(ret==1) {
1554         boolMsg* bmsg = meshMod[minChunk].verifyLock(index,sharedIdx,0);
1555         delete bmsg;
1556       }
1557       return ret;
1558     }
1559   }
1560   else if(FEM_Is_ghost_index(nodeId)) {
1561     //find the index of the least chunk that owns it & try to lock it
1562     int index = m->getfmMM()->getIdx();
1563     int numchunks;
1564     IDXL_Share **chunks1;
1565     m->getfmMM()->getfmUtil()->getChunkNos(0,nodeId,&numchunks,&chunks1);
1566     int minChunk = MAX_CHUNK;
1567     for(int j=0; j<numchunks; j++) {
1568       int pchk = chunks1[j]->chk;
1569       if(pchk == index) continue;
1570       if(pchk < minChunk) minChunk = pchk;
1571     }
1572     for(int j=0; j<numchunks; j++) {
1573       delete chunks1[j];
1574     }
1575     if(numchunks!=0) free(chunks1);
1576     if(minChunk==MAX_CHUNK) return -1;
1577     int sharedIdx = m->getfmMM()->getfmUtil()->exists_in_IDXL(m,nodeId,minChunk,2);
1578     if(sharedIdx < 0) return -1;
1579     intMsg* imsg;
1580     if(readlock) {
1581       imsg = meshMod[minChunk].lockRemoteNode(sharedIdx, index, 1, 1);
1582     } else {
1583       imsg = meshMod[minChunk].lockRemoteNode(sharedIdx, index, 1, 0);
1584     }
1585     ret = imsg->i;
1586     delete imsg;
1587     /*if(ret==1) {
1588       if(nodeId==-8 && index==4) {
1589         CkPrintf("Locking node %d on chunk %d\n",nodeId, minChunk);
1590       }
1591       boolMsg* bmsg = meshMod[minChunk].verifyLock(index, sharedIdx, 1);
1592       delete bmsg;
1593       }*/
1594     return ret;
1595   }
1596   else {
1597     if(readlock) {
1598       ret = m->getfmMM()->fmLockN[nodeId].rlock();
1599     } else {
1600       int index = m->getfmMM()->getIdx();
1601       ret = m->getfmMM()->fmLockN[nodeId].wlock(index);
1602     }
1603     /*if(ret==1) {
1604       m->getfmMM()->fmLockN[nodeId].verifyLock();
1605       }*/
1606     return ret;
1607   }
1608   return -1; //should not reach here
1609 }
1610
1611 /** Unlock this node
1612     find the chunk that holds the lock and unlock it on that chunk
1613     the type of the lock is also considered while unlocking a lock
1614 */
1615 int FEM_Modify_UnlockN(FEM_Mesh *m, int nodeId, int readlock) {
1616   if(is_shared(m,nodeId)) {
1617     //find the index of the least chunk it is shared on and try to unlock it
1618     int index = m->getfmMM()->getIdx();
1619     int numchunks;
1620     IDXL_Share **chunks1;
1621     m->getfmMM()->getfmUtil()->getChunkNos(0,nodeId,&numchunks,&chunks1);
1622     int minChunk = MAX_CHUNK;
1623     for(int j=0; j<numchunks; j++) {
1624       int pchk = chunks1[j]->chk;
1625       if(pchk < minChunk) minChunk = pchk;
1626     }
1627     for(int j=0; j<numchunks; j++) {
1628       delete chunks1[j];
1629     }
1630     if(numchunks!=0) free(chunks1);
1631     if(minChunk==index) {
1632       if(readlock) {
1633         return m->getfmMM()->fmLockN[nodeId].runlock();
1634       } else {
1635         return m->getfmMM()->fmLockN[nodeId].wunlock(index);
1636       }
1637     }
1638     else {
1639       CkAssert(minChunk!=MAX_CHUNK);
1640       int sharedIdx = m->getfmMM()->getfmUtil()->exists_in_IDXL(m,nodeId,minChunk,0);
1641       intMsg* imsg;
1642       if(readlock) {
1643         imsg = meshMod[minChunk].unlockRemoteNode(sharedIdx, index, 0, 1);
1644       } else {
1645         imsg = meshMod[minChunk].unlockRemoteNode(sharedIdx, index, 0, 0);
1646       }
1647       int ret = imsg->i;
1648       delete imsg;
1649       return ret;
1650     }
1651   }
1652   else if(FEM_Is_ghost_index(nodeId)) {
1653     //find the index of the least chunk that owns it & try to unlock it
1654     int index = m->getfmMM()->getIdx();
1655     int numchunks;
1656     IDXL_Share **chunks1;
1657     m->getfmMM()->getfmUtil()->getChunkNos(0,nodeId,&numchunks,&chunks1);
1658     int minChunk = MAX_CHUNK;
1659     for(int j=0; j<numchunks; j++) {
1660       int pchk = chunks1[j]->chk;
1661       if(pchk == index) continue;
1662       if(pchk < minChunk) minChunk = pchk;
1663     }
1664     for(int j=0; j<numchunks; j++) {
1665       delete chunks1[j];
1666     }
1667     if(numchunks!=0) free(chunks1);
1668     if(minChunk==MAX_CHUNK) return -1;
1669     int sharedIdx = m->getfmMM()->getfmUtil()->exists_in_IDXL(m,nodeId,minChunk,2);
1670     intMsg* imsg;
1671     if(readlock) {
1672       imsg = meshMod[minChunk].unlockRemoteNode(sharedIdx, index, 1, 1);
1673     } else {
1674       imsg = meshMod[minChunk].unlockRemoteNode(sharedIdx, index, 1, 0);
1675     }
1676     int ret = imsg->i;
1677     delete imsg;
1678     return ret;
1679   }
1680   else {
1681     if(readlock) {
1682       return m->getfmMM()->fmLockN[nodeId].runlock();
1683     } else {
1684       int index = m->getfmMM()->getIdx();
1685       return m->getfmMM()->fmLockN[nodeId].wunlock(index);
1686     }
1687   }
1688   return -1; //should not reach here
1689 }
1690
1691
1692
1693 /** When a chunk is losing a node, the lock might need to be reassigned
1694     lock it on the minchunk other than 'this chunk'
1695     should only be called when 'this chunk' is losing the node as local/shared
1696 */
1697 void FEM_Modify_LockAll(FEM_Mesh*m, int nodeId, bool lockall) {
1698   int index = m->getfmMM()->getIdx();
1699   int numchunks;
1700   const IDXL_Rec *irec = m->node.shared.getRec(nodeId);
1701   //if it is a corner node, it might not have been a shared node, so can't lock it on anything else
1702   if(irec) {
1703     numchunks = irec->getShared();
1704     if(!lockall) {
1705       int minchunk=MAX_CHUNK;
1706       int sharedIdx = -1;
1707       for(int i=0; i<numchunks; i++) {
1708         int pchk = irec->getChk(i); 
1709         if(pchk<minchunk) {
1710           minchunk = pchk;
1711           sharedIdx = irec->getIdx(i);
1712         }
1713       }
1714       CkAssert(minchunk!=MAX_CHUNK && sharedIdx!=-1);
1715       //lock it on this chunk, if not already locked
1716       intMsg* imsg = meshMod[minchunk].lockRemoteNode(sharedIdx, index, 0, 0);
1717       int done = imsg->i;
1718       delete imsg;
1719       //if done=-1, then it is already locked, otherwise we just locked it
1720     }
1721     else {
1722       for(int i=0; i<numchunks; i++) {
1723         int pchk = irec->getChk(i); 
1724         int sharedIdx = irec->getIdx(i);
1725         intMsg* imsg = meshMod[pchk].lockRemoteNode(sharedIdx, index, 0, 0);
1726         int done = imsg->i;
1727         delete imsg;
1728       }
1729       m->getfmMM()->fmLockN[nodeId].wlock(index);
1730     }
1731   }
1732   return;
1733 }
1734
1735 /**Deprecated: Update the lock on this node (by locking the newly formed node
1736    The locks will be cleared later
1737    must be a local node, lock it & then unlock it if needed
1738  */
1739 void FEM_Modify_LockUpdate(FEM_Mesh*m, int nodeId, bool lockall) {
1740   int index = m->getfmMM()->getIdx();
1741   int numchunks;
1742   const IDXL_Rec *irec = m->node.shared.getRec(nodeId);
1743   //if it is a corner, the new node might not be shared
1744   //nothing was locked, hence nothing needs to be unlocked too
1745   if(irec) {
1746     numchunks = irec->getShared();
1747     int minchunk=MAX_CHUNK;
1748     int minI=-1;
1749     for(int i=0; i<numchunks; i++) {
1750       int pchk = irec->getChk(i); 
1751       if(pchk<minchunk) {
1752         minchunk = pchk;
1753         minI=i;
1754       }
1755     }
1756     if(!lockall) {
1757       if(minchunk>index) {
1758         int prevminchunk=minchunk;
1759         minchunk=index;
1760         int sharedIdx = irec->getIdx(minI);
1761         CkAssert(prevminchunk!=MAX_CHUNK && sharedIdx!=-1);
1762         intMsg* imsg = meshMod[prevminchunk].unlockRemoteNode(sharedIdx, index, 0, 0);
1763         delete imsg;
1764       }
1765       else if(minchunk < index) {
1766         //unlock the previously acquired lock
1767         int sharedIdx = irec->getIdx(minI);
1768         intMsg* imsg = meshMod[minchunk].lockRemoteNode(sharedIdx, index, 0, 0);
1769         delete imsg;
1770         m->getfmMM()->fmLockN[nodeId].wunlock(index);
1771       }
1772     }
1773     else {
1774       if(minchunk>index) minchunk=index;
1775       if(minchunk!=index) {
1776         m->getfmMM()->fmLockN[nodeId].wunlock(index);
1777       }
1778       for(int i=0; i<numchunks; i++) {
1779         int pchk = irec->getChk(i);
1780         if(pchk!=minchunk) {
1781           int sharedIdx = irec->getIdx(i);
1782           intMsg* imsg = meshMod[pchk].unlockRemoteNode(sharedIdx, index, 0, 0);
1783           delete imsg;
1784         }
1785       }
1786     }
1787   }
1788   return;
1789 }
1790
1791 /** Deprecated: For the newly acquired node, correct the lock by removing superfluous locks
1792     Should always be called on the new node to correct the locks
1793     should make sure that it is called before unlocking, i.e. before the entire operation completes.
1794     gets rid of the extra safety lock & makes the locks consistent
1795 */
1796 void FEM_Modify_correctLockN(FEM_Mesh *m, int nodeId) {
1797   int index = m->getfmMM()->getIdx();
1798   int numchunks;
1799   IDXL_Share **chunks1;
1800   m->getfmMM()->getfmUtil()->getChunkNos(0,nodeId,&numchunks,&chunks1);
1801   int minChunk = MAX_CHUNK;
1802   int owner = -1;
1803   for(int j=0; j<numchunks; j++) {
1804     int pchk = chunks1[j]->chk;
1805     if(pchk < minChunk) minChunk = pchk;
1806   }
1807   if(is_shared(m,nodeId)) {
1808     for(int j=0; j<numchunks; j++) {
1809       int pchk = chunks1[j]->chk;
1810       if(pchk == index) owner = m->getfmMM()->fmLockN[nodeId].lockOwner();
1811       else {
1812         int sharedIdx = m->getfmMM()->getfmUtil()->exists_in_IDXL(m,nodeId,pchk,0);
1813         intMsg* imsg = meshMod[pchk].hasLockRemoteNode(sharedIdx, index, 0);
1814         owner = imsg->i;
1815         delete imsg;
1816       }
1817       if(owner != -1) { //this node is locked
1818         if(pchk == minChunk) {
1819           //the lock is good
1820           break;
1821         }
1822         else {
1823           //unlock the node on pchk & lock it on minChunk
1824           int locknodes = nodeId;
1825           int gotlocks = 1;
1826           int done = -1;
1827           if(pchk==index) m->getfmMM()->fmLockN[nodeId].wunlock(index);
1828           else {
1829             int sharedIdx = m->getfmMM()->getfmUtil()->exists_in_IDXL(m,nodeId,pchk,0);
1830             intMsg* imsg = meshMod[pchk].unlockRemoteNode(sharedIdx, index, 0, 0);
1831             delete imsg;
1832           }
1833           if(minChunk==index) done = m->getfmMM()->fmLockN[nodeId].wlock(index);
1834           else {
1835             CkAssert(minChunk!=MAX_CHUNK);
1836             int sharedIdx = m->getfmMM()->getfmUtil()->exists_in_IDXL(m,nodeId,minChunk,0);
1837             intMsg* imsg = meshMod[minChunk].lockRemoteNode(sharedIdx, index, 0, 0);
1838             done = imsg->i;
1839             delete imsg;
1840           }
1841           break;
1842         }
1843       }
1844     }
1845   }
1846   else if(FEM_Is_ghost_index(nodeId)) {
1847     //someone must have the lock
1848     for(int j=0; j<numchunks; j++) {
1849       int pchk = chunks1[j]->chk;
1850       int sharedIdx = m->getfmMM()->getfmUtil()->exists_in_IDXL(m,nodeId,pchk,2);
1851       intMsg* imsg = meshMod[pchk].hasLockRemoteNode(sharedIdx, index, 1);
1852       owner = imsg->i;
1853       delete imsg;
1854       if(owner != -1) { //this node is locked
1855         if(pchk == minChunk) {
1856           //the lock is good
1857           break;
1858         }
1859         else {
1860           //unlock the node on pchk & lock it on minChunk
1861           int locknodes = nodeId;
1862           int gotlocks = 1;
1863           int done = -1;
1864           int sharedIdx = m->getfmMM()->getfmUtil()->exists_in_IDXL(m,nodeId,pchk,2);
1865           intMsg* imsg = meshMod[pchk].unlockRemoteNode(sharedIdx, index, 1, 0);
1866           delete imsg;
1867           gotlocks=-1;
1868           while(done==-1) {
1869             done = m->getfmMM()->fmAdaptL->lockNodes(&gotlocks, &locknodes, 0, &locknodes, 1);
1870           }
1871           break;
1872         }
1873       }
1874     }
1875     if(owner==-1) {
1876       int locknodes = nodeId;
1877       int done = -1;
1878       int gotlocks=-1;
1879       while(done==-1) {
1880         done = m->getfmMM()->fmAdaptL->lockNodes(&gotlocks, &locknodes, 0, &locknodes, 1);
1881       }
1882     }
1883   }
1884   for(int j=0; j<numchunks; j++) {
1885     delete chunks1[j];
1886   }
1887   if(numchunks!=0) free(chunks1);
1888   return;
1889 }
1890 //========================Locking operations=======================
1891
1892
1893
1894
1895 //=======================Read/Write operations=======================
1896 /* The following are functions to read/write data from/to the mesh */
1897 /** wrapper function that takes a 'mesh pointer' instead of meshId, and gets/sets the 
1898     data corresponding to the this 'attr' of this 'entity'
1899     'firstItem' is the starting point in the list, 'length' is the number of entries
1900     'datatype' is the type of data while 'width' is the size of each data entry
1901 */
1902 void FEM_Mesh_dataP(FEM_Mesh *fem_mesh,int entity,int attr,void *data, int firstItem,int length, int datatype,int width) {
1903   IDXL_Layout lo(datatype,width);
1904   FEM_Mesh_data_layoutP(fem_mesh,entity,attr,data,firstItem,length,lo);
1905 }
1906
1907 /** similar to the above function, but takes a 'layout' instead of datatype and width
1908  */
1909 void FEM_Mesh_data_layoutP(FEM_Mesh *fem_mesh,int entity,int attr,void *data, int firstItem,int length, IDXL_Layout_t layout) {
1910   const char *caller="FEM_Mesh_data_layout";
1911   FEM_Mesh_data_layoutP(fem_mesh,entity,attr,data,firstItem,length,IDXL_Layout_List::get().get(layout,caller));
1912 }
1913
1914 /** The helper function that both the above functions call
1915     does the actual setting/getting of the data in/from the mesh
1916 */
1917 void FEM_Mesh_data_layoutP(FEM_Mesh *m,int entity,int attr,void *data, int firstItem,int length, const IDXL_Layout &layout) {
1918   const char *caller="FEM_Mesh_data";
1919   //FEMAPI(caller);
1920   //FEM_Mesh *m=FEM_Mesh_lookup(fem_mesh,caller);
1921   FEM_Attribute *a=m->lookup(entity,caller)->lookup(attr,caller);
1922   if (m->isSetting()) 
1923     a->set(data,firstItem,length,layout,caller);
1924   else /* m->isGetting()*/
1925     a->get(data,firstItem,length,layout,caller);
1926 }
1927
1928
1929
1930 /** Copy the essential attributes for this ghost node from remote chunks
1931     This would be only done for attributes of a node, currently for coords & boundary
1932  */
1933 void FEM_Ghost_Essential_attributes(FEM_Mesh *m, int coord_attr, int bc_attr, int nodeid) {
1934   femMeshModify *theMod = m->getfmMM();
1935   FEM_MUtil *util = theMod->getfmUtil();
1936   int index = theMod->idx;
1937   if(FEM_Is_ghost_index(nodeid)) {
1938     //build up the list of chunks it is ghost to
1939     int numchunks;
1940     IDXL_Share **chunks1;
1941     util->getChunkNos(0,nodeid,&numchunks,&chunks1);
1942     for(int j=0; j<numchunks; j++) {
1943       int chk = chunks1[j]->chk;
1944       if(chk==index) continue;
1945       int ghostidx = util->exists_in_IDXL(m,nodeid,chk,2);
1946       double2Msg *d = meshMod[chk].getRemoteCoord(index,ghostidx);
1947       intMsg *im = meshMod[chk].getRemoteBound(index,ghostidx);
1948       double *coord = new double[2];
1949       coord[0] = d->i; coord[1] = d->j;
1950       int bound = im->i;
1951       CkVec<FEM_Attribute *>*attrs = m->node.ghost->getAttrVec();
1952       for (int i=0; i<attrs->size(); i++) {
1953         FEM_Attribute *a = (FEM_Attribute *)(*attrs)[i];
1954         if (a->getAttr() == theMod->fmAdaptAlgs->coord_attr) {
1955           FEM_DataAttribute *d = (FEM_DataAttribute *)a;
1956           d->getDouble().setRow(FEM_From_ghost_index(nodeid),coord,0);
1957         }
1958         else if(a->getAttr() == FEM_BOUNDARY) {
1959           FEM_DataAttribute *d = (FEM_DataAttribute *)a;
1960           d->getInt().setRow(FEM_From_ghost_index(nodeid),bound);
1961         }
1962       }
1963       delete [] coord;
1964       for(int j=0; j<numchunks; j++) {
1965         delete chunks1[j];
1966       }
1967       if(numchunks != 0) free(chunks1);
1968       delete im;
1969       delete d;
1970       break;
1971     }
1972   }
1973   return;
1974 }
1975 //=======================Read/Write operations=======================
1976
1977
1978
1979
1980
1981
1982
1983 //=======================femMeshModify: shadow array=======================
1984 femMeshModify::femMeshModify(femMeshModMsg *fm) {
1985   numChunks = fm->numChunks;
1986   idx = fm->myChunk;
1987   fmLock = new FEM_lock(idx, this);
1988   fmUtil = new FEM_MUtil(idx, this);
1989   fmAdapt = NULL;
1990   fmAdaptL = NULL;
1991   fmAdaptAlgs = NULL;
1992   fmInp = NULL;
1993   fmMesh = NULL;
1994   delete fm;
1995 }
1996
1997 femMeshModify::femMeshModify(CkMigrateMessage *m) {
1998   tc = NULL;
1999   fmMesh = NULL;
2000   fmLock = new FEM_lock(this);
2001   fmUtil = new FEM_MUtil(this);
2002   fmInp = new FEM_Interpolate(this);
2003   fmAdapt = new FEM_Adapt(this);
2004   fmAdaptL = new FEM_AdaptL(this);
2005   fmAdaptAlgs = new FEM_Adapt_Algs(this);
2006 }
2007
2008 femMeshModify::~femMeshModify() {
2009   if(fmLock != NULL) {
2010     delete fmLock;
2011   }
2012   if(fmUtil != NULL) {
2013     delete fmUtil;
2014   }
2015 }
2016
2017
2018
2019 void femMeshModify::pup(PUP::er &p) {
2020   p|numChunks;
2021   p|idx;
2022   p|tproxy;
2023   fmLock->pup(p);
2024   p|fmLockN;
2025   p|fmIdxlLock;
2026   p|fmfixedNodes;
2027   fmUtil->pup(p);
2028   fmInp->pup(p);
2029   fmAdapt->pup(p);
2030   fmAdaptL->pup(p);
2031   fmAdaptAlgs->pup(p);
2032 }
2033
2034 enum {FEM_globalID=33};
2035 void femMeshModify::ckJustMigrated(void) {
2036   ArrayElement1D::ckJustMigrated();
2037   //set the pointer to fmMM
2038   tc = tproxy[idx].ckLocal();
2039   CkVec<TCharm::UserData> &v=tc->sud;
2040   FEM_chunk *c = (FEM_chunk*)(v[FEM_globalID].getData());
2041   fmMesh = c->getMesh("ckJustMigrated");
2042   fmMesh->fmMM = this;
2043   setPointersAfterMigrate(fmMesh);
2044 }
2045
2046 void femMeshModify::setPointersAfterMigrate(FEM_Mesh *m) {
2047   fmMesh = m;
2048   fmInp->FEM_InterpolateSetMesh(fmMesh);
2049   fmAdapt->FEM_AdaptSetMesh(fmMesh);
2050   fmAdaptL->FEM_AdaptLSetMesh(fmMesh);
2051   fmAdaptAlgs->FEM_AdaptAlgsSetMesh(fmMesh);
2052   for(int i=0; i<fmLockN.size(); i++) fmLockN[i].setMeshModify(this);
2053 }
2054
2055
2056
2057 /** Part of the initialization phase, create all the new objects
2058     Populate all data structures, include all locks
2059     It also computes all the fixed nodes and populates a data structure with it
2060  */
2061 void femMeshModify::setFemMesh(FEMMeshMsg *fm) {
2062   fmMesh = fm->m;
2063   tc = fm->t;
2064   tproxy = tc->getProxy();
2065   fmMesh->setFemMeshModify(this);
2066   fmAdapt = new FEM_Adapt(fmMesh, this);
2067   fmAdaptL = new FEM_AdaptL(fmMesh, this);
2068   fmAdaptAlgs = new FEM_Adapt_Algs(fmMesh, this);
2069   fmInp = new FEM_Interpolate(fmMesh, this);
2070   //populate the node locks
2071   int nsize = fmMesh->node.size();
2072   for(int i=0; i<nsize; i++) {
2073     fmLockN.push_back(FEM_lockN(i,this));
2074   }
2075   /*int gsize = fmMesh->node.ghost->size();
2076     for(int i=0; i<gsize; i++) {
2077     fmgLockN.push_back(new FEM_lockN(FEM_To_ghost_index(i),this));
2078     }*/
2079   for(int i=0; i<numChunks; i++) {
2080     fmIdxlLock.push_back(false);
2081   }
2082   //compute all the fixed nodes
2083   for(int i=0; i<nsize; i++) {
2084     if(fmAdaptL->isCorner(i)) {
2085       fmfixedNodes.push_back(i);
2086     }
2087   }
2088   delete fm;
2089   return;
2090 }
2091
2092
2093
2094 /* Coarse chunk locks */
2095 intMsg *femMeshModify::lockRemoteChunk(int2Msg *msg) {
2096   CtvAccess(_curTCharm) = tc;
2097   intMsg *imsg = new intMsg(0);
2098   int ret = fmLock->lock(msg->i, msg->j);
2099   imsg->i = ret;
2100   delete msg;
2101   return imsg;
2102 }
2103
2104 intMsg *femMeshModify::unlockRemoteChunk(int2Msg *msg) {
2105   CtvAccess(_curTCharm) = tc;
2106   intMsg *imsg = new intMsg(0);
2107   int ret = fmLock->unlock(msg->i, msg->j);
2108   imsg->i = ret;
2109   delete msg;
2110   return imsg;
2111 }
2112
2113 /* node locks */
2114 intMsg *femMeshModify::lockRemoteNode(int sharedIdx, int fromChk, int isGhost, int readLock) {
2115   CtvAccess(_curTCharm) = tc;
2116   int localIdx;
2117   if(isGhost) {
2118     localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 1);
2119   } else {
2120     localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 0);
2121   }
2122   //CkAssert(localIdx != -1);
2123   intMsg *imsg = new intMsg(0);
2124   int ret;
2125   if(localIdx == -1) {
2126     ret = -1;
2127   }
2128   else {
2129     if(readLock) {
2130       ret = fmLockN[localIdx].rlock();
2131     } else {
2132       ret = fmLockN[localIdx].wlock(fromChk);
2133     }
2134   }
2135   imsg->i = ret;
2136   return imsg;
2137 }
2138
2139 intMsg *femMeshModify::unlockRemoteNode(int sharedIdx, int fromChk, int isGhost, int readLock) {
2140   CtvAccess(_curTCharm) = tc;
2141   int localIdx;
2142   if(isGhost) {
2143     localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 1);
2144   } else {
2145     localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 0);
2146   }
2147   CkAssert(localIdx != -1);
2148   intMsg *imsg = new intMsg(0);
2149   int ret;
2150   if(readLock) {
2151     ret = fmLockN[localIdx].runlock();
2152   } else {
2153     ret = fmLockN[localIdx].wunlock(fromChk);
2154   }
2155   imsg->i = ret;
2156   return imsg;
2157 }
2158
2159
2160
2161 chunkListMsg *femMeshModify::getChunksSharingGhostNode(int2Msg *i2m) {
2162   CtvAccess(_curTCharm) = tc;
2163   chunkListMsg *clm;
2164   clm = fmUtil->getChunksSharingGhostNodeRemote(fmMesh, i2m->i, i2m->j);
2165   delete i2m;
2166   return clm;
2167 }
2168
2169
2170
2171 intMsg *femMeshModify::addNodeRemote(addNodeMsg *msg) {
2172   CtvAccess(_curTCharm) = tc;
2173   intMsg *imsg = new intMsg(-1);
2174   //translate the indices
2175   int *localIndices = (int*)malloc(msg->nBetween*sizeof(int));
2176   int *chunks = (int*)malloc(msg->numChunks*sizeof(int));
2177   CkAssert(msg->numChunks==1);
2178   chunks[0] = msg->chunks[0];
2179   for(int i=0; i<msg->nBetween; i++) {
2180     localIndices[i] = fmUtil->lookup_in_IDXL(fmMesh, msg->between[i], chunks[0], 0);
2181     CkAssert(localIndices[i] != -1);
2182   }
2183   int ret = FEM_add_node(fmMesh, localIndices, msg->nBetween, chunks, msg->numChunks, msg->forceShared);
2184   //this is a ghost on that chunk,
2185   //add it to the idxl & update that guys idxl list
2186   fmMesh->node.ghostSend.addNode(ret,chunks[0]);
2187   imsg->i = fmUtil->exists_in_IDXL(fmMesh,ret,chunks[0],1);
2188   free(localIndices);
2189   delete msg;
2190   return imsg;
2191 }
2192
2193 void femMeshModify::addSharedNodeRemote(sharedNodeMsg *fm) {
2194   CtvAccess(_curTCharm) = tc;
2195   FEM_add_shared_node_remote(fmMesh, fm->chk, fm->nBetween, fm->between);
2196   delete fm;
2197   return;
2198 }
2199
2200 void femMeshModify::removeSharedNodeRemote(removeSharedNodeMsg *fm) {
2201   CtvAccess(_curTCharm) = tc;
2202   fmUtil->removeNodeRemote(fmMesh, fm->chk, fm->index);
2203   delete fm;
2204   return;
2205 }
2206
2207 void femMeshModify::removeGhostNode(int fromChk, int sharedIdx) {
2208   CtvAccess(_curTCharm) = tc;
2209   fmUtil->removeGhostNodeRemote(fmMesh, fromChk, sharedIdx);
2210   return;
2211 }
2212
2213
2214
2215 void femMeshModify::addGhostElem(addGhostElemMsg *fm) {
2216   CtvAccess(_curTCharm) = tc;
2217   fmUtil->addGhostElementRemote(fmMesh, fm->chk, fm->elemType, fm->indices, fm->typeOfIndex, fm->connSize);
2218   delete fm;
2219   return;
2220 }
2221
2222 intMsg *femMeshModify::addElementRemote(addElemMsg *fm) {
2223   CtvAccess(_curTCharm) = tc;
2224   int newEl = fmUtil->addElemRemote(fmMesh, fm->chk, fm->elemtype, fm->connSize, fm->conn, fm->numGhostIndex, fm->ghostIndices);
2225   //expecting that this element will always be local to this chunk
2226   CkAssert(!FEM_Is_ghost_index(newEl));
2227   int shidx = fmUtil->exists_in_IDXL(fmMesh,newEl,fm->chk,3);
2228   intMsg *imsg = new intMsg(shidx);
2229   delete fm;
2230   return imsg;
2231 }
2232
2233 void femMeshModify::removeGhostElem(removeGhostElemMsg *fm) {
2234   CtvAccess(_curTCharm) = tc;
2235   fmUtil->removeGhostElementRemote(fmMesh, fm->chk, fm->elementid, fm->elemtype, fm->numGhostIndex, fm->ghostIndices, fm->numGhostRNIndex, fm->ghostRNIndices, fm->numGhostREIndex, fm->ghostREIndices, fm->numSharedIndex, fm->sharedIndices);
2236   delete fm;
2237   return;
2238 }
2239
2240 void femMeshModify::removeElementRemote(removeElemMsg *fm) {
2241   CtvAccess(_curTCharm) = tc;
2242   fmUtil->removeElemRemote(fmMesh, fm->chk, fm->elementid, fm->elemtype, fm->permanent);
2243   delete fm;
2244   return;
2245 }
2246
2247
2248
2249 intMsg *femMeshModify::eatIntoElement(int fromChk, int sharedIdx) {
2250   CtvAccess(_curTCharm) = tc;
2251   intMsg *imsg = new intMsg(-1);
2252   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 4);
2253   if(localIdx==-1) return imsg;
2254   int newEl = fmUtil->eatIntoElement(FEM_To_ghost_index(localIdx));
2255   int returnIdx = fmUtil->exists_in_IDXL(fmMesh, newEl, fromChk, 3);
2256   //returnIdx=-1 means that it was a discontinuous part
2257   if(returnIdx==-1) {
2258     //the other chunk has lost the entire region, unlock all the nodes involved
2259     //if there is a node 'fromChk' doesn't know about, 
2260     //but it is holding a lock on it, unlock it
2261     // nodesPerEl should be the number of nodes that can be adjacent to this element
2262     const int nodesPerEl = fmMesh->elem[0].getConn().width();
2263     int *adjnodes = new int[nodesPerEl];
2264     int numLocks = nodesPerEl + 1; //for 2D.. will be different for 3D
2265     int *lockednodes = new int[numLocks];
2266     fmMesh->e2n_getAll(newEl, adjnodes, 0);
2267     int newNode = -1;
2268     bool foundNewNode = false;
2269     //get the other node, which will be on an element that is an e2e of newEl
2270     //and which fromChk doesn't know abt, but has a lock on
2271     /*int *adjelems = new int[nodesPerEl];
2272     fmMesh->e2e_getAll(newEl, adjelems, 0);
2273     int *nnds = new int[nodesPerEl];
2274     for(int i=0; i<nodesPerEl; i++) {
2275       if(adjelems[i]!=-1) {
2276         fmMesh->e2n_getAll(adjelems[i], nnds, 0);
2277         for(int j=0; j<nodesPerEl; j++) {
2278           bool istarget = true;
2279           for(int k=0; k<nodesPerEl; k++) {
2280             if(nnds[j]==adjnodes[k]) {
2281               istarget = false;
2282             }
2283           }
2284           //get the owner of the lock
2285           if(istarget) {
2286             int lockowner = fmUtil->getLockOwner(nnds[j]);
2287             if(lockowner==fromChk) {
2288               bool knows = fmUtil->knowsAbtNode(fromChk,nnds[j]);
2289               //fromChk doesn't know abt this node
2290               if(!knows) {
2291                 newNode = nnds[j];
2292                 foundNewNode = true;
2293                 break;
2294               }
2295             }
2296           }
2297           if(foundNewNode) break;
2298         }
2299       }
2300       if(foundNewNode) break;
2301     }
2302     delete[] adjelems;
2303     delete[] nnds;
2304     */
2305     int *gotlocks = new int[numLocks];
2306     for(int i=0; i<nodesPerEl; i++) {
2307       lockednodes[i] = adjnodes[i];
2308       gotlocks[i] = 1;
2309     }
2310     if(foundNewNode) {
2311       lockednodes[nodesPerEl] = newNode;
2312       gotlocks[nodesPerEl] = 1;
2313     }
2314     else numLocks--;
2315     fmAdaptL->unlockNodes(gotlocks, lockednodes, 0, lockednodes, numLocks);
2316     delete[] lockednodes;
2317     delete[] gotlocks;
2318     delete[] adjnodes;
2319   }
2320   imsg->i = returnIdx;
2321   return imsg;
2322 }
2323
2324
2325
2326 intMsg *femMeshModify::getLockOwner(int fromChk, int sharedIdx) {
2327   CtvAccess(_curTCharm) = tc;
2328   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 0);
2329   int ret = fmUtil->getLockOwner(localIdx);
2330   intMsg *imsg = new intMsg(ret);
2331   return imsg;
2332 }
2333
2334 boolMsg *femMeshModify::knowsAbtNode(int fromChk, int toChk, int sharedIdx) {
2335   CtvAccess(_curTCharm) = tc;
2336   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 1);
2337   bool ret = fmUtil->knowsAbtNode(toChk,localIdx);
2338   boolMsg *bmsg = new boolMsg(ret);
2339   return bmsg;
2340 }
2341
2342
2343
2344 void femMeshModify::refine_flip_element_leb(int fromChk, int propElemT, int propNodeT,
2345                                             int newNodeT, int nbrOpNodeT,
2346                                             int nbrghost, double longEdgeLen)
2347 {
2348   CtvAccess(_curTCharm) = tc;
2349   int propElem, propNode, newNode, nbrOpNode;
2350   propElem = getfmUtil()->lookup_in_IDXL(fmMesh, propElemT, fromChk, 3, 0);
2351   propNode = getfmUtil()->lookup_in_IDXL(fmMesh, propNodeT, fromChk, 0, -1);
2352   newNode = getfmUtil()->lookup_in_IDXL(fmMesh, newNodeT, fromChk, 0, -1);
2353   if(nbrghost) {
2354     nbrOpNode = getfmUtil()->lookup_in_IDXL(fmMesh, nbrOpNodeT, fromChk, 1,-1); 
2355   }
2356   else {
2357     nbrOpNode = getfmUtil()->lookup_in_IDXL(fmMesh, nbrOpNodeT, fromChk, 0,-1); 
2358   }
2359   fmAdaptAlgs->refine_flip_element_leb(propElem, propNode, newNode, nbrOpNode, longEdgeLen);
2360   return;
2361 }
2362
2363
2364
2365 void femMeshModify::addToSharedList(int fromChk, int sharedIdx) {
2366   CtvAccess(_curTCharm) = tc;
2367   fmUtil->addToSharedList(fmMesh, fromChk, sharedIdx);
2368
2369
2370 void femMeshModify::updateAttrs(updateAttrsMsg* umsg) {
2371   CtvAccess(_curTCharm) = tc;
2372   int localIdx = -1;
2373   if(!umsg->isnode) {
2374     localIdx = fmUtil->lookup_in_IDXL(fmMesh, umsg->sharedIdx, umsg->fromChk, 3);
2375   }
2376   else {
2377     if(!umsg->isGhost) {
2378       localIdx = fmUtil->lookup_in_IDXL(fmMesh, umsg->sharedIdx, umsg->fromChk, 0);
2379     }
2380     else localIdx = fmUtil->lookup_in_IDXL(fmMesh, umsg->sharedIdx, umsg->fromChk, 1);
2381   }
2382   CkAssert(localIdx!=-1);
2383   fmUtil->updateAttrs(umsg->data,umsg->datasize,localIdx,umsg->isnode,umsg->elemType);
2384   delete umsg;
2385   return;
2386 }
2387
2388 void femMeshModify::updateghostsend(verifyghostsendMsg *vmsg) {
2389   CtvAccess(_curTCharm) = tc;
2390   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, vmsg->sharedIdx, vmsg->fromChk, 0);
2391   fmUtil->UpdateGhostSend(localIdx, vmsg->chunks, vmsg->numchks);
2392   delete vmsg;
2393 }
2394
2395 findgsMsg *femMeshModify::findghostsend(int fromChk, int sharedIdx) {
2396   CtvAccess(_curTCharm) = tc;
2397   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 0);
2398   int *chkl, numchkl=0;
2399   fmUtil->findGhostSend(localIdx, chkl, numchkl);
2400   findgsMsg *fmsg = new(numchkl)findgsMsg();
2401   fmsg->numchks = numchkl;
2402   for(int i=0; i<numchkl; i++) fmsg->chunks[i] = chkl[i];
2403   if(numchkl>0) delete[] chkl;
2404   return fmsg;
2405 }
2406
2407 intMsg *femMeshModify::getIdxGhostSend(int fromChk, int idxshared, int toChk) {
2408   CtvAccess(_curTCharm) = tc;
2409   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, idxshared, fromChk, 0);
2410   int idxghostsend = -1;
2411   if(localIdx != -1) { 
2412     const IDXL_Rec *irec = fmMesh->node.ghostSend.getRec(localIdx);
2413     if(irec) {
2414       for(int i=0; i<irec->getShared(); i++) {
2415         if(irec->getChk(i) == toChk) {
2416           idxghostsend = fmUtil->exists_in_IDXL(fmMesh, localIdx, toChk, 1);
2417           break;
2418         }
2419       }
2420     }
2421   }
2422   intMsg *d = new intMsg(idxghostsend);
2423   return d;
2424 }
2425
2426
2427
2428 double2Msg *femMeshModify::getRemoteCoord(int fromChk, int ghostIdx) {
2429   CtvAccess(_curTCharm) = tc;
2430   //int localIdx = fmUtil->lookup_in_IDXL(fmMesh, ghostIdx, fromChk, 1);
2431   if(ghostIdx == fmMesh->node.ghostSend.addList(fromChk).size()) {
2432     double2Msg *d = new double2Msg(-2.0,-2.0);
2433     return d;
2434   }
2435   else {
2436     int localIdx = fmMesh->node.ghostSend.addList(fromChk)[ghostIdx];
2437     double coord[2];
2438     FEM_Mesh_dataP(fmMesh, FEM_NODE, fmAdaptAlgs->coord_attr, coord, localIdx, 1, FEM_DOUBLE, 2);
2439     double2Msg *d = new double2Msg(coord[0], coord[1]);
2440     return d;
2441   }
2442 }
2443
2444 intMsg *femMeshModify::getRemoteBound(int fromChk, int ghostIdx) {
2445   CtvAccess(_curTCharm) = tc;
2446   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, ghostIdx, fromChk, 1);
2447   int bound;
2448   FEM_Mesh_dataP(fmMesh, FEM_NODE, FEM_BOUNDARY, &bound, localIdx, 1, FEM_INT, 1);
2449   intMsg *d = new intMsg(bound);
2450   return d;
2451 }
2452
2453
2454
2455 void femMeshModify::updateIdxlList(int fromChk, int idxTrans, int transChk) {
2456   CtvAccess(_curTCharm) = tc;
2457   int idxghostrecv = fmUtil->lookup_in_IDXL(fmMesh, idxTrans, transChk, 2);
2458   CkAssert(idxghostrecv != -1);
2459   fmMesh->node.ghost->ghostRecv.addNode(idxghostrecv,fromChk);
2460   return;
2461 }
2462
2463 void femMeshModify::removeIDXLRemote(int fromChk, int sharedIdx, int type) {
2464   CtvAccess(_curTCharm) = tc;
2465   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, type);
2466   //CkAssert(localIdx != -1);
2467   if(localIdx!=-1) {
2468     fmMesh->node.ghostSend.removeNode(localIdx,fromChk);
2469   }
2470   return;
2471 }
2472
2473 void femMeshModify::addTransIDXLRemote(int fromChk, int sharedIdx, int transChk) {
2474   CtvAccess(_curTCharm) = tc;
2475   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, transChk, 0);
2476   CkAssert(localIdx != -1);
2477   fmMesh->node.ghostSend.addNode(localIdx,fromChk);
2478   return;
2479 }
2480
2481 void femMeshModify::verifyIdxlList(int fromChk, int size, int type) {
2482   CtvAccess(_curTCharm) = tc;
2483   fmUtil->verifyIdxlListRemote(fmMesh, fromChk, size, type);
2484   return;
2485 }
2486
2487
2488
2489 void femMeshModify::idxllockRemote(int fromChk, int type) {
2490   CtvAccess(_curTCharm) = tc;
2491   if(type==1) type = 2;
2492   else if(type ==2) type = 1;
2493   else if(type ==3) type = 4;
2494   else if(type ==4) type = 3;
2495   fmUtil->idxllockLocal(fmMesh, fromChk, type);
2496   return;
2497 }
2498
2499 void femMeshModify::idxlunlockRemote(int fromChk, int type) {
2500   CtvAccess(_curTCharm) = tc;
2501   if(type==1) type = 2;
2502   else if(type ==2) type = 1;
2503   else if(type ==3) type = 4;
2504   else if(type ==4) type = 3;
2505   fmUtil->idxlunlockLocal(fmMesh, fromChk, type);
2506   return;
2507 }
2508
2509
2510
2511 intMsg *femMeshModify::hasLockRemoteNode(int sharedIdx, int fromChk, int isGhost) {
2512   CtvAccess(_curTCharm) = tc;
2513   int localIdx;
2514   if(isGhost) {
2515     localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 1);
2516   } else {
2517     localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 0);
2518   }
2519   CkAssert(localIdx != -1);
2520   intMsg *imsg = new intMsg(0);
2521   int ret = fmLockN[localIdx].lockOwner();
2522   imsg->i = ret;
2523   return imsg;
2524 }
2525
2526 void femMeshModify::modifyLockAll(int fromChk, int sharedIdx) {
2527   CtvAccess(_curTCharm) = tc;
2528   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 1);
2529   FEM_Modify_LockAll(fmMesh, localIdx);
2530   return;
2531 }
2532
2533 boolMsg *femMeshModify::verifyLock(int fromChk, int sharedIdx, int isGhost) {
2534   CtvAccess(_curTCharm) = tc;
2535   int localIdx;
2536   if(isGhost) {
2537     localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 1);
2538   } else {
2539     localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 0);
2540   }
2541   //CkAssert(localIdx != -1);
2542   boolMsg *bmsg = new boolMsg(0);
2543   bool ret;
2544   if(localIdx == -1) {
2545     ret = false;
2546   }
2547   else {
2548     ret = fmLockN[localIdx].verifyLock();
2549   }
2550   bmsg->b = ret;
2551   return bmsg;
2552 }
2553
2554 void femMeshModify::verifyghostsend(verifyghostsendMsg *vmsg) {
2555   CtvAccess(_curTCharm) = tc;
2556   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, vmsg->sharedIdx, vmsg->fromChk, 1);
2557   const IDXL_Rec *irec = fmMesh->node.shared.getRec(localIdx);
2558   if (irec!=NULL) {
2559     int numsh = irec->getShared();
2560     CkAssert(numsh==vmsg->numchks-1);
2561     for(int i=0; i<numsh; i++) {
2562       int ckl = irec->getChk(i);
2563       bool found = false;
2564       for(int j=0; j<numsh+1; j++) {
2565         if(vmsg->chunks[j]==ckl) {
2566           found = true; break;
2567         }
2568       }
2569       CkAssert(found);
2570     }
2571   }
2572   delete vmsg;
2573 }
2574
2575 boolMsg *femMeshModify::shouldLoseGhost(int fromChk, int sharedIdx, int toChk) {
2576   CtvAccess(_curTCharm) = tc;
2577   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 0);
2578   int *elems, numElems;
2579   fmMesh->n2e_getAll(localIdx, elems, numElems);
2580   bool shouldBeDeleted = true;
2581   for(int k=0; k<numElems; k++) {
2582     if(elems[k]>=0) {
2583       if(fmMesh->getfmMM()->fmUtil->exists_in_IDXL(fmMesh,elems[k],toChk,3)!=-1) {
2584         shouldBeDeleted = false;
2585         break;
2586       }
2587     }
2588   }
2589   if(numElems>0) delete[] elems;
2590   boolMsg *bmsg = new boolMsg(shouldBeDeleted);
2591   return bmsg;
2592 }
2593
2594
2595
2596 /** The node index is the sharedIdx with 'fromChk'
2597     This local node index is then to be send as a ghost to 'toChk'
2598 */
2599 void femMeshModify::addghostsendl(int fromChk, int sharedIdx, int toChk, int transIdx) {
2600   CtvAccess(_curTCharm) = tc;
2601   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 0);
2602   int sharedghost = fmUtil->exists_in_IDXL(fmMesh,localIdx,toChk,1);
2603   if(sharedghost==-1) { //it needs to be added from this chunk
2604     //lock idxl
2605     fmMesh->node.ghostSend.addNode(localIdx,toChk);
2606     meshMod[toChk].addghostsendl1(idx,fromChk,transIdx);
2607     //unlock idxl
2608   }
2609 }
2610
2611 /** The node index is obtained as a the ghost received from 'transChk' at 'transIdx'
2612  */
2613 void femMeshModify::addghostsendl1(int fromChk, int transChk, int transIdx) {
2614   CtvAccess(_curTCharm) = tc;
2615   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, transIdx, transChk, 2);
2616   fmMesh->node.ghost->ghostRecv.addNode(localIdx,fromChk);
2617 }
2618
2619 /** The node is the obtained on the ghost recv list from 'fromChk' at 'sharedIdx' index.
2620  */
2621 void femMeshModify::addghostsendr(int fromChk, int sharedIdx, int toChk, int transIdx) {
2622   CtvAccess(_curTCharm) = tc;
2623   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 2);
2624   int sharedghost = fmUtil->exists_in_IDXL(fmMesh,FEM_To_ghost_index(localIdx),toChk,2);
2625   if(sharedghost==-1) { //it needs to be added from this chunk
2626     //lock idxl
2627     fmUtil->idxllock(fmMesh,toChk,0);
2628     fmMesh->node.ghost->ghostRecv.addNode(localIdx,toChk);
2629     meshMod[toChk].addghostsendr1(idx,fromChk,transIdx);
2630     fmUtil->idxlunlock(fmMesh,toChk,0);
2631     //unlock idxl
2632   }
2633 }
2634
2635 /** The local node index is obtained on the shared list of 'transChk' at 'transIdx'
2636  */
2637 void femMeshModify::addghostsendr1(int fromChk, int transChk, int transIdx) {
2638   CtvAccess(_curTCharm) = tc;
2639   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, transIdx, transChk, 0);
2640   fmMesh->node.ghostSend.addNode(localIdx,fromChk);
2641 }
2642
2643 boolMsg *femMeshModify::willItLose(int fromChk, int sharedIdx) {
2644   CtvAccess(_curTCharm) = tc;
2645   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 3);
2646   int nnbrs[3];
2647   fmMesh->e2n_getAll(localIdx,nnbrs,0);
2648   //if it loses any node if it loses this element, then it should let fromChk, acquire this elem
2649   bool willlose = true;
2650   for(int i=0; i<3; i++) {
2651     int *enbrs, numenbrs;
2652     fmMesh->n2e_getAll(nnbrs[i],enbrs,numenbrs);
2653     willlose = true;
2654     for(int j=0; j<numenbrs; j++) {
2655       if(enbrs[j]>=0 && enbrs[j]!=localIdx) {
2656         willlose = false;
2657         break;
2658       }
2659     }
2660     if(numenbrs>0) delete [] enbrs;
2661     if(willlose) break;
2662   }
2663   boolMsg *bmsg = new boolMsg(willlose);
2664   return bmsg;
2665 }
2666
2667
2668
2669 /** The two local elements indices are obtained from the ghost send list to 'fromChk'
2670     Data is copied from the first element to the second
2671 */
2672 void femMeshModify::interpolateElemCopy(int fromChk, int sharedIdx1, int sharedIdx2) {
2673   CtvAccess(_curTCharm) = tc;
2674   int localIdx1 = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx1, fromChk, 3);
2675   int localIdx2 = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx2, fromChk, 3);
2676   CkAssert(localIdx1!=-1 && localIdx2!=-1);
2677   fmUtil->copyElemData(0,localIdx1,localIdx2);
2678 }
2679
2680 void femMeshModify::cleanupIDXL(int fromChk, int sharedIdx) {
2681   CtvAccess(_curTCharm) = tc;
2682   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 4);
2683   CkAssert(fmUtil->exists_in_IDXL(fmMesh,FEM_To_ghost_index(localIdx),fromChk,4)!=-1);
2684   fmMesh->elem[0].ghost->ghostRecv.removeNode(localIdx, fromChk);
2685   fmMesh->elem[0].getGhost()->set_invalid(localIdx,false);
2686 }
2687
2688 void femMeshModify::purgeElement(int fromChk, int sharedIdx) {
2689   CtvAccess(_curTCharm) = tc;
2690   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 3);
2691   CkAssert(localIdx!=-1);
2692   FEM_purge_element(fmMesh,localIdx,0);
2693 }
2694
2695 entDataMsg *femMeshModify::packEntData(int fromChk, int sharedIdx, bool isnode, int elemType) {
2696   CtvAccess(_curTCharm) = tc;
2697   int localIdx=-1;
2698   if(isnode) {
2699     localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 1);
2700   }
2701   else {
2702     localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 3);
2703   }
2704   CkAssert(localIdx!=-1);
2705   char *data; 
2706   int size, count;
2707   fmUtil->packEntData(&data, &size, &count, localIdx, isnode, elemType);
2708   entDataMsg *edm = new (size, 0) entDataMsg(count,size);
2709   memcpy(edm->data, data, size);
2710   free(data);
2711   return edm;
2712 }
2713
2714 boolMsg *femMeshModify::isFixedNodeRemote(int fromChk, int sharedIdx) {
2715   CtvAccess(_curTCharm) = tc;
2716   int localIdx = fmUtil->lookup_in_IDXL(fmMesh, sharedIdx, fromChk, 0);
2717   for(int i=0; i<fmfixedNodes.size(); i++) {
2718     if(fmfixedNodes[i]==localIdx) {
2719       return new boolMsg(true);
2720     }
2721   }
2722   return new boolMsg(false);
2723 }
2724
2725
2726
2727 //need these functions for debugging
2728 void femMeshModify::finish(void) {
2729   CkPrintf("[%d]finished this program!\n",idx);
2730   return;
2731 }
2732
2733 void femMeshModify::finish1(void) {
2734   meshMod[0].finish();
2735   return;
2736 }
2737
2738 #include "ParFUM_Adapt.def.h"
2739 /*@}*/