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