00001 
00002 
00003 
00004 
00005 
00006 #include "ParFUM.h"
00007 #include "ParFUM_internals.h"
00008 
00009 extern void splitEntity(IDXL_Side &c, int localIdx, int nBetween, int *between, int idxbase);
00010 
00011 
00012 FEM_MUtil::FEM_MUtil(int i, femMeshModify *m) {
00013   idx = i;
00014   mmod = m;
00015 }
00016 
00017 FEM_MUtil::FEM_MUtil(femMeshModify *m) {
00018   mmod = m;
00019 }
00020 
00021 FEM_MUtil::~FEM_MUtil() {
00022   outStandingMappings.removeAll();
00023 }
00024 
00025 
00026 
00030 int FEM_MUtil::getRemoteIdx(FEM_Mesh *m, int elementid, int elemtype) {
00031   CkAssert(elementid < -1);
00032   int ghostid = FEM_To_ghost_index(elementid);
00033   const IDXL_Rec *irec = m->elem[elemtype].ghost->ghostRecv.getRec(ghostid);
00034   int size = irec->getShared();
00035   CkAssert(size == 1);
00036   int remoteChunk = irec->getChk(0);
00037   int sharedIdx = irec->getIdx(0);
00038   return remoteChunk;
00039 }
00040 
00045 bool FEM_MUtil::isShared(int index) {
00046   const IDXL_Rec *irec = mmod->fmMesh->node.shared.getRec(index);
00047   if(irec != NULL) return true;
00048   return false;
00049 }
00050 
00055 int FEM_MUtil::exists_in_IDXL(FEM_Mesh *m, int localIdx, int chk, int type, int elemType) {
00056   int exists  = -1;
00057   IDXL_List ll;
00058   if(type == 0) { 
00059     ll = m->node.shared.addList(chk);
00060   }
00061   else if(type == 1) { 
00062     ll = m->node.ghostSend.addList(chk);
00063   }
00064   else if(type == 2) { 
00065     ll = m->node.ghost->ghostRecv.addList(chk);
00066     localIdx = FEM_To_ghost_index(localIdx);
00067   }
00068   else if(type == 3) { 
00069     ll = m->elem[elemType].ghostSend.addList(chk);
00070   }
00071   else if(type == 4) { 
00072     ll = m->elem[elemType].ghost->ghostRecv.addList(chk);
00073     localIdx = FEM_To_ghost_index(localIdx);
00074   }
00075   for(int w2=0; w2<ll.size(); w2++) {
00076     if(ll[w2] == localIdx) {
00077       exists = w2;
00078       break;
00079     }
00080   }
00081   return exists;
00082 }
00083 
00087 int FEM_MUtil::lookup_in_IDXL(FEM_Mesh *m, int sharedIdx, int chk, int type, int elemType) {
00088   int localIdx  = -1;
00089   IDXL_List ll;
00090   if(type == 0) { 
00091     ll = m->node.shared.addList(chk);
00092   }
00093   else if(type == 1) { 
00094     ll = m->node.ghostSend.addList(chk);
00095   }
00096   else if(type == 2) { 
00097     ll = m->node.ghost->ghostRecv.addList(chk);
00098   }
00099   else if(type == 3) { 
00100     ll = m->elem[elemType].ghostSend.addList(chk);
00101   }
00102   else if(type == 4) { 
00103     ll = m->elem[elemType].ghost->ghostRecv.addList(chk);
00104   }
00105   CkAssert(sharedIdx < ll.size());
00106   localIdx = ll[sharedIdx];
00107   return localIdx;
00108 }
00109 
00110 
00111 
00123 void FEM_MUtil::getChunkNos(int entType, int entNo, int *numChunks, IDXL_Share ***chunks, int elemType) {
00124   int type = 0; 
00125   if(entType == 0) { 
00126     
00127     if(FEM_Is_ghost_index(entNo)) type = 2;
00128     else if(isShared(entNo)) type = 1;
00129     else type = 0;
00130     if(type == 2) {
00131       int ghostid = FEM_To_ghost_index(entNo);
00132       int noShared = 0;
00133       const IDXL_Rec *irec = mmod->fmMesh->node.ghost->ghostRecv.getRec(ghostid);
00134       if(irec) {
00135     noShared = irec->getShared(); 
00136     
00137     *numChunks = noShared;
00138     *chunks = (IDXL_Share**)malloc((*numChunks)*sizeof(IDXL_Share*));
00139     int i=0;
00140     for(i=0; i<*numChunks; i++) {
00141       int chk = irec->getChk(i); 
00142       int index = -1; 
00143       (*chunks)[i] = new IDXL_Share(chk, index);
00144     }
00145     
00146       }
00147       else { 
00148     *numChunks = 0;
00149       }
00150     }
00151     else if(type == 1) {
00152       const IDXL_Rec *irec = mmod->fmMesh->node.shared.getRec(entNo);
00153       if(irec) {
00154     *numChunks = irec->getShared() + 1; 
00155     *chunks = (IDXL_Share**)malloc((*numChunks)*sizeof(IDXL_Share*));
00156     int i=0;
00157     for(i=0; i<*numChunks-1; i++) {
00158       int chk = irec->getChk(i);
00159       int index = irec->getIdx(i);
00160       (*chunks)[i] = new IDXL_Share(chk, index);
00161     }
00162     (*chunks)[i] = new IDXL_Share(idx, -1);
00163       }
00164       else { 
00165     *numChunks = 0;
00166       }
00167     }
00168     else if(type == 0) {
00169       *numChunks = 0;
00170       *chunks = (IDXL_Share**)malloc((*numChunks)*sizeof(IDXL_Share*));
00171       for(int i=0; i<*numChunks; i++) {
00172     int chk = idx; 
00173     int index = entNo;
00174     (*chunks)[i] = new IDXL_Share(chk, index);
00175       }
00176     }
00177   }
00178   else if(entType == 1) { 
00179     
00180     if(FEM_Is_ghost_index(entNo)) type = 2;
00181     else type = 0;
00182 
00183     if(type == 2) {
00184       int ghostid = FEM_To_ghost_index(entNo);
00185       const IDXL_Rec *irec = mmod->fmMesh->elem[elemType].ghost->ghostRecv.getRec(ghostid);
00186       if(irec) {
00187     *numChunks = irec->getShared(); 
00188     *chunks = (IDXL_Share**)malloc((*numChunks)*sizeof(IDXL_Share*));
00189     for(int i=0; i<*numChunks; i++) {
00190       int chk = irec->getChk(i);
00191       int index = irec->getIdx(i);
00192       (*chunks)[i] = new IDXL_Share(chk, index);
00193     }
00194       }
00195       else { 
00196     *numChunks = 0;
00197       }
00198     }
00199     else if(type == 0) {
00200       *numChunks = 0;
00201       *chunks = (IDXL_Share**)malloc((*numChunks)*sizeof(IDXL_Share*));
00202       for(int i=0; i<*numChunks; i++) {
00203     int chk = idx; 
00204     int index = entNo;
00205     (*chunks)[i] = new IDXL_Share(chk, index);
00206       }
00207     }
00208   }
00209   return;
00210 }
00211 
00221 void FEM_MUtil::buildChunkToNodeTable(int *nodetype, int sharedcount, int ghostcount, int localcount, int *conn, int connSize, CkVec<int> ***allShared, int *numSharedChunks, CkVec<int> **allChunks, int ***sharedConn) {
00222   if((sharedcount > 0 && ghostcount == 0) || (ghostcount > 0 && localcount == 0)) { 
00223     *allShared = (CkVec<int> **)malloc(connSize*sizeof(CkVec<int> *));
00224     for(int i=0; i<connSize; i++) {
00225       (*allShared)[i] = new CkVec<int>; 
00226       int numchunks;
00227       IDXL_Share **chunks1;
00228       if(nodetype[i] == 1) { 
00229     
00230     getChunkNos(0,conn[i],&numchunks,&chunks1);
00231       }
00232       else if(nodetype[i] == 2) { 
00233     
00234     getChunkNos(0,conn[i],&numchunks,&chunks1);
00235       }
00236       else if(nodetype[i] == 0) {
00237     numchunks = 1;
00238     (*allShared)[i]->push_back(getIdx());
00239       }
00240       if((nodetype[i] == 1) || (nodetype[i] == 2)) {
00241     for(int j=0; j<numchunks; j++) {
00242       (*allShared)[i]->push_back(chunks1[j]->chk);
00243     }
00244       }
00245       if(nodetype[i]==1 || nodetype[i]==2) {
00246     for(int j=0; j<numchunks; j++) {
00247       delete chunks1[j];
00248     }
00249     if(numchunks!=0) free(chunks1);
00250       }
00251     }
00252     
00253     *allChunks = new CkVec<int>;
00254     for(int i=0; i<connSize; i++) {
00255       for(int j=0; j<(*allShared)[i]->size(); j++) {
00256     int exists = 0;
00257     for(int k=0; k<(*allChunks)->size(); k++) {
00258       if((*(*allChunks))[k]==(*(*allShared)[i])[j]) {
00259         exists = 1;
00260         break;
00261       }
00262     }
00263     if(!exists) {
00264       (*allChunks)->push_back((*(*allShared)[i])[j]);
00265       (*numSharedChunks)++;
00266     }
00267       }
00268     }
00269     *sharedConn = (int **)malloc((*numSharedChunks)*sizeof(int *));
00270     for(int j=0; j<*numSharedChunks; j++) {
00271       (*sharedConn)[j] = (int *)malloc(connSize*sizeof(int));
00272     }
00273     int index = getIdx();
00274     for(int i=0; i<connSize; i++) {
00275       
00276       int chkindex = -1;
00277       if((nodetype[i] == 1) || (nodetype[i] == 2)) {
00278     for(int j=0; j<(*numSharedChunks); j++) {
00279       (*sharedConn)[j][i] = -1;
00280     }
00281     for(int j=0; j<(*allShared)[i]->size(); j++) {
00282       for(int k=0; k<*numSharedChunks; k++) {
00283         if((*(*allShared)[i])[j] == (*(*allChunks))[k]) chkindex = k;
00284       }
00285       (*sharedConn)[chkindex][i] = 1; 
00286     }
00287     if(nodetype[i] == 2) {
00288       for(int k=0; k<*numSharedChunks; k++) {
00289         if(index == (*(*allChunks))[k]) chkindex = k;
00290       }
00291       (*sharedConn)[chkindex][i] = 2;
00292       if((*allShared)[i]->size()==1) {
00293         for(int k=0; k<*numSharedChunks; k++) {
00294           if((*(*allShared)[i])[0] == (*(*allChunks))[k]) chkindex = k;
00295         }
00296         (*sharedConn)[chkindex][i] = 0;
00297       }
00298     }
00299       }
00300       else {
00301     
00302     for(int j=0; j<(*numSharedChunks); j++) {
00303       (*sharedConn)[j][i] = -1; 
00304     }
00305     for(int k=0; k<*numSharedChunks; k++) {
00306       if(index == (*(*allChunks))[k]) chkindex = k;
00307     }
00308     (*sharedConn)[chkindex][i] = 0;
00309       }
00310     }
00311   }
00312   CkAssert(*numSharedChunks>0);
00313   return;
00314 }
00315 
00321 chunkListMsg *FEM_MUtil::getChunksSharingGhostNodeRemote(FEM_Mesh *m, int chk, int sharedIdx) {
00322   const IDXL_List ll = m->node.ghostSend.addList(chk);
00323   int localIdx = ll[sharedIdx];
00324   int numChunkList = 0;
00325   const IDXL_Rec *tween = m->node.shared.getRec(localIdx);
00326   if(tween) {
00327     numChunkList = tween->getShared();
00328   }
00329   chunkListMsg *clm = new (numChunkList, numChunkList, 0) chunkListMsg;
00330   clm->numChunkList = numChunkList;
00331   for(int i=0; i<numChunkList; i++) {
00332     clm->chunkList[i] = tween->getChk(i);
00333     clm->indexList[i] = tween->getIdx(i);
00334   }
00335   
00336   for(int i=0; i<numChunkList; i++) {
00337     for(int j=i+1; j<numChunkList; j++) {
00338       if(clm->chunkList[i]> clm->chunkList[j]) {
00339     int tmp = clm->chunkList[i];
00340     clm->chunkList[i] = clm->chunkList[j];
00341     clm->chunkList[j] = tmp;
00342     tmp = clm->indexList[i];
00343     clm->indexList[i] = clm->indexList[j];
00344     clm->indexList[j] = tmp;
00345       }
00346     }
00347   }
00348   return clm;
00349 }
00350 
00351 
00352 
00357 void FEM_MUtil::splitEntityAll(FEM_Mesh *m, int localIdx, int nBetween, int *between)
00358 {
00359   
00360   IDXL_Side *c = &(m->node.shared);
00361   const IDXL_Rec **tween = (const IDXL_Rec **)malloc(nBetween*sizeof(IDXL_Rec *));
00362   
00363   
00364   tween[0] = c->getRec(between[0]);
00365   for (int zs=tween[0]->getShared()-1; zs>=0; zs--) {
00366     for (int w1=0; w1<nBetween; w1++) {
00367       tween[w1] = c->getRec(between[w1]);
00368     }
00369     int chk = tween[0]->getChk(zs);
00370     
00371     int w = 0;
00372     for (w=0; w<nBetween; w++) {
00373       if (!tween[w]->hasChk(chk)) {
00374     break;
00375       }
00376     }
00377     if (w == nBetween) {
00378       idxllock(m,chk,0); 
00379       c->addNode(localIdx,chk); 
00380       
00381       int *sharedIndices = (int *)malloc(nBetween*sizeof(int));
00382       const IDXL_List ll = m->node.shared.addList(chk);
00383       for(int w1=0; w1<nBetween; w1++) {
00384     for(int w2=0; w2<ll.size(); w2++) {
00385       if(ll[w2] == between[w1]) {
00386         sharedIndices[w1] = w2;
00387         break;
00388       }
00389     }
00390       }
00391       sharedNodeMsg *fm = new (nBetween, 0) sharedNodeMsg;
00392       fm->chk = mmod->idx;
00393       fm->nBetween = nBetween;
00394       for(int j=0; j<nBetween; j++) {
00395     fm->between[j] = sharedIndices[j];
00396       }
00397       meshMod[chk].addSharedNodeRemote(fm);
00398       idxlunlock(m,chk,0); 
00399       free(sharedIndices);
00400     }
00401   }
00402   free(tween);
00403   return;
00404 }
00405 
00411 void FEM_MUtil::splitEntitySharing(FEM_Mesh *m, int localIdx, int nBetween, int *between, int numChunks, int *chunks)
00412 {
00413   for(int i=0; i<numChunks; i++) {
00414     int chk = chunks[i];
00415     if(chk==idx) continue;
00416     sharedNodeMsg *fm = new (nBetween, 0) sharedNodeMsg;
00417     fm->chk = idx;
00418     fm->nBetween = nBetween;
00419     for(int j=0; j<nBetween; j++) {
00420       fm->between[j] = exists_in_IDXL(m,between[j],chk,0);
00421       CkAssert(fm->between[j]!=-1);
00422     }
00423     idxllock(m,chk,0); 
00424     m->node.shared.addNode(localIdx,chk); 
00425     meshMod[chk].addSharedNodeRemote(fm);
00426     idxlunlock(m,chk,0); 
00427   }
00428   return;
00429 }
00430 
00436 void FEM_MUtil::splitEntityRemote(FEM_Mesh *m, int chk, int localIdx, int nBetween, int *between)
00437 {
00438   
00439   int *localIndices = (int *)malloc(nBetween*sizeof(int));
00440   const IDXL_List ll = m->node.shared.addList(chk);
00441   for(int i=0; i<nBetween; i++) {
00442     localIndices[i] = ll[between[i]];
00443   }
00444   FEM_Interpolate *inp = m->getfmMM()->getfmInp();
00445   FEM_Interpolate::NodalArgs nm;
00446   nm.n = localIdx;
00447   for(int i=0; i<nBetween; i++) {
00448     nm.nodes[i] = localIndices[i];
00449   }
00450   nm.frac = 0.5;
00451   nm.addNode = true;
00452   inp->FEM_InterpolateNodeOnEdge(nm);
00453   m->node.shared.addNode(localIdx,chk);
00454   free(localIndices);
00455   return;
00456 }
00457 
00458 
00459 
00465 void FEM_MUtil::removeNodeAll(FEM_Mesh *m, int localIdx)
00466 {
00467   IDXL_Side *c = &(m->node.shared);
00468   const IDXL_Rec *tween = c->getRec(localIdx);
00469   int size = 0;
00470   if(tween) 
00471     size = tween->getShared();
00472   if(size>0) {
00473     int *schunks = (int*)malloc(size*sizeof(int));
00474     int *sidx = (int*)malloc(size*sizeof(int));
00475     for(int i=0; i<size; i++) {
00476       schunks[i] = tween->getChk(i);
00477       sidx[i] = tween->getIdx(i);
00478     }
00479     for(int i=0; i<size; i++) {
00480       removeSharedNodeMsg *fm = new removeSharedNodeMsg;
00481       fm->chk = mmod->idx;
00482       fm->index = sidx[i];
00483       meshMod[schunks[i]].removeSharedNodeRemote(fm);
00484       m->node.shared.removeNode(localIdx, schunks[i]);      
00485     }
00486     free(schunks);
00487     free(sidx);
00488   }
00489   
00490   FEM_remove_node_local(m,localIdx);
00491   return;
00492 }
00493 
00498 void FEM_MUtil::removeNodeRemote(FEM_Mesh *m, int chk, int sharedIdx) {
00499   int localIdx;
00500   const IDXL_List ll = m->node.shared.addList(chk);
00501   localIdx = ll[sharedIdx];
00502   m->node.shared.removeNode(localIdx, chk);
00503   FEM_remove_node_local(m,localIdx);
00504   return;
00505 }
00506 
00512 void FEM_MUtil::removeGhostNodeRemote(FEM_Mesh *m, int fromChk, int sharedIdx) {
00513   int localIdx = lookup_in_IDXL(m,sharedIdx,fromChk,2); 
00514   if(localIdx >= 0) {
00515     m->node.ghost->ghostRecv.removeNode(localIdx, fromChk);
00516     if(m->node.ghost->ghostRecv.getRec(localIdx)==NULL) {
00517       int ghostid = FEM_To_ghost_index(localIdx);
00518       int numAdjNodes, numAdjElts;
00519       int *adjNodes, *adjElts;
00520       m->n2n_getAll(ghostid, adjNodes, numAdjNodes);
00521       m->n2e_getAll(ghostid, adjElts, numAdjElts);
00522       
00523       if(!((numAdjNodes==0) && (numAdjElts==0))) {
00524     CkPrintf("Error: Node %d cannot be removed, it is connected to :\n",ghostid);
00525     FEM_Print_n2e(m,ghostid);
00526     FEM_Print_n2n(m,ghostid);
00527     
00528       }
00529       
00530       m->node.ghost->set_invalid(localIdx,true);
00531     }
00532     
00533   }
00534   return;
00535 }
00536 
00537 
00538 
00545 int FEM_MUtil::addElemRemote(FEM_Mesh *m, int chk, int elemtype, int connSize, int *conn, int numGhostIndex, int *ghostIndices) {
00546   
00547   
00548   
00549   const IDXL_List ll1 = m->node.ghostSend.addList(chk);
00550   const IDXL_List ll2 = m->node.shared.addList(chk);
00551   int *localIndices = (int *)malloc(connSize*sizeof(int));
00552   int j=0;
00553   int ghostsRemaining = numGhostIndex;
00554   for(int i=0; i<connSize; i++) {
00555     if(ghostsRemaining > 0) {
00556       if(ghostIndices[j] == i) {
00557     localIndices[i] = ll1[conn[i]];
00558     ghostsRemaining--;
00559     j++;
00560       }
00561       else {
00562     localIndices[i] = ll2[conn[i]];
00563       }
00564     }
00565     else {
00566       localIndices[i] = ll2[conn[i]];
00567     }
00568   }
00569   int newEl = FEM_add_element(m, localIndices, connSize, elemtype, idx);
00570   free(localIndices);
00571   return newEl;
00572 }
00573 
00585 void FEM_MUtil::addGhostElementRemote(FEM_Mesh *m, int chk, int elemType, int *indices, int *typeOfIndex, int connSize) {
00586   
00587   const IDXL_List ll1 = m->node.ghost->ghostRecv.addList(chk);
00588   
00589   const IDXL_List ll2 = m->node.shared.addList(chk);
00590   int *conn = (int *)malloc(connSize*sizeof(int));
00591   for(int i=0; i<connSize; i++) {
00592     if(typeOfIndex[i]==-1) {
00593       
00594       int newGhostNode = FEM_add_node_local(m, 1);
00595       m->node.ghost->ghostRecv.addNode(newGhostNode,chk);
00596       
00597       int sharedIdx = exists_in_IDXL(m,FEM_To_ghost_index(newGhostNode),chk,2);
00598       int2Msg *i2 = new int2Msg(idx, sharedIdx);
00599       chunkListMsg *clm = meshMod[chk].getChunksSharingGhostNode(i2);
00600       conn[i] = FEM_To_ghost_index(newGhostNode);
00601       for(int j=0; j<clm->numChunkList; j++) {
00602     if(clm->chunkList[j]==idx) continue;
00603     int chk1 = clm->chunkList[j];
00604     int sharedIdx1 = clm->indexList[j];
00605     idxllock(m,chk1,0);
00606     m->node.ghost->ghostRecv.addNode(newGhostNode,chk1);
00607     
00608     meshMod[chk1].addTransIDXLRemote(idx,sharedIdx1,chk);
00609     idxlunlock(m,chk1,0);
00610       }
00611       delete clm;
00612       FEM_Ghost_Essential_attributes(m, mmod->fmAdaptAlgs->coord_attr, FEM_BOUNDARY, conn[i]);
00613     }
00614     else if(typeOfIndex[i]==1) {
00615       conn[i] = FEM_To_ghost_index(ll1[indices[i]]);
00616     }
00617     else if(typeOfIndex[i]==0) {
00618       conn[i] = ll2[indices[i]];
00619     }
00620   }
00621   int newGhostElement = FEM_add_element_local(m, conn, connSize, elemType, 1);
00622   m->elem[elemType].ghost->ghostRecv.addNode(FEM_To_ghost_index(newGhostElement),chk);
00623   free(conn);
00624   return;
00625 }
00626 
00629 void FEM_MUtil::removeElemRemote(FEM_Mesh *m, int chk, int elementid, int elemtype, int permanent,
00630         bool aggressive_node_removal) {
00631   const IDXL_List ll = m->elem[elemtype].ghostSend.addList(chk);
00632   int localIdx = ll[elementid];
00633   FEM_remove_element(m, localIdx, elemtype, permanent, aggressive_node_removal);
00634   return;
00635 }
00636 
00647 void FEM_MUtil::removeGhostElementRemote(FEM_Mesh *m, int chk, int elementid, int elemtype, int numGhostIndex, int *ghostIndices, int numGhostRNIndex, int *ghostRNIndices, int numGhostREIndex, int *ghostREIndices, int numSharedIndex, int *sharedIndices) {
00648   
00649   
00650   int localIdx;
00651   const IDXL_List lgre = m->elem[elemtype].ghost->ghostRecv.addList(chk);
00652   localIdx = lgre[elementid];
00653   if(localIdx == -1) {
00654 #ifndef FEM_SILENT
00655     CkPrintf("Ghost element at shared index %d, already removed\n",elementid);
00656 #endif
00657     return;
00658   }
00659   
00660   
00661   FEM_remove_element_local(m, FEM_To_ghost_index(localIdx), elemtype);
00662   
00663   if(numGhostIndex > 0) {
00664     const IDXL_List lgrn = m->node.ghost->ghostRecv.addList(chk);
00665     for(int i=0; i<numGhostIndex; i++) {
00666       localIdx = lgrn[ghostIndices[i]];
00667       m->node.ghost->ghostRecv.removeNode(localIdx, chk); 
00668       if(m->node.ghost->ghostRecv.getRec(localIdx)) {
00669     
00670       }
00671       else {
00672     FEM_remove_node_local(m, FEM_To_ghost_index(localIdx));
00673       }
00674     }
00675   }
00676   if(numGhostREIndex > 0) {
00677     const IDXL_List lgse = m->elem[elemtype].ghostSend.addList(chk);
00678     for(int i=0; i<numGhostREIndex; i++) {
00679       localIdx = lgse[ghostREIndices[i]];
00680       m->elem[elemtype].ghostSend.removeNode(localIdx, chk); 
00681     }
00682   }
00683   if(numGhostRNIndex > 0) {
00684     const IDXL_List lgsn = m->node.ghostSend.addList(chk);
00685     for(int i=0; i<numGhostRNIndex; i++) {
00686       localIdx = lgsn[ghostRNIndices[i]];
00687       m->node.ghostSend.removeNode(localIdx, chk); 
00688     }
00689   }
00690   if(numSharedIndex > 0) { 
00691     const IDXL_List lsn = m->node.shared.addList(chk);
00692     for(int i=0; i<numSharedIndex; i++) {
00693       bool flag1 = true;
00694       if(sharedIndices[i]<-500000000) {
00695     
00696     const IDXL_List lgrn = m->node.ghost->ghostRecv.addList(chk);
00697     if(sharedIndices[i]<-1000000000) {
00698       sharedIndices[i] += 500000000;
00699       flag1 = false; 
00700     }
00701     sharedIndices[i] += 1000000000;
00702     
00703     localIdx = lgrn[sharedIndices[i]]; 
00704     
00705     int newN = Replace_node_local(m, FEM_To_ghost_index(localIdx), -1);
00706     
00707     
00708     if(flag1) m->node.ghostSend.addNode(newN, chk);
00709     m->node.ghost->ghostRecv.removeNode(localIdx, chk);
00710     outStandingMappings.push_back(tuple(FEM_To_ghost_index(localIdx),newN));
00711     FEM_remove_node_local(m,FEM_To_ghost_index(localIdx));
00712     
00713 #ifdef DEBUG_1
00714     CkPrintf("Corner node %d converted from ghost to local\n");
00715 #endif
00716     mmod->fmfixedNodes.push_back(newN);
00717       }
00718       else {
00719     if(sharedIndices[i]<0 && sharedIndices[i]>=-500000000) {
00720       sharedIndices[i] += 500000000;
00721       flag1 = false;
00722     }
00723     localIdx = lsn[sharedIndices[i]];
00724     m->node.shared.removeNode(localIdx, chk);
00725     if(flag1) m->node.ghostSend.addNode(localIdx, chk);
00726       }
00727     }
00728   }
00729   return;
00730 }
00731 
00732 
00733 
00739 int FEM_MUtil::Replace_node_local(FEM_Mesh *m, int oldIdx, int newIdx) {
00740 #ifdef CPSD
00741   bool dropLock = false;
00742 #endif
00743 
00744   if(newIdx==-1) {
00745     
00746     
00747     newIdx = m->node.size();
00748     m->node.setLength(newIdx+1); 
00749     m->node.set_valid(newIdx,true);   
00750     m->n2e_removeAll(newIdx);    
00751     m->n2n_removeAll(newIdx);    
00752     mmod->fmLockN.push_back(FEM_lockN(newIdx,mmod));
00753     mmod->fmLockN[newIdx].wlock(idx); 
00754 #ifdef CPSD
00755     dropLock = true;
00756 #endif
00757   }
00758   
00759   FEM_Interpolate *inp = mmod->getfmInp();
00760   inp->FEM_InterpolateCopyAttributes(oldIdx,newIdx);
00761   
00762   int *nnbrs;
00763   int nsize;
00764   m->n2n_getAll(oldIdx, nnbrs, nsize);
00765   m->n2n_removeAll(newIdx);
00766   for(int i=0; i<nsize; i++) {
00767     m->n2n_add(newIdx, nnbrs[i]);
00768     m->n2n_replace(nnbrs[i], oldIdx, newIdx);
00769   }
00770   int *enbrs;
00771   int esize;
00772   m->n2e_getAll(oldIdx, enbrs, esize);
00773   m->n2e_removeAll(newIdx);
00774   for(int i=0; i<esize; i++) {
00775     m->n2e_add(newIdx, enbrs[i]);
00776     m->e2n_replace(enbrs[i], oldIdx, newIdx, 0);
00777   }
00778   
00779   m->n2n_removeAll(oldIdx);
00780   m->n2e_removeAll(oldIdx);
00781   if(nsize>0) delete[] nnbrs;
00782   if(esize>0) delete[] enbrs;
00783 #ifdef CPSD
00784   if (dropLock)
00785     mmod->fmLockN[newIdx].wunlock(idx); 
00786 #endif
00787   return newIdx;  
00788 }
00789 
00793 void FEM_MUtil::addToSharedList(FEM_Mesh *m, int fromChk, int sharedIdx) {
00794   int elemType = 0;
00795   int connSize = m->elem[elemType].getConn().width();
00796   
00797   int localIdx = mmod->fmUtil->lookup_in_IDXL(m,sharedIdx,fromChk,1); 
00798   
00799   m->node.shared.addNode(localIdx, fromChk);
00800   m->node.ghostSend.removeNode(localIdx, fromChk);
00801   
00802   const IDXL_Rec *irec2 = m->node.ghostSend.getRec(localIdx);
00803   if(irec2!=NULL) {
00804     for(int i=0; i<irec2->getShared(); i++) {
00805       
00806       int pchk = irec2->getChk(i);
00807       int shidx = irec2->getIdx(i);
00808       
00809       meshMod[pchk].addghostsendr(idx,shidx,fromChk,exists_in_IDXL(m,localIdx,fromChk,0));
00810       
00811     }
00812   }
00813   int *enbrs;
00814   int esize;
00815   m->n2e_getAll(localIdx, enbrs, esize);
00816   for(int i=0; i<esize; i++) {
00817     if(enbrs[i] >= 0) { 
00818       
00819       int exists = mmod->fmUtil->exists_in_IDXL(m, enbrs[i], fromChk, 3);
00820       if(exists == -1) {
00821     m->elem[elemType].ghostSend.addNode(enbrs[i], fromChk);
00822     
00823     int *indices = (int*)malloc(connSize*sizeof(int));
00824     int *typeOfIndex = (int*)malloc(connSize*sizeof(int));
00825     int *nnbrs = (int*)malloc(connSize*sizeof(int));
00826     m->e2n_getAll(enbrs[i], nnbrs, elemType);
00827     for(int j=0; j<connSize; j++) {
00828       int sharedNode = mmod->fmUtil->exists_in_IDXL(m,nnbrs[j],fromChk,0);
00829       if(sharedNode == -1) {
00830         
00831         int sharedGhost = mmod->fmUtil->exists_in_IDXL(m,nnbrs[j],fromChk,1);
00832         if( sharedGhost == -1) {
00833           
00834           
00835           const IDXL_Rec *irec = m->node.shared.getRec(nnbrs[j]);
00836           if(irec) {
00837         int noShared = irec->getShared();
00838         for(int sharedck=0; sharedck<noShared; sharedck++) {
00839           int ckshared = irec->getChk(sharedck);
00840           int idxshared = irec->getIdx(sharedck);
00841           if(ckshared == fromChk) continue;
00842           CkAssert(fromChk!=idx && fromChk!=ckshared && ckshared!=idx);
00843           intMsg* imsg = meshMod[ckshared].getIdxGhostSend(idx,idxshared,fromChk);
00844           int idxghostsend = imsg->i;
00845           delete imsg;
00846           if(idxghostsend != -1) {
00847             m->node.ghostSend.addNode(nnbrs[j],fromChk);
00848             meshMod[fromChk].updateIdxlList(idx,idxghostsend,ckshared);
00849             sharedGhost = exists_in_IDXL(m,nnbrs[j],fromChk,1);
00850             CkAssert(sharedGhost != -1);
00851             break; 
00852           }
00853           
00854         }
00855           }
00856           
00857         }
00858         if( sharedGhost == -1) {
00859           
00860           indices[j] = nnbrs[j];
00861           typeOfIndex[j] = -1;
00862         }
00863         else {
00864           
00865           indices[j] = sharedGhost;
00866           typeOfIndex[j] = 1;
00867         }
00868       }
00869       else {
00870         
00871         indices[j] = sharedNode;
00872         typeOfIndex[j] = 0;
00873       }
00874     }
00875     
00876     addGhostElemMsg *fm = new (connSize, connSize, 0)addGhostElemMsg;
00877     fm->chk = getIdx();
00878     fm->elemType = elemType;
00879     for(int j=0; j<connSize; j++) {
00880       fm->indices[j] = indices[j];
00881       fm->typeOfIndex[j] = typeOfIndex[j];
00882       if(typeOfIndex[j]==-1) {
00883         m->node.ghostSend.addNode(indices[j],fromChk);
00884       }
00885     }
00886     fm->connSize = connSize;
00887     meshMod[fromChk].addGhostElem(fm); 
00888     for(int j=0; j<connSize; j++) {
00889       
00890       
00891       if(typeOfIndex[j]==-1) {
00892         const IDXL_Rec *irec1 = m->node.shared.getRec(indices[j]);
00893         if(irec1!=NULL) {
00894           for(int sh=0; sh<irec1->getShared(); sh++) {
00895         int transIdx = exists_in_IDXL(m,indices[j],fromChk,1);
00896         meshMod[irec1->getChk(sh)].addghostsendl(idx,irec1->getIdx(sh),fromChk,transIdx);
00897           }
00898         }
00899       }
00900     }
00901     free(indices);
00902     free(typeOfIndex);
00903     free(nnbrs);
00904       }
00905     }
00906   }
00907   delete[] enbrs;
00908   return;
00909 }
00910 
00917 int FEM_MUtil::eatIntoElement(int localIdx, bool aggressive_node_removal) {
00918   CkAssert(FEM_Is_ghost_index(localIdx));
00919   int nodesPerEl = mmod->fmMesh->elem[0].getConn().width();
00920   int *adjnodes = new int[nodesPerEl];
00921   int *oldnodes = new int[nodesPerEl];
00922   int elemtype = 0;
00923   mmod->fmMesh->e2n_getAll(localIdx, adjnodes, elemtype);
00924   CkPrintf("Chunk %d eating elem %d(%d,%d,%d) %d\n",idx,localIdx,adjnodes[0],adjnodes[1],adjnodes[2], aggressive_node_removal);
00925 #ifdef DEBUG_1
00926   CkPrintf("Chunk %d eating elem %d(%d,%d,%d)\n",idx,localIdx,adjnodes[0],adjnodes[1],adjnodes[2]);
00927 #endif
00928   
00929   int remChk = mmod->fmMesh->elem[0].ghost->ghostRecv.getRec(FEM_From_ghost_index(localIdx))->getChk(0);
00930 #ifndef CPSD
00931   for(int i=0; i<nodesPerEl; i++) {
00932     if(FEM_Is_ghost_index(adjnodes[i])) { 
00933       
00934       int sharedIdx = exists_in_IDXL(mmod->fmMesh,adjnodes[i],remChk,2);
00935       meshMod[remChk].modifyLockAll(idx,sharedIdx);
00936     }
00937   }
00938 #endif
00939 #ifdef DEBUG_1
00940   CkPrintf("eatIntoElement::remove\n");
00941 #endif
00942   FEM_remove_element(mmod->fmMesh,localIdx,elemtype,idx,aggressive_node_removal);
00943 #ifdef DEBUG_1
00944   CkPrintf("eatIntoElement::done removing\n");
00945 #endif
00946   for(int j=0; j<nodesPerEl; j++) oldnodes[j] = adjnodes[j];
00947   
00948   
00949   bool flag1 = false;
00950   for(int i=0; i<outStandingMappings.size(); i++) {
00951     for(int j=0; j<nodesPerEl; j++) {
00952       if(outStandingMappings[i].oldIdx==adjnodes[j]) {
00953     adjnodes[j] = outStandingMappings[i].newIdx;
00954     outStandingMappings.remove(i);
00955     flag1=true;
00956     break;
00957       }
00958     }
00959     if(flag1) break;
00960   }
00961 #ifdef DEBUG_1
00962   CkPrintf("eatIntoElement::add\n");
00963 #endif
00964   int newEl = FEM_add_element(mmod->fmMesh, adjnodes, nodesPerEl, elemtype, idx);
00965 #ifdef DEBUG_1
00966   CkPrintf("eatIntoElement::done adding\n");
00967 #endif
00968   copyElemData(0,localIdx,newEl); 
00969   FEM_purge_element(mmod->fmMesh,localIdx,elemtype);
00970 
00971   for(int i=0; i<nodesPerEl; i++) {
00972     if(adjnodes[i]!=oldnodes[i]) {
00973       
00974       FEM_Modify_LockUpdate(mmod->fmMesh,adjnodes[i]);
00975     }
00976   }
00977 
00978   delete [] adjnodes;
00979   delete [] oldnodes;
00980   return newEl;
00981 }
00982 
00983 
00984 
00988 bool FEM_MUtil::knowsAbtNode(int chk, int nodeId) {
00989   if(nodeId >= 0) {
00990     const IDXL_Rec *irec = mmod->fmMesh->node.ghostSend.getRec(nodeId);
00991     if(irec) {
00992       for(int i=0; i<irec->getShared(); i++) {
00993     if(irec->getChk(i) == chk) return true;
00994       }
00995     }
00996     irec = mmod->fmMesh->node.shared.getRec(nodeId);
00997     if(irec) {
00998       for(int i=0; i<irec->getShared(); i++) {
00999     if(irec->getChk(i) == chk) return true;
01000       }
01001     }
01002   }
01003   else { 
01004     
01005     int owner = mmod->fmMesh->node.ghost->ghostRecv.getRec(FEM_From_ghost_index(nodeId))->getChk(0);
01006     int sharedIdx = mmod->fmMesh->node.ghost->ghostRecv.getRec(FEM_From_ghost_index(nodeId))->getIdx(0);
01007     
01008     boolMsg* bmsg = meshMod[owner].knowsAbtNode(idx, chk, sharedIdx);
01009     bool b1 = bmsg->b;
01010     delete bmsg;
01011     return b1;
01012   }
01013   return false;
01014 }
01015 
01021 void FEM_MUtil::UpdateGhostSend(int nodeId, int *chunkl, int numchunkl) {
01022   if(nodeId<0) return;
01023   const IDXL_Rec *irec = mmod->fmMesh->node.ghostSend.getRec(nodeId);
01024   int numchunks=0;
01025   int *chunks1, *inds1;
01026   if(irec) {
01027     numchunks = irec->getShared();
01028     chunks1 = new int[numchunks];
01029     inds1 = new int[numchunks];
01030     for(int i=0; i<numchunks; i++) {
01031       chunks1[i] = irec->getChk(i);
01032       inds1[i] = irec->getIdx(i);
01033     }
01034   }
01035   for(int i=0; i<numchunks; i++) {
01036     bool shouldbeSent = false;
01037     for(int j=0; j<numchunkl; j++) { 
01038       if(chunks1[i]==chunkl[j]) {
01039     shouldbeSent=true;
01040     break;
01041       }
01042     }
01043     if(!shouldbeSent) {
01044       meshMod[chunks1[i]].removeGhostNode(idx, inds1[i]);
01045       mmod->fmMesh->node.ghostSend.removeNode(nodeId,chunks1[i]);
01046     }
01047   }
01048   if(numchunks>0) {
01049     delete[] chunks1;
01050     delete[] inds1;
01051   }
01052   return;
01053 }
01054 
01060 void FEM_MUtil::findGhostSend(int nodeId, int *&chunkl, int &numchunkl) {
01061   if(nodeId<0) return;
01062   CkVec<int> chkl;
01063   int *adjnds, numadjnds=0;
01064   mmod->fmMesh->n2n_getAll(nodeId, adjnds, numadjnds);
01065   for(int j=0; j<numadjnds; j++) {
01066     int node1 = adjnds[j];
01067     if(node1>=0) {
01068       const IDXL_Rec *irec = mmod->fmMesh->node.shared.getRec(node1);
01069       if(irec) {
01070     for(int k=0; k<irec->getShared(); k++) {
01071       int pchk = irec->getChk(k);
01072       
01073       
01074       if(exists_in_IDXL(mmod->fmMesh,nodeId,pchk,0)!=-1) continue;
01075       bool shouldadd=true;
01076       for(int l=0; l<chkl.size(); l++) {
01077         if(chkl[l]==pchk) {
01078           shouldadd=false;
01079           break;
01080         }
01081       }
01082       if(shouldadd) chkl.push_back(pchk);
01083     }
01084       }
01085     }
01086   }
01087   if(numadjnds>0) delete[] adjnds;
01088   numchunkl = chkl.size();
01089   if(numchunkl>0) 
01090     chunkl = new int[numchunkl];
01091   for(int i=0;i<numchunkl;i++) {
01092     chunkl[i] = chkl[i];
01093   }
01094   
01095   chkl.free();
01096   return;
01097 }
01098 
01099 
01100 
01104 void FEM_MUtil::copyElemData(int etype, int elemid, int newEl) {
01105   FEM_Interpolate *inp = mmod->getfmInp();
01106   FEM_Interpolate::ElementArgs em;
01107   em.e = newEl;
01108   em.oldElement = elemid;
01109   em.elType = etype;
01110   inp->FEM_InterpolateElementCopy(em);
01111 }
01112 
01119 void FEM_MUtil::packEntData(char **data, int *size, int *cnt, int localIdx, bool isnode, int elemType){
01120   CkVec<FEM_Attribute *>*entattrs;
01121   if(isnode) {
01122     entattrs = (mmod->fmMesh->node).getAttrVec();
01123   }
01124   else {
01125     entattrs = (mmod->fmMesh->elem[elemType]).getAttrVec();
01126   }
01127   int count = 0;
01128   PUP::sizer psizer;
01129   for(int j=0;j<entattrs->size();j++){
01130     FEM_Attribute *attr = (FEM_Attribute *)(*entattrs)[j];
01131     if(attr->getAttr() < FEM_ATTRIB_FIRST){ 
01132       
01133       attr->pupSingle(psizer, localIdx);
01134       count++;
01135     }
01136     else if(attr->getAttr()==FEM_MESH_SIZING || attr->getAttr()==FEM_BOUNDARY) {
01137       
01138       attr->pupSingle(psizer, localIdx);
01139       count++;
01140     }
01141   }
01142   *cnt = count;
01143   *size = psizer.size();
01144   *data = (char*)malloc((*size)*sizeof(char));
01145   PUP::toMem pmem(*data);
01146   for(int j=0;j<entattrs->size();j++){
01147     FEM_Attribute *attr = (FEM_Attribute *)(*entattrs)[j];
01148     if(attr->getAttr() < FEM_ATTRIB_FIRST){ 
01149       
01150       attr->pupSingle(pmem, localIdx);
01151     }
01152     else if(attr->getAttr()==FEM_MESH_SIZING || attr->getAttr()==FEM_BOUNDARY) {
01153       
01154       attr->pupSingle(pmem, localIdx);
01155     }
01156   }
01157   return;
01158 }
01159 
01163 void FEM_MUtil::updateAttrs(char *data, int size, int newIndex, bool isnode, int elemType) {
01164   PUP::fromMem pmem(data);
01165   int count=0;
01166   CkVec<FEM_Attribute *>*attrs;
01167   if(isnode) {
01168     attrs = (mmod->fmMesh->node).getAttrVec();
01169   }
01170   else {
01171     attrs = (mmod->fmMesh->elem[elemType]).getAttrVec();
01172   }
01173   for(int j=0;j<attrs->size();j++){
01174     FEM_Attribute *attr = (FEM_Attribute *)(*attrs)[j];
01175     if(attr->getAttr() < FEM_ATTRIB_FIRST){
01176       
01177       attr->pupSingle(pmem, newIndex);
01178       count++;
01179     }
01180     else if(attr->getAttr()==FEM_MESH_SIZING || attr->getAttr()==FEM_BOUNDARY) {
01181       
01182       attr->pupSingle(pmem,newIndex);
01183       count++;
01184     }
01185   }
01186   CkAssert(size==count);
01187   return;
01188 }
01189 
01190 
01191 
01197 int FEM_MUtil::getLockOwner(int nodeId) {
01198   int owner = -1;
01199   if(nodeId>=0) {
01200     const IDXL_Rec *irec = mmod->fmMesh->node.shared.getRec(nodeId);
01201     
01202     int minchunk = MAX_CHUNK;
01203     if(irec) {
01204       for(int i=0; i<irec->getShared(); i++) {
01205     int pchk = irec->getChk(i);
01206     if(pchk<minchunk) minchunk = pchk;
01207       }
01208     }
01209     else minchunk = idx;
01210     if(minchunk == idx) owner = mmod->fmLockN[nodeId].lockOwner();
01211     else {
01212       CkAssert(minchunk!=MAX_CHUNK);
01213       int sharedIdx = mmod->getfmUtil()->exists_in_IDXL(mmod->fmMesh,nodeId,minchunk,0);
01214       intMsg* imsg = meshMod[minchunk].hasLockRemoteNode(sharedIdx, idx, 0);
01215       owner = imsg->i;
01216       delete imsg;
01217     }
01218   }
01219   else {
01220     int otherchk = mmod->fmMesh->node.ghost->ghostRecv.getRec(FEM_From_ghost_index(nodeId))->getChk(0);
01221     int sharedIdx = mmod->fmMesh->node.ghost->ghostRecv.getRec(FEM_From_ghost_index(nodeId))->getIdx(0);
01222     intMsg* imsg1 = meshMod[otherchk].getLockOwner(idx, sharedIdx);
01223     owner = imsg1->i;
01224     delete imsg1;
01225   }
01226   
01227   return owner;
01228 }
01229 
01230 
01231 
01235 void FEM_MUtil::idxllock(FEM_Mesh *m, int chk, int type) {
01236   if(idx < chk) {
01237     idxllockLocal(m,chk,type);
01238   } else {
01239     meshMod[chk].idxllockRemote(idx,type);
01240   }
01241   return;
01242 }
01243 
01246 void FEM_MUtil::idxlunlock(FEM_Mesh *m, int chk, int type) {
01247   if(idx < chk) {
01248     idxlunlockLocal(m,chk,type);
01249   } else {
01250     meshMod[chk].idxlunlockRemote(idx,type);
01251   }
01252   return;
01253 }
01254 
01258 void FEM_MUtil::idxllockLocal(FEM_Mesh *m, int toChk, int type) {
01259   CkAssert(toChk>=0 && toChk<mmod->numChunks && toChk!=idx && type==0);
01260   while(mmod->fmIdxlLock[toChk + type] == true) {
01261     
01262     CthYield();
01263   }
01264   
01265   mmod->fmIdxlLock[toChk + type] = true;
01266 #ifdef DEBUG_IDXLLock
01267   CkPrintf("[%d]locked idxl lock %d!\n",idx,toChk+type);
01268 #endif
01269   return;
01270 }
01271 
01274 void FEM_MUtil::idxlunlockLocal(FEM_Mesh *m, int toChk, int type) {
01275   CkAssert(toChk>=0 && toChk<mmod->numChunks && toChk!=idx && type==0);
01276   CkAssert(mmod->fmIdxlLock[toChk + type] == true);
01277   
01278   mmod->fmIdxlLock[toChk + type] = false;
01279 #ifdef DEBUG_IDXLLock
01280   CkPrintf("[%d]unlocked idxl lock %d!\n",idx,toChk+type);
01281 #endif
01282   return;
01283 }
01284 
01285 
01286 
01289 void FEM_MUtil::FEM_Print_n2n(FEM_Mesh *m, int nodeid){
01290   CkPrintf("node %d is adjacent to nodes:", nodeid);
01291   int *adjnodes;
01292   int sz;
01293   m->n2n_getAll(nodeid, adjnodes, sz); 
01294   for(int i=0;i<sz;i++)
01295     CkPrintf(" %d", adjnodes[i]);
01296   if(sz!=0) delete[] adjnodes;  
01297   CkPrintf("\n");
01298 }
01299 
01302 void FEM_MUtil::FEM_Print_n2e(FEM_Mesh *m, int eid){
01303   CkPrintf("node %d is adjacent to elements:", eid);
01304   int *adjes;
01305   int sz;
01306   m->n2e_getAll(eid, adjes, sz);
01307   for(int i=0;i<sz;i++)
01308     CkPrintf(" %d", adjes[i]);
01309   if(sz!=0) delete[] adjes;
01310   CkPrintf("\n");
01311 }
01312 
01315 void FEM_MUtil::FEM_Print_e2n(FEM_Mesh *m, int eid){
01316   CkPrintf("element %d is adjacent to nodes:", eid);
01317   int consz = m->elem[0].getConn().width();
01318   int *adjns = new int[consz];
01319   m->e2n_getAll(eid, adjns, 0); 
01320   for(int i=0;i<consz;i++)
01321     CkPrintf(" %d", adjns[i]);
01322   CkPrintf("\n");
01323   delete [] adjns;
01324 }
01325 
01328 void FEM_MUtil::FEM_Print_e2e(FEM_Mesh *m, int eid){
01329   CkPrintf("element %d is adjacent to elements:", eid);
01330   int consz = m->elem[0].getConn().width();
01331   int *adjes = new int[consz];
01332   m->e2e_getAll(eid, adjes, 0); 
01333   for(int i=0;i<consz;i++)
01334     CkPrintf(" %d", adjes[i]);
01335   CkPrintf("\n");
01336   delete [] adjes;
01337 }
01338 
01343 void FEM_MUtil::FEM_Print_coords(FEM_Mesh *m, int nodeid) {
01344   double crds[2];
01345   int bound;
01346   if(!FEM_Is_ghost_index(nodeid)) {
01347     FEM_Mesh_dataP(m, FEM_NODE, mmod->fmAdaptAlgs->coord_attr, (void *)crds, nodeid, 1, FEM_DOUBLE, 2);
01348     FEM_Mesh_dataP(m, FEM_NODE, FEM_BOUNDARY, (void *)&bound, nodeid, 1, FEM_INT, 1);
01349   }
01350   else {
01351     int numchunks;
01352     IDXL_Share **chunks1;
01353     getChunkNos(0,nodeid,&numchunks,&chunks1);
01354     int index = mmod->idx;
01355     for(int j=0; j<numchunks; j++) {
01356       int chk = chunks1[j]->chk;
01357       if(chk==index) continue;
01358       int ghostidx = exists_in_IDXL(m,nodeid,chk,2);
01359       double2Msg *d = meshMod[chk].getRemoteCoord(index,ghostidx);
01360       intMsg *im = meshMod[chk].getRemoteBound(index,ghostidx);
01361       crds[0] = d->i;
01362       crds[1] = d->j;
01363       bound = im->i;
01364       for(int j=0; j<numchunks; j++) {
01365     delete chunks1[j];
01366       }
01367       if(numchunks != 0) free(chunks1);
01368       delete im;
01369       delete d;
01370       break;
01371     }
01372   }
01373 #ifndef FEM_SILENT  
01374   CkPrintf("node %d (%f,%f) and boundary %d\n",nodeid,crds[0],crds[1],bound);
01375 #endif
01376 }
01377 
01378 
01379 
01384 void FEM_MUtil::StructureTest(FEM_Mesh *m) {
01385     int noNodes = m->node.size();
01386     int noEle = m->elem[0].size();
01387     int noGhostEle = m->elem[0].ghost->size();
01388     int noGhostNodes = m->node.ghost->size();
01389     int wdt = m->elem[0].getConn().width();
01390     int *e2n = (int*)malloc(wdt*sizeof(int));
01391     for(int i=0; i<noEle; i++) {
01392         if(m->elem[0].is_valid(i)) {
01393             m->e2n_getAll(i,e2n,0);
01394             
01395             if(e2n[0]==e2n[1] || e2n[1]==e2n[2] || e2n[2]==e2n[0]) {
01396                 CkPrintf("ERROR: element %d, has connectivity (%d,%d,%d)\n",i,e2n[0],e2n[1],e2n[2]);
01397                 CkAssert(false);
01398             }
01399             
01400             if(e2n[0]<0 || e2n[1]<0 || e2n[2]<0) {
01401                 CkPrintf("ERROR: element %d, has connectivity (%d,%d,%d)\n",i,e2n[0],e2n[1],e2n[2]);
01402                 CkAssert(false);
01403             }
01404             for(int j=0; j<3; j++) {
01405                 
01406                 CkAssert(m->node.is_valid(e2n[j]));
01407             }
01408             int e2e[3];
01409             m->e2e_getAll(i,e2e,0);
01410             if((e2e[0]==e2e[1]) && (e2e[0]==e2e[2]) && (e2e[0]!=-1)) {
01411                 CkPrintf("ERROR: element %d, has e2e (%d,%d,%d)\n",i,e2e[0],e2e[1],e2e[2]);
01412                 CkAssert(false);
01413             }
01414         }
01415     }
01416     for(int i=0; i<noGhostEle; i++) {
01417         if(m->elem[0].ghost->is_valid(i)) {
01418             int ghostIndex = FEM_To_ghost_index(i);
01419             m->e2n_getAll(ghostIndex,e2n,0);
01420             if(e2n[0]==e2n[1] || e2n[1]==e2n[2] || e2n[2]==e2n[0]) {
01421                 
01422                 CkPrintf("ERROR: element %d, has connectivity (%d,%d,%d)\n",ghostIndex,e2n[0],e2n[1],e2n[2]);
01423                 CkAssert(false);
01424             }
01425             if(!(e2n[0]>=0 || e2n[1]>=0 || e2n[2]>=0)) {
01426                 
01427                 CkPrintf("ERROR: element %d, has connectivity (%d,%d,%d)\n",ghostIndex,e2n[0],e2n[1],e2n[2]);
01428                 CkAssert(false);
01429             }
01430             for(int j=0; j<3; j++) {
01431                 
01432                 if(e2n[j]>=0) CkAssert(m->node.is_valid(e2n[j]));
01433                 else CkAssert(m->node.ghost->is_valid(FEM_From_ghost_index(e2n[j])));
01434             }
01435             int e2e[3];
01436             m->e2e_getAll(ghostIndex,e2e,0);
01437             if((e2e[0]==e2e[1]) && (e2e[0]==e2e[2]) && (e2e[0]!=-1)) {
01438                 CkPrintf("ERROR: element %d, has e2e (%d,%d,%d)\n",i,e2e[0],e2e[1],e2e[2]);
01439                 CkAssert(false);
01440             }
01441         }
01442     }
01443     int *n2e, n2esize=0;
01444     int *n2n, n2nsize=0;
01445     for(int i=0; i<noNodes; i++) {
01446         if(m->node.is_valid(i)) {
01447             m->n2e_getAll(i,n2e,n2esize);
01448             m->n2n_getAll(i,n2n,n2nsize);
01449             if(n2esize>n2nsize) {
01450                 FEM_Print_coords(m,i);
01451                 FEM_Print_n2e(m,i);
01452                 FEM_Print_n2n(m,i);
01453                 CkPrintf("ERROR: local node %d, with inconsistent adjacency list\n",i);
01454                 CkAssert(false);
01455             }
01456         }
01457         else {
01458             continue;
01459         }
01460         if(n2esize > 0) {
01461             for(int j=0; j<n2esize; j++) {
01462                 CkAssert(n2e[j]!=-1);
01463                 if(FEM_Is_ghost_index(n2e[j])) CkAssert(m->elem[0].ghost->is_valid(FEM_From_ghost_index(n2e[j]))==1);
01464                 else CkAssert(m->elem[0].is_valid(n2e[j])==1);
01465             }
01466             m->e2n_getAll(n2e[0],e2n,0);
01467             
01468             bool done = false;
01469             for(int j=0; j<n2esize; j++) {
01470                 if(n2e[j] >= 0) {
01471                     done = true; 
01472                     break;
01473                 }
01474             }
01475             if(!done) {
01476                 FEM_Print_coords(m,i);
01477                 FEM_Print_n2e(m,i);
01478                 FEM_Print_n2n(m,i);
01479                 CkPrintf("ERROR: isolated local node %d, with no local element connectivity\n",i);
01480                 CkAssert(false);
01481             }
01482             
01483             int testnode = i;
01484             int startnode = (e2n[0]==testnode) ? e2n[1] : e2n[0];
01485             int othernode = (e2n[2]==testnode) ? e2n[1] : e2n[2];
01486             int previousnode = startnode;
01487             int nextnode = -1;
01488             int numdeadends = 0;
01489             int numunused = n2esize-1;
01490             n2e[0] = -1;
01491             for(int j=0; j<n2esize-1; j++) {
01492                 nextnode = -1;
01493                 for(int k=1; k<n2esize; k++) {
01494                     if(n2e[k]==-1) continue;
01495                     m->e2n_getAll(n2e[k],e2n,0);
01496                     if(e2n[0]==previousnode || e2n[1]==previousnode || e2n[2]==previousnode) {
01497                         nextnode = (e2n[0]==previousnode) ? ((e2n[1]==testnode)? e2n[2]:e2n[1]) : ((e2n[1]==previousnode)? ((e2n[0]==testnode)? e2n[2]:e2n[0]):((e2n[1]==testnode)? e2n[0]:e2n[1]));
01498                         previousnode = nextnode;
01499                         n2e[k] = -1;
01500                         numunused--;
01501                     }
01502                 }
01503                 if(nextnode==othernode && othernode!=-1) {
01504                     
01505                     break;
01506                 }
01507                 else if(nextnode==-1) {
01508                     
01509                     numdeadends++;
01510                     previousnode = othernode;
01511                     othernode = -1;
01512                 }
01513                 if(numdeadends>2 && numunused!=0) {
01514                     FEM_Print_coords(m,i);
01515                     FEM_Print_n2e(m,i);
01516                     FEM_Print_n2n(m,i);
01517                     CkPrintf("ERROR: cloud connectivity of node %d is discontinuous\n",i);
01518                     CkAssert(false);
01519                 }
01520             }
01521             if(n2esize>0) delete[] n2e; n2esize=0;
01522             
01523             int n2n1size = 0;
01524             int *n2n1 = (int*)malloc(n2nsize*sizeof(int));
01525             int n2n1Count = 0;
01526             m->n2e_getAll(i,n2e,n2esize);
01527             for(int j=0; j<n2esize; j++) {
01528                 CkAssert(n2e[j]!=-1);
01529                 
01530                 int e2n1[3];
01531                 m->e2n_getAll(n2e[j],e2n1,0);
01532                 if(e2n1[0]!=i && e2n1[1]!=i && e2n1[2]!=i) {
01533                     FEM_Print_coords(m,i);
01534                     FEM_Print_n2e(m,i);
01535                     FEM_Print_e2n(m,n2e[j]);
01536                     CkPrintf("ERROR: ghost elem %d & ghost node %d have inconsistent adjacency list\n",n2e[j],i);
01537                     CkAssert(false);
01538                 }
01539                 for(int k=0; k<3;k++) {
01540                     if(e2n1[k] == i) continue;
01541                     bool flag1 = true;
01542                     for(int l=0; l<n2n1Count; l++) {
01543                         if(e2n1[k] == n2n1[l]) flag1 = false;
01544                     }
01545                     if(flag1 && n2n1Count<n2nsize) { 
01546                         n2n1[n2n1Count] = e2n1[k];
01547                         n2n1Count++;
01548                     }
01549                 }
01550             }
01551             
01552             bool flag2 = true;
01553             if(n2n1Count!=n2nsize) flag2 = false;
01554             for(int j=0; j<n2n1Count; j++) {
01555                 bool flag1 = false;
01556                 for(int k=0; k<n2nsize; k++) {
01557                     if(n2n[k]==n2n1[j]) flag1 = true;
01558                 }
01559                 if(!flag1) {
01560                     flag2 = false;
01561                     break;
01562                 }
01563             }
01564             if(!flag2) {
01565                 FEM_Print_coords(m,i);
01566                 FEM_Print_n2n(m,i);
01567                 for(int l=0; l<n2esize; l++) FEM_Print_e2n(m,n2e[l]);
01568                 CkPrintf("ERROR: ghost node %d has inconsistent adjacency list\n",i);
01569                 CkAssert(false);
01570             }
01571 
01572             free(n2n1);
01573             if(n2esize>0) delete [] n2e; n2esize=0;
01574             if(n2nsize>0) delete [] n2n; n2nsize=0;
01575         }
01576     }
01577     if(n2esize>0) delete [] n2e; n2esize=0;
01578     if(n2nsize>0) delete [] n2n; n2nsize=0;
01579     for(int i=0; i<noGhostNodes; i++) {
01580         int ghostidx = FEM_To_ghost_index(i);
01581         if(m->node.ghost->is_valid(i)) {
01582             m->n2e_getAll(ghostidx,n2e,n2esize);
01583             m->n2n_getAll(ghostidx,n2n,n2nsize);
01584             bool done = false;
01585             for(int k=0;k<n2nsize;k++) {
01586                 if(n2n[k]>=0) {
01587                     done = true;
01588                     break;
01589                 }
01590             }
01591             if(n2esize>n2nsize || !done) {
01592                 FEM_Print_coords(m,ghostidx);
01593                 FEM_Print_n2e(m,ghostidx);
01594                 FEM_Print_n2n(m,ghostidx);
01595                 CkPrintf("ERROR: ghost node %d, with inconsistent adjacency list\n",ghostidx);
01596                 CkAssert(false);
01597             }
01598             if(n2esize > 0) {
01599                 
01600                 int n2n1size = 0;
01601                 int *n2n1 = (int*)malloc(n2nsize*sizeof(int));
01602                 int n2n1Count = 0;
01603                 for(int j=0; j<n2esize; j++) {
01604                     CkAssert(n2e[j]!=-1);
01605                     if(FEM_Is_ghost_index(n2e[j])) {
01606                         CkAssert(m->elem[0].ghost->is_valid(FEM_From_ghost_index(n2e[j]))==1);
01607                     }
01608                     else {
01609                         CkAssert(m->elem[0].is_valid(n2e[j])==1);
01610                     }
01611                     
01612                     int e2n1[3];
01613                     m->e2n_getAll(n2e[j],e2n1,0);
01614                     if(e2n1[0]!=ghostidx && e2n1[1]!=ghostidx && e2n1[2]!=ghostidx) {
01615                         FEM_Print_coords(m,ghostidx);
01616                         FEM_Print_n2e(m,ghostidx);
01617                         FEM_Print_e2n(m,n2e[j]);
01618                         CkPrintf("ERROR: ghost elem %d & ghost node %d have inconsistent adjacency list\n",n2e[j],ghostidx);
01619                         CkAssert(false);
01620                     }
01621                     for(int k=0; k<3;k++) {
01622                         if(e2n1[k] == ghostidx) continue;
01623                         bool flag1 = true;
01624                         for(int l=0; l<n2n1Count; l++) {
01625                             if(e2n1[k] == n2n1[l]) flag1 = false;
01626                         }
01627                         if(flag1 && n2n1Count<n2nsize) { 
01628                             n2n1[n2n1Count] = e2n1[k];
01629                             n2n1Count++;
01630                         }
01631                     }
01632                 }
01633                 
01634                 bool flag2 = true;
01635                 if(n2n1Count!=n2nsize) flag2 = false;
01636                 for(int j=0; j<n2n1Count; j++) {
01637                     bool flag1 = false;
01638                     for(int k=0; k<n2nsize; k++) {
01639                         if(n2n[k]==n2n1[j]) flag1 = true;
01640                     }
01641                     if(!flag1) {
01642                         flag2 = false;
01643                         break;
01644                     }
01645                 }
01646                 if(!flag2) {
01647                     FEM_Print_coords(m,ghostidx);
01648                     FEM_Print_n2n(m,ghostidx);
01649                     for(int l=0; l<n2esize; l++) FEM_Print_e2n(m,n2e[l]);
01650                     CkPrintf("ERROR: ghost node %d has inconsistent adjacency list\n",ghostidx);
01651                     CkAssert(false);
01652                 }
01653                 free(n2n1);
01654                 delete[] n2e;
01655             }
01656             if(n2nsize > 0) {
01657                 for(int j=0; j<n2nsize; j++) {
01658                     CkAssert(n2n[j]!=-1);
01659                 }
01660                 delete[] n2n;
01661             }
01662             
01663             const IDXL_Rec *irec = m->node.ghost->ghostRecv.getRec(i);
01664             
01665             int remChk = irec->getChk(0);
01666             int sharedIdx = irec->getIdx(0);
01667             int numsh = irec->getShared();
01668             verifyghostsendMsg *vmsg = new(numsh,0)verifyghostsendMsg();
01669             vmsg->fromChk = idx;
01670             vmsg->sharedIdx = sharedIdx;
01671             vmsg->numchks = numsh;
01672             for(int k=0; k<numsh; k++) vmsg->chunks[k] = irec->getChk(k);
01673             meshMod[remChk].verifyghostsend(vmsg);
01674         }
01675         else {
01676             continue;
01677         }
01678     }
01679     free(e2n);
01680     return;
01681 }
01682 
01685 int FEM_MUtil::AreaTest(FEM_Mesh *m) {
01686   int noEle = m->elem[0].size();
01687   int wdt = m->elem[0].getConn().width();
01688   int con[3];
01689   double n1_coord[2], n2_coord[2], n3_coord[3];
01690   double smallestarea = 1.0, smallestedge = 1.0, smallestalt = 1.0, largestQ=1.0;
01691   int worstEl = 0;
01692   for(int i=0; i<noEle; i++) {
01693     if(m->elem[0].is_valid(i)) {
01694       m->e2n_getAll(i,con,0);
01695       mmod->fmAdaptAlgs->getCoord(con[0],n1_coord);
01696       mmod->fmAdaptAlgs->getCoord(con[1],n2_coord);
01697       mmod->fmAdaptAlgs->getCoord(con[2],n3_coord);
01698       double area = mmod->fmAdaptAlgs->getSignedArea(n1_coord,n2_coord,n3_coord);
01699       double len1 = mmod->fmAdaptAlgs->length(n1_coord,n2_coord);
01700       double len2 = mmod->fmAdaptAlgs->length(n2_coord,n3_coord);
01701       double len3 = mmod->fmAdaptAlgs->length(n3_coord,n1_coord);
01702       double min = len1, max = len1;
01703       if(len2>max) max=len2;
01704       if(len3>max) max=len3;
01705       if(len2<min) min=len2;
01706       if(len3<min) min=len3;
01707       double shAl = fabs(area/max);
01708       double larR = max/shAl;
01709       if(fabs(area)<smallestarea) smallestarea = fabs(area);
01710       if(min<smallestedge) smallestedge = min;
01711       if(shAl<smallestalt) smallestalt = shAl;
01712       if(larR>largestQ) {
01713     largestQ = larR;
01714     worstEl = i;
01715       }
01716 #ifdef FEM_ELEMSORDERED
01717       if(-area < SLIVERAREA || larR>100.0) {
01718     CkAssert(false);
01719     return -1;
01720       }
01721 #else
01722       if(fabs(area) < SLIVERAREA || larR>100.0) {
01723     CkAssert(false);
01724     return -1;
01725       }
01726 #endif
01727     }
01728   }
01729 #ifdef DEBUG_QUALITY
01730   CkPrintf("SmallestArea %lf, SmallestEdge %lf, SmallestAlt %lf worstQuality %lf\n",smallestarea,smallestedge,smallestalt,largestQ);  
01731   m->e2n_getAll(worstEl,con,0);
01732   CkPrintf("WorstEl %d\n",worstEl);
01733   FEM_Print_coords(m,con[0]);
01734   FEM_Print_coords(m,con[1]);
01735   FEM_Print_coords(m,con[2]);
01736 #endif
01737   return 1;
01738 }
01739 
01743 int FEM_MUtil::IdxlListTest(FEM_Mesh *m) {
01744   for(int type=0; type <5; type++) {
01745     int listsize = 0;
01746     if(type==0) listsize = m->node.shared.size();
01747     else if(type==1) listsize = m->node.ghostSend.size();
01748     else if(type==2) listsize = m->node.ghost->ghostRecv.size();
01749     else if(type==3) listsize = m->elem[0].ghostSend.size();
01750     else if(type==4) listsize = m->elem[0].ghost->ghostRecv.size();
01751     for(int i=0; i<listsize; i++) {
01752       IDXL_List il;
01753       if(type==0) il = m->node.shared.getLocalList(i);
01754       else if(type==1) il = m->node.ghostSend.getLocalList(i);
01755       else if(type==2) il = m->node.ghost->ghostRecv.getLocalList(i);
01756       else if(type==3) il = m->elem[0].ghostSend.getLocalList(i);
01757       else if(type==4) il = m->elem[0].ghost->ghostRecv.getLocalList(i);
01758       int shck = il.getDest();
01759       int size = il.size();
01760       
01761       meshMod[shck].verifyIdxlList(idx,size,type);
01762       
01763       for(int j=0; j<size; j++) {
01764     CkAssert(il[j] >= -1);
01765       }
01766     }
01767   }
01768   return 1;
01769 }
01770 
01773 void FEM_MUtil::verifyIdxlListRemote(FEM_Mesh *m, int fromChk, int fsize, int type) {
01774   IDXL_Side is;
01775   IDXL_List il;
01776   if(type==0) il = m->node.shared.addList(fromChk);
01777   else if(type==1) il = m->node.ghost->ghostRecv.addList(fromChk);
01778   else if(type==2) il = m->node.ghostSend.addList(fromChk);
01779   else if(type==3) il = m->elem[0].ghost->ghostRecv.addList(fromChk);
01780   else if(type==4) il = m->elem[0].ghostSend.addList(fromChk);
01781   int size = il.size();
01782   CkAssert(fsize == size);
01783   
01784   for(int j=0; j<size; j++) {
01785     CkAssert(il[j] >= -1);
01786   }
01787   return;
01788 }
01789 
01793 int FEM_MUtil::residualLockTest(FEM_Mesh *m) {
01794   int noNodes = m->node.size();
01795   for(int i=0; i<noNodes; i++) {
01796     if(m->node.is_valid(i)) {
01797       if (mmod->fmLockN[i].haslocks()) {
01798           CkPrintf("[%d] Node %d has a residual lock\n", FEM_My_partition(), i);
01799           CkAssert(false);
01800       }
01801     }
01802   }
01803   for(int i=0; i<mmod->numChunks; i++) {
01804     CkAssert(mmod->fmIdxlLock[i]==false);
01805   }
01806   return 1;
01807 }
01808 
01813 void FEM_MUtil::unlockAll(FEM_Mesh *m) {
01814   int noNodes = m->node.size();
01815   for(int i=0; i<noNodes; i++) {
01816     if(m->node.is_valid(i)) {
01817       if (mmod->fmLockN[i].haslocks()) {
01818           mmod->fmLockN[i].reset(i, mmod);
01819       }
01820     }
01821   }
01822   for(int i=0; i<mmod->numChunks; i++) {
01823     mmod->fmIdxlLock[i] = false;
01824   }
01825 }
01826