Charm++: Charm Source Code Documentation

00001 /*
00002 Simple sample implementation of refinement interface.
00003 Orion Sky Lawlor, olawlor@acm.org, 4/9/2002
00004 Modified by Terry Wilmarth, wilmarth@cse.uiuc.edu, 4/16/2002
00005  */
00006 #include <stdlib.h>
00007 #include <unistd.h>
00008 #include <stdio.h>
00009 #include <math.h>
00010 #include "charm++.h"
00011 #include "charm-api.h"
00012 #include "tcharm.h"
00013 #include "mpi.h"
00014 #include "tri.h"
00015 #include "refine.h"
00016 
00017 /********************* Attach *****************/
00018 CLINKAGE void REFINE2D_Init(void) {
00019   TCHARM_API_TRACE("REFINE2D_Init", "refine");
00020   TCharm *tc=TCharm::get();
00021   
00022   // Create the refinement array
00023   MPI_Comm comm=MPI_COMM_WORLD; /* fixme: make user pass in communicator */
00024   int rank; MPI_Comm_rank(comm,&rank);
00025   CkArrayID refArrayID;
00026   if (rank==0) 
00027   { /* Master creates the (empty) array and broadcasts arrayID: */
00028     CkArrayOptions opts;
00029     opts.bindTo(tc->getProxy());
00030     refArrayID=CProxy_chunk::ckNew(new chunkMsg, opts);
00031   }
00032   MPI_Bcast(&refArrayID,sizeof(refArrayID),MPI_BYTE, 0,comm);
00033   mesh=refArrayID; /* set up readonly proxy */
00034   
00035   // Now everybody inserts their element into the array
00036   chunkMsg *cm = new chunkMsg;
00037   cm->nChunks = tc->getNumElements();
00038   cm->myThreads = tc->getProxy();
00039   mesh[rank].insert(cm);
00040   tc->suspend(); /* will resume from chunk constructor */
00041 }
00042 FLINKAGE void FTN_NAME(REFINE2D_INIT,refine2d_init)(void)
00043 {
00044   REFINE2D_Init();
00045 }
00046 
00047 /******************** NewMesh *******************/
00048 CLINKAGE void REFINE2D_NewMesh(int nEl,int nGhost,const int *conn,const int *gid)
00049 {
00050   TCHARM_API_TRACE("REFINE2D_NewMesh", "refine");
00051   if (!CtvAccess(_refineChunk))
00052     CkAbort("Forgot to call REFINE_Attach!!\n");
00053 
00054 /* FIXME: global barrier is a silly way to avoid early edge numbering messages */
00055   MPI_Barrier(MPI_COMM_WORLD);
00056   CtvAccess(_refineChunk)->newMesh(nEl,nGhost,conn, gid, 0);
00057   CkWaitQD(); //Wait for all edge numbering messages to finish
00058 }
00059 FLINKAGE void FTN_NAME(REFINE2D_NEWMESH,refine2d_newmesh)
00060 (int *nEl,int *nGhost,const int *conn,const int *gid)
00061 {
00062   TCHARM_API_TRACE("REFINE2D_NewMesh", "refine");
00063   if (!CtvAccess(_refineChunk))
00064     CkAbort("Forgot to call REFINE_Attach!!\n"); 
00065   
00066   CtvAccess(_refineChunk)->newMesh(*nEl, *nGhost,conn, gid, 1);
00067   CkWaitQD(); //Wait for all edge numbering messages to finish
00068 }
00069 
00070 /********************** Splitting ******************/
00071 class refineResults {
00072     int nResults;
00073     class resRec {
00074     public:
00075         int t,s,n;
00076         double f;
00077         int flag;
00078         resRec(int t_,int s_,int n_,double f_) 
00079             :t(t_), s(s_), n(n_), f(f_) {flag =0;}
00080         resRec(int t_,int s_,int n_,double f_,int flag_) 
00081             :t(t_), s(s_), n(n_), f(f_),flag(flag_) {}
00082     };
00083     std::vector<resRec> res;
00084     
00085   //Return anything on [0,2] than a or b
00086   int otherThan(int a,int b) {
00087     if (a==b) CkAbort("Opposite node is moving!");
00088     for (int i=0;i<3;i++)
00089         if (i!=a && i!=b) return i;
00090     CkAbort("Logic error in refine.C::otherThan");
00091     return -1;
00092   }
00093 public:
00094     refineResults(void) {nResults=0;}
00095     void add(int tri_,int side_,int n_,double frac_) {
00096         nResults++;
00097         res.push_back(resRec(tri_,side_,n_,frac_));
00098     }
00099     void add(int tri_,int side_,int n_,double frac_,int flag) {
00100         nResults++;
00101         res.push_back(resRec(tri_,side_,n_,frac_,flag));
00102     }
00103 
00104     int countResults(void) const {return nResults;}
00105     void extract(int i,const int *conn,int *triDest,int *A,int *B,int *C,double *fracDest,int idxBase,int *flags) {
00106         if ((i<0) || (i>=(int)res.size()))
00107             CkAbort("Invalid index in REFINE2D_Get_Splits");
00108         
00109         int tri=res[i].t;
00110         *triDest=tri+idxBase;
00111         int edgeOfTri=res[i].s;
00112         int movingNode=res[i].n;
00113         
00114         int c=(edgeOfTri+2)%3; //==opnode
00115         *A=conn[3*tri+movingNode]; //==othernode
00116         *B=conn[3*tri+otherThan(c,movingNode)];
00117         *C=conn[3*tri+c];
00118         *fracDest=res[i].f;
00119         *flags = res[i].flag;
00120         if (i==(int)res.size()-1) {
00121           delete this;
00122           chunk *C = CtvAccess(_refineChunk);
00123           C->refineResultsStorage=NULL;
00124         }
00125     }
00126 };
00127 
00128 class resultsRefineClient : public refineClient {
00129   refineResults *res;
00130 public:
00131   resultsRefineClient(refineResults *res_) :res(res_) {}
00132   void split(int tri, int side, int node, double frac) {
00133   #if 0
00134     //Convert from "tri.C edges" to sensible edges
00135     if (side==1) side=2;
00136     else if (side==2) side=1;
00137   #endif
00138     res->add(tri, side, node, frac);
00139   }
00140     void split(int tri, int side, int node, double frac,int flag) {
00141     res->add(tri, side, node, frac,flag);
00142   }
00143 
00144 };
00145 
00146 // this function should be called from a thread
00147 CLINKAGE void REFINE2D_Split(int nNode,double *coord,int nEl,double *desiredArea)
00148 {
00149   TCHARM_API_TRACE("REFINE2D_Split", "refine");
00150   chunk *C = CtvAccess(_refineChunk);
00151   if (!C)
00152     CkAbort("REFINE2D_Split failed> Did you forget to call REFINE2D_Attach?");
00153   C->refineResultsStorage=new refineResults;
00154   resultsRefineClient client(C->refineResultsStorage);
00155 
00156   C->updateNodeCoords(nNode, coord, nEl);
00157   C->multipleRefine(desiredArea, &client);
00158   CkWaitQD();
00159 }
00160 FLINKAGE void FTN_NAME(REFINE2D_SPLIT,refine2d_split)
00161    (int *nNode,double *coord,int *nEl,double *desiredArea)
00162 {
00163   REFINE2D_Split(*nNode,coord,*nEl,desiredArea);
00164 }
00165 
00166 static refineResults *getResults(void) {
00167   chunk *C = CtvAccess(_refineChunk);
00168   if (!C)
00169     CkAbort("Did you forget to call REFINE2D_Init?");
00170   refineResults *ret=C->refineResultsStorage;
00171   if (ret==NULL)
00172     CkAbort("Did you forget to call REFINE2D_Begin?");
00173   return ret;
00174 }
00175 
00176 CLINKAGE int REFINE2D_Get_Split_Length(void)
00177 {
00178   TCHARM_API_TRACE("REFINE2D_Get_Split_Length", "refine");
00179   return getResults()->countResults();
00180 }
00181 FLINKAGE int FTN_NAME(REFINE2D_GET_SPLIT_LENGTH,refine2d_get_split_length)(void)
00182 {
00183   return REFINE2D_Get_Split_Length();
00184 }
00185 
00186 CLINKAGE void REFINE2D_Get_Split
00187     (int splitNo,const int *conn,int *triDest,int *A,int *B,int *C,double *fracDest,int *flags)
00188 {
00189   TCHARM_API_TRACE("REFINE2D_Get_Split", "refine");
00190   refineResults *r=getResults();
00191   r->extract(splitNo,conn,triDest,A,B,C,fracDest,0,flags);
00192 }
00193 FLINKAGE void FTN_NAME(REFINE2D_GET_SPLIT,refine2d_get_split)
00194     (int *splitNo,const int *conn,int *triDest,int *A,int *B,int *C,double *fracDest, int *flags)
00195 {
00196   TCHARM_API_TRACE("REFINE2D_Get_Split", "refine");
00197   refineResults *r=getResults();
00198   r->extract(*splitNo-1,conn,triDest,A,B,C,fracDest,1,flags);
00199 }
00200 
00201 /********************* Check *****************/
00202 static int checkElement(chunk *C,const element &e,const int *uc,int idxBase)
00203 {
00204   int nMismatch=0;
00205   //Make sure the nodes are equal:
00206   /* This check no longer applies...
00207   for (int j=0;j<3;j++) {
00208     if (e.nodes[j].cid!=C->cid)
00209         CkAbort("REFINE2D_Check> Triangle has non-local node!\n");
00210     if (uc[j]-idxBase!=e.nodes[j].idx) {
00211       //Very bad-- mismatching indices
00212       nMismatch++;
00213     }
00214   }
00215   */
00216   return nMismatch;
00217 }
00218 
00219 static void checkConn(int nEl,const int *conn,int idxBase,int nNode)
00220 {
00221   chunk *C = CtvAccess(_refineChunk);
00222   if (!C) CkAbort("Did you forget to call REFINE2D_Attach?");
00223   C->sanityCheck();
00224   //This check stolen from updateNodeCoords:
00225   if (nEl != C->numElements || nNode != C->numNodes) {
00226     CkPrintf("ERROR: inconsistency in REFINE2D_Check on chunk %d:\n"
00227        "  your nEl (%d); my numElements (%d)\n"
00228        "  your nNode (%d); my numNodes (%d)\n",
00229        C->cid, nEl, C->numElements, nNode,C->numNodes);
00230     CkAbort("User code/library numbering inconsistency in REFINE2D");
00231   }
00232   int i;
00233   int nErrs=0;
00234   for (i=0;i<nEl;i++) 
00235   {
00236     const element &e=C->theElements[i];
00237     const int *uc=&conn[3*i];
00238     int elErrs=checkElement(C,e,uc,idxBase);
00239     nErrs+=elErrs;
00240     if (elErrs!=0) 
00241     { //A bad element-- print out debugging info:
00242       /* This check is no longer valid...
00243       CkError("REFINE2D Chunk %d>   Your triangle %d: %d %d %d\n"
00244               "REFINE2D Chunk %d> but my triangle %d: %d %d %d\n",
00245           C->cid+idxBase,i+idxBase, uc[0], uc[1], uc[2],
00246           C->cid+idxBase,i+idxBase,
00247           e.nodes[0].idx+idxBase,e.nodes[1].idx+idxBase,e.nodes[2].idx+idxBase);
00248       */
00249     }
00250   }
00251   if (nErrs!=0) {
00252     CkAbort("REFINE2D_Check> Major errors found. Exiting.");
00253   }
00254 }
00255 
00256 CLINKAGE void REFINE2D_Check(int nEl,const int *conn,int nNodes) {
00257   checkConn(nEl,conn,0,nNodes);
00258 }
00259 FLINKAGE void FTN_NAME(REFINE2D_CHECK,refine2d_check)
00260   (int *nEl,const int *conn,int *nNodes)
00261 {
00262   checkConn(*nEl,conn,1,*nNodes);
00263 }
00264 
00265 
00266
libs/ck-libs/tmr/refine.C
