AMPI/CUDA: Synchronous and asynchronous invocation, some cleanup
[charm.git] / src / libs / ck-libs / ampi / ampiimpl.h
1 #ifndef _AMPIIMPL_H
2 #define _AMPIIMPL_H
3
4 #include <string.h> /* for strlen */
5
6 #include "ampi.h"
7 #include "charm++.h"
8
9 #if AMPI_COMLIB
10 //#warning COMPILING IN UNTESTED AMPI COMLIB SUPPORT
11 #include "StreamingStrategy.h"
12 #include "EachToManyMulticastStrategy.h" /* for ComlibManager Strategy*/
13 #include "BroadcastStrategy.h"
14 #else
15 #define ComlibInstanceHandle int
16 #endif
17
18 #if 0
19 #define AMPI_DEBUG CkPrintf
20 #else
21 #define AMPI_DEBUG /* empty */
22 #endif
23
24 #ifdef AMPIMSGLOG
25
26 static int msgLogRank;
27 static int msgLogWrite;
28 static int msgLogRead;
29 static char *msgLogFilename;
30
31 #if CMK_PROJECTIONS_USE_ZLIB
32 #include <zlib.h>
33 namespace PUP{
34 class zdisk : public er {
35  protected:
36   gzFile F;//Disk file to read from/write to
37   zdisk(unsigned int type,gzFile f):er(type),F(f) {}
38   zdisk(const zdisk &p);                        //You don't want to copy
39   void operator=(const zdisk &p);       // You don't want to copy
40
41   //For seeking (pack/unpack in different orders)
42   virtual void impl_startSeek(seekBlock &s); /*Begin a seeking block*/
43   virtual int impl_tell(seekBlock &s); /*Give the current offset*/
44   virtual void impl_seek(seekBlock &s,int off); /*Seek to the given offset*/
45 };
46
47 //For packing to a disk file
48 class tozDisk : public zdisk {
49  protected:
50   //Generic bottleneck: pack n items of size itemSize from p.
51   virtual void bytes(void *p,int n,size_t itemSize,dataType t);
52  public:
53   //Write data to the given file pointer
54   // (must be opened for binary write)
55   // You must close the file yourself when done.
56   tozDisk(gzFile f):zdisk(IS_PACKING,f) {}
57 };
58
59 //For unpacking from a disk file
60 class fromzDisk : public zdisk {
61  protected:
62   //Generic bottleneck: unpack n items of size itemSize from p.
63   virtual void bytes(void *p,int n,size_t itemSize,dataType t);
64  public:
65   //Write data to the given file pointer 
66   // (must be opened for binary read)
67   // You must close the file yourself when done.
68   fromzDisk(gzFile f):zdisk(IS_UNPACKING,f) {}
69 };
70 }; // namespace PUP
71 #endif
72 #endif // AMPIMSGLOG
73
74 #define AMPI_COUNTER 0
75
76 #define AMPI_ALLTOALL_SHORT_MSG   32
77 #if CMK_CONVERSE_LAPI ||  CMK_BLUEGENE_CHARM
78 #define AMPI_ALLTOALL_MEDIUM_MSG   4194304
79 #else
80 #define AMPI_ALLTOALL_MEDIUM_MSG   32768
81 #endif
82
83 #if AMPI_COUNTER
84 class AmpiCounters{
85 public:
86         int send,recv,isend,irecv,barrier,bcast,gather,scatter,allgather,alltoall,reduce,allreduce,scan;
87         AmpiCounters(){
88                 send=0;recv=0;isend=0;irecv=0;barrier=0;bcast=0;gather=0;scatter=0;allgather=0;alltoall=0;reduce=0;allreduce=0;scan=0;
89         }
90         void output(int idx){
91                 printf("[%d]send=%d;recv=%d;isend=%d;irecv=%d;barrier=%d;bcast=%d;gather=%d;scatter=%d;allgather=%d;alltoall=%d;reduce=%d;allreduce=%d;scan=%d\n",idx,send,recv,isend,irecv,barrier,bcast,gather,scatter,allgather,alltoall,reduce,allreduce,scan);
92         }
93 };
94 #endif
95
96
97
98 void applyOp(MPI_Datatype datatype, MPI_Op op, int count, void* invec, void* inoutvec);
99 PUPfunctionpointer(MPI_Op)
100 class AmpiOpHeader {
101 public:
102   MPI_User_function* func;
103   MPI_Datatype dtype;  
104   int len;
105   int szdata;
106   AmpiOpHeader(MPI_User_function* f,MPI_Datatype d,int l,int szd):
107     func(f),dtype(d),len(l),szdata(szd) { }
108 };
109
110 //------------------- added by YAN for one-sided communication -----------
111 /* the index is unique within a communicator */
112 class WinStruct{
113 public:
114   MPI_Comm comm;
115   int index;
116   WinStruct(void):comm(MPI_COMM_NULL),index(-1){ }
117   WinStruct(MPI_Comm comm_, int index_):comm(comm_),index(index_){ }
118   void pup(PUP::er &p){ p|comm; p|index; }
119 };
120
121 class ampi;
122 class lockQueueEntry {
123 public:
124   int requestRank;
125   int pe_src; 
126   int ftHandle; 
127   int lock_type;
128   lockQueueEntry (int _requestRank, int _pe_src, int _ftHandle, int _lock_type) 
129         : requestRank(_requestRank), pe_src(_pe_src), ftHandle(_ftHandle),  lock_type(_lock_type) {}
130   lockQueueEntry () {}
131 };
132
133 typedef CkQ<lockQueueEntry *> LockQueue;
134
135 class win_obj {
136  public:
137   char* winName;
138   int winNameLeng;
139   int initflag; 
140   
141   void *baseAddr;
142   MPI_Aint winSize;
143   int disp_unit;
144   MPI_Comm comm;
145
146   int owner;   // Rank of owner of the lock, -1 if not locked
147   LockQueue lockQueue;  // queue of waiting processors for the lock
148                      // top of queue is the one holding the lock
149                      // queue is empty if lock is not applied
150                          
151   void setName(const char *src,int len);
152   void getName(char *src,int *len);
153   
154  public:
155   void pup(PUP::er &p); 
156
157   win_obj();
158   win_obj(char *name, void *base, MPI_Aint size, int disp_unit, MPI_Comm comm);
159   ~win_obj();
160   
161   int create(char *name, void *base, MPI_Aint size, int disp_unit, 
162              MPI_Comm comm);
163   int free();
164   
165   int put(void *orgaddr, int orgcnt, int orgunit, 
166           MPI_Aint targdisp, int targcnt, int targunit);
167
168   int get(void *orgaddr, int orgcnt, int orgunit, 
169           MPI_Aint targdisp, int targcnt, int targunit);
170   int accumulate(void *orgaddr, int orgcnt, MPI_Datatype orgtype, MPI_Aint targdisp, int targcnt, 
171                   MPI_Datatype targtype, MPI_Op op);
172
173   int iget(int orgcnt, MPI_Datatype orgtype,
174           MPI_Aint targdisp, int targcnt, MPI_Datatype targtype);
175   int igetWait(MPI_Request *req, MPI_Status *status);
176   int igetFree(MPI_Request *req, MPI_Status *status);
177
178   int fence();
179           
180   int lock(int requestRank, int pe_src, int ftHandle, int lock_type);
181   int unlock(int requestRank, int pe_src, int ftHandle);
182  
183   int wait();
184   int post();
185   int start();
186   int complete();
187
188   void lockTopQueue();
189   void enqueue(int requestRank, int pe_src, int ftHandle, int lock_type);       
190   void dequeue();
191   bool emptyQueue();
192 };
193 //-----------------------End of code by YAN ----------------------
194
195 class KeyvalPair{
196 protected:
197   int klen, vlen;
198   char* key;
199   char* val;
200 public:
201   KeyvalPair(void){ }
202   KeyvalPair(char* k, char* v);
203   ~KeyvalPair(void);
204   void pup(PUP::er& p){
205     p|klen; p|vlen;
206     if(p.isUnpacking()){
207       key=new char[klen];
208       val=new char[vlen];
209     }
210     p(key,klen);
211     p(val,vlen);
212   }
213 friend class InfoStruct;
214 };
215
216 class InfoStruct{
217   CkPupPtrVec<KeyvalPair> nodes;
218   bool valid;
219 public:
220   InfoStruct(void):valid(true) { }
221   void setvalid(bool valid_){ valid = valid_; }
222   bool getvalid(void){ return valid; }
223   void set(char* k, char* v);
224   void dup(InfoStruct& src);
225   int get(char* k, int vl, char*& v); // return flag
226   int deletek(char* k); // return -1 when not found
227   int get_valuelen(char* k, int* vl); // return flag
228   int get_nkeys(void) { return nodes.size(); }
229   int get_nthkey(int n,char* k);
230   void myfree(void);
231   void pup(PUP::er& p);
232 };
233
234 class CProxy_ampi;
235 class CProxyElement_ampi;
236
237 //Describes an AMPI communicator
238 class ampiCommStruct {
239         MPI_Comm comm; //Communicator
240         CkArrayID ampiID; //ID of corresponding ampi array
241         int size; //Number of processes in communicator
242         int isWorld; //1 if ranks are 0..size-1?
243         int isInter; // 0: intra-communicator; 1: inter-communicator
244         CkVec<int> indices;  //indices[r] gives the array index for rank r
245         CkVec<int> remoteIndices;  // remote group for inter-communicator
246         // cartesian virtual topology parameters
247         int ndims;
248         CkVec<int> dims;
249         CkVec<int> periods;
250         
251         // For communicator attributes (MPI_Attr_get): indexed by keyval
252         CkVec<void *> keyvals;
253
254         // graph virtual topology parameters
255         int nvertices;
256         CkVec<int> index;
257         CkVec<int> edges;
258
259         // Lazily fill world communicator indices
260         void makeWorldIndices(void) const {
261                 // cast away constness of "index" list
262           CkVec<int> *ind=(CkVec<int> *)&indices;  // changed by Isaac (as a guess to fix a bug). Was "index" not "indices"
263                 for (int i=0;i<size;i++) ind->push_back(i);
264         }
265 public:
266         ampiCommStruct(int ignored=0) {size=-1;isWorld=-1;isInter=0;}
267         ampiCommStruct(MPI_Comm comm_,const CkArrayID &id_,int size_)
268                 :comm(comm_), ampiID(id_),size(size_), isWorld(1), isInter(0) {}
269         ampiCommStruct(MPI_Comm comm_,const CkArrayID &id_,
270                 int size_,const CkVec<int> &indices_)
271                 :comm(comm_), ampiID(id_),size(size_),isInter(0),
272                  isWorld(0), indices(indices_) {}
273         ampiCommStruct(MPI_Comm comm_,const CkArrayID &id_,
274                 int size_,const CkVec<int> &indices_,
275                 const CkVec<int> &remoteIndices_)
276                 :comm(comm_),ampiID(id_),size(size_),isWorld(0),isInter(1),
277                 indices(indices_),remoteIndices(remoteIndices_) {}
278         void setArrayID(const CkArrayID &nID) {ampiID=nID;}
279
280         MPI_Comm getComm(void) const {return comm;}
281         const CkVec<int> &getIndices(void) const {
282                 if (isWorld && indices.size()!=size) makeWorldIndices();
283                 return indices;
284         }
285         const CkVec<int> &getRemoteIndices(void) const {return remoteIndices;}
286         CkVec<void *> &getKeyvals(void) {return keyvals;}
287
288         //Get the proxy for the entire array
289         CProxy_ampi getProxy(void) const;
290
291         //Get the array index for rank r in this communicator
292         int getIndexForRank(int r) const {
293 #ifndef CMK_OPTIMIZE
294                 if (r>=size) CkAbort("AMPI> You passed in an out-of-bounds process rank!");
295 #endif
296                 if (isWorld) return r;
297                 else return indices[r];
298         }
299         int getIndexForRemoteRank(int r) const {
300 #ifndef CMK_OPTIMIZE
301                 if (r>=remoteIndices.size()) CkAbort("AMPI> You passed in an out-of-bounds process rank!");
302 #endif
303                 if (isWorld) return r;
304                 else return remoteIndices[r];
305         }
306         //Get the rank for this array index (Warning: linear time)
307         int getRankForIndex(int i) const {
308                 if (isWorld) return i;
309                 else {
310                   for (int r=0;r<indices.size();r++)
311                     if (indices[r]==i) return r;
312                   return -1; /*That index isn't in this communicator*/
313                 }
314         }
315
316         int getSize(void) const {return size;}
317
318         inline int isinter(void) const { return isInter; }
319         inline const CkVec<int> &getindices() const {
320                 if (isWorld && indices.size()!=size) makeWorldIndices();
321                 return indices;
322         }
323         inline const CkVec<int> &getdims() const {return dims;}
324         inline const CkVec<int> &getperiods() const {return periods;}
325
326         inline int getndims() {return ndims;}
327         inline void setndims(int ndims_) {ndims = ndims_; }
328         inline void setdims(const CkVec<int> &dims_) { dims = dims_; }
329         inline void setperiods(const CkVec<int> &periods_) { periods = periods_; }
330
331         /* Similar hack for graph vt */
332         inline int getnvertices() {return nvertices;}
333         inline const CkVec<int> &getindex() const {return index;}
334         inline const CkVec<int> &getedges() const {return edges;}
335
336         inline void setnvertices(int nvertices_) {nvertices = nvertices_; }
337         inline void setindex(const CkVec<int> &index_) { index = index_; }
338         inline void setedges(const CkVec<int> &edges_) { edges = edges_; }
339
340         void pup(PUP::er &p) {
341                 p|comm;
342                 p|ampiID;
343                 p|size;
344                 p|isWorld;
345                 p|isInter;
346                 p|indices;
347                 p|remoteIndices;
348                 p|ndims;
349                 p|dims;
350                 p|periods;
351                 p|nvertices;
352                 p|index;
353                 p|edges;
354         }
355 };
356 PUPmarshall(ampiCommStruct)
357
358 struct mpi_comm_world
359 {
360   mpi_comm_world(const mpi_comm_world &m); //DO NOT USE
361   void operator=(const mpi_comm_world &m);
362   char *name; //new'd human-readable zero-terminated string name, or NULL
363 public:
364   ampiCommStruct comm;
365   mpi_comm_world() {
366     name=NULL;
367   }
368   ~mpi_comm_world() {
369           if (name) { delete[] name; name=0; }
370   }
371   void setName(const char *src) {
372     setName(src,strlen(src));
373   }
374   void setName(const char *src,int len) {
375         name=new char[len+1];
376         memcpy(name,src,len);
377         name[len] = '\0';
378   }
379   const char *getName(void) const { return name; }
380   void pup(PUP::er &p) {
381     p|comm;
382     int len=0;
383     if (name!=NULL) len=strlen(name)+1;
384     p|len;
385     if (p.isUnpacking()) name=new char[len];
386     p(name,len);
387   }
388 };
389 class mpi_comm_worlds {
390         mpi_comm_world s[MPI_MAX_COMM_WORLDS];
391 public:
392         mpi_comm_world &operator[](int i) {return s[i];}
393         void pup(PUP::er &p) {
394                 for (int i=0;i<MPI_MAX_COMM_WORLDS;i++)
395                         s[i].pup(p);
396         }
397 };
398
399 typedef CkVec<int> groupStruct;
400 // groupStructure operations
401 inline void outputOp(groupStruct vec){
402   if(vec.size()>50){
403     CkPrintf("vector too large to output!\n");
404     return;
405   }
406   CkPrintf("output vector: size=%d  {",vec.size());
407   for(int i=0;i<vec.size();i++)
408     CkPrintf(" %d ",vec[i]);
409   CkPrintf("}\n");
410 }
411 inline int getPosOp(int idx, groupStruct vec){
412   for (int r=0;r<vec.size();r++)
413     if (vec[r]==idx) return r;
414   return MPI_UNDEFINED;
415 }
416 inline groupStruct unionOp(groupStruct vec1, groupStruct vec2){
417   groupStruct newvec(vec1);
418   for(int i=0;i<vec2.size();i++){
419     if(getPosOp(vec2[i],vec1)==MPI_UNDEFINED)
420       newvec.push_back(vec2[i]);
421   }
422   return newvec;
423 }
424 inline groupStruct intersectOp(groupStruct vec1, groupStruct vec2){
425   groupStruct newvec;
426   for(int i=0;i<vec1.size();i++){
427     if(getPosOp(vec1[i],vec2)!=MPI_UNDEFINED)
428       newvec.push_back(vec1[i]);
429   }
430   return newvec;
431 }
432 inline groupStruct diffOp(groupStruct vec1, groupStruct vec2){
433   groupStruct newvec;
434   for(int i=0;i<vec1.size();i++){
435     if(getPosOp(vec1[i],vec2)==MPI_UNDEFINED)
436       newvec.push_back(vec1[i]);
437   }
438   return newvec;
439 }
440 inline int* translateRanksOp(int n,groupStruct vec1,int* ranks1,groupStruct vec2){
441   int* ret = new int[n];
442   for(int i=0;i<n;i++){
443     ret[i] = getPosOp(vec1[ranks1[i]],vec2);
444   }
445   return ret;
446 }
447 inline int compareVecOp(groupStruct vec1,groupStruct vec2){
448   int i,pos,ret = MPI_IDENT;
449   if(vec1.size() != vec2.size()) return MPI_UNEQUAL;
450   for(i=0;i<vec1.size();i++){
451     pos = getPosOp(vec1[i],vec2);
452     if(pos == MPI_UNDEFINED) return MPI_UNEQUAL;
453     if(pos != i)   ret = MPI_SIMILAR;
454   }
455   return ret;
456 }
457 inline groupStruct inclOp(int n,int* ranks,groupStruct vec){
458   groupStruct retvec;
459   for(int i=0;i<n;i++){
460     retvec.push_back(vec[ranks[i]]);
461   }
462   return retvec;
463 }
464 inline groupStruct exclOp(int n,int* ranks,groupStruct vec){
465   groupStruct retvec;
466   int add=1;
467   for(int j=0;j<vec.size();j++){
468     for(int i=0;i<n;i++)
469       if(j==ranks[i]){ add=0; break; }
470     if(add==1)  retvec.push_back(vec[j]);
471     else add=1;
472   }
473   return retvec;
474 }
475 inline groupStruct rangeInclOp(int n, int ranges[][3], groupStruct vec){
476   groupStruct retvec;
477   int first,last,stride;
478   for(int i=0;i<n;i++){
479     first = ranges[i][0];
480     last = ranges[i][1];
481     stride = ranges[i][2];
482     for(int j=0;j<=(last-first)/stride;j++)
483       retvec.push_back(vec[first+stride*j]);
484   }
485   return retvec;
486 }
487 inline groupStruct rangeExclOp(int n, int ranges[][3], groupStruct vec){
488   groupStruct retvec;
489   CkVec<int> ranksvec;
490   int first,last,stride;
491   int *ranks,cnt;
492   int i,j;
493   for(i=0;i<n;i++){
494     first = ranges[i][0];
495     last = ranges[i][1];
496     stride = ranges[i][2];
497     for(j=0;j<=(last-first)/stride;j++)
498       ranksvec.push_back(first+stride*j);
499   }
500   cnt=ranksvec.size();
501   ranks=new int[cnt];
502   for(i=0;i<cnt;i++)
503     ranks[i]=ranksvec[i];
504   return exclOp(cnt,ranks,vec);
505 }
506
507 #include "tcharm.h"
508 #include "tcharmc.h"
509
510 #include "ampi.decl.h"
511 #include "ddt.h"
512 #include "charm-api.h"
513 #include <sys/stat.h> // for mkdir
514
515 extern int _mpi_nworlds;
516
517 #define MPI_BCAST_TAG   MPI_TAG_UB_VALUE+10
518 #define MPI_BARR_TAG    MPI_TAG_UB_VALUE+11
519 #define MPI_REDUCE_TAG  MPI_TAG_UB_VALUE+12
520 #define MPI_GATHER_TAG  MPI_TAG_UB_VALUE+13
521 #define MPI_SCATTER_TAG MPI_TAG_UB_VALUE+14
522 #define MPI_SCAN_TAG    MPI_TAG_UB_VALUE+15
523 #define MPI_ATA_TAG     MPI_TAG_UB_VALUE+16
524
525 #if 0
526 // This is currently not used.
527 class BlockMap : public CkArrayMap {
528  public:
529   BlockMap(void) {}
530   BlockMap(CkMigrateMessage *m) {}
531   int registerArray(CkArrayMapRegisterMessage *m) {
532     delete m;
533     return 0;
534   }
535   int procNum(int /*arrayHdl*/,const CkArrayIndex &idx) {
536     int elem=*(int *)idx.data();
537     int penum =  (elem/(32/CkNumPes()));
538     CkPrintf("%d mapped to %d proc\n", elem, penum);
539     return penum;
540   }
541 };
542 #endif
543
544 #define MyAlign8(x) (((x)+7)&(~7))
545
546 /**
547 Represents an MPI request that has been initiated
548 using Isend, Irecv, Ialltoall, Send_init, etc.
549 */
550 class AmpiRequest {
551 public:
552         void *buf;
553         int count;
554         int type;
555         int tag;            // the order must match MPI_Status
556         int src;
557         int comm;
558
559 #if CMK_BLUEGENE_CHARM
560 public:
561         void *event;    // the event point that corresponding to this message
562 #endif
563 protected:
564         bool isvalid;
565 public:
566         AmpiRequest(){ }
567         /// Close this request (used by free and cancel)
568         virtual ~AmpiRequest(){ }
569
570         /// Activate this persistent request.
571         ///  Only meaningful for persistent requests,
572         ///  other requests just abort.
573         virtual int start(void){ return -1; }
574
575         /// Return true if this request is finished (progress).
576         virtual CmiBool test(MPI_Status *sts) =0;
577
578         /// Completes the operation hanging on the request
579         virtual void complete(MPI_Status *sts) =0;
580
581         /// Block until this request is finished,
582         ///  returning a valid MPI error code.
583         virtual int wait(MPI_Status *sts) =0;
584
585         virtual void receive(ampi *ptr, AmpiMsg *msg) = 0; 
586
587         /// Frees up the request: invalidate it
588         virtual void free(void){ isvalid=false; }
589         inline bool isValid(void){ return isvalid; }
590
591         /// Returns the type of request: 1-PersReq, 2-IReq, 3-ATAReq,
592         /// 4-SReq, 5-GPUReq
593         virtual int getType(void) =0;
594
595         virtual void pup(PUP::er &p) {
596                 p((char *)&buf,sizeof(void *));  //supposed to work only with isomalloc
597                 p(count);
598                 p(type);
599                 p(src);
600                 p(tag);
601                 p(comm);
602                 p(isvalid);
603 #if CMK_BLUEGENE_CHARM
604                 //needed for bigsim out-of-core emulation
605                 //as the "log" is not moved from memory, this pointer is safe
606                 //to be reused
607                 p((char *)&event, sizeof(void *));
608 #endif
609         }
610         
611         //added due to BIGSIM_OOC DEBUGGING
612         virtual void print();
613 };
614
615 class PersReq : public AmpiRequest {
616         int sndrcv; // 1 if send , 2 if recv, 3 if ssend
617 public:
618         PersReq(void *buf_, int count_, int type_, int src_, int tag_, 
619                 MPI_Comm comm_, int sndrcv_) 
620         {
621                 buf=buf_;  count=count_;  type=type_;  src=src_;  tag=tag_; 
622                 comm=comm_;  sndrcv=sndrcv_;  isvalid=true; 
623         }
624         PersReq(){};
625         ~PersReq(){ }
626         int start();
627         CmiBool test(MPI_Status *sts);
628         void complete(MPI_Status *sts);
629         int wait(MPI_Status *sts);
630         void receive(ampi *ptr, AmpiMsg *msg) {}
631         inline int getType(void){ return 1; }
632         virtual void pup(PUP::er &p){
633                 AmpiRequest::pup(p);
634                 p(sndrcv);
635         }
636         //added due to BIGSIM_OOC DEBUGGING
637         virtual void print();
638 };
639
640 class IReq : public AmpiRequest {
641 public:
642         bool statusIreq;
643         int length;     // recv'ed length
644         IReq(void *buf_, int count_, int type_, int src_, int tag_, MPI_Comm comm_)
645         {
646                 buf=buf_;  count=count_;  type=type_;  src=src_;  tag=tag_; 
647                 comm=comm_;  isvalid=true; statusIreq=false; length=0;
648         }
649         IReq(): statusIreq(false){};
650         ~IReq(){ }
651         CmiBool test(MPI_Status *sts);
652         void complete(MPI_Status *sts);
653         int wait(MPI_Status *sts);
654         inline int getType(void){ return 2; }
655         void receive(ampi *ptr, AmpiMsg *msg); 
656         virtual void pup(PUP::er &p){
657                 AmpiRequest::pup(p);
658                 p|statusIreq;  p|length;
659         }
660         //added due to BIGSIM_OOC DEBUGGING
661         virtual void print();
662 };
663
664 class ATAReq : public AmpiRequest {
665         class Request {
666         protected:
667                 void *buf;
668                 int count;
669                 int type;
670                 int src;
671                 int tag;
672                 int comm;
673 #if CMK_BLUEGENE_CHARM
674                 void *event;             // event buffered for the request
675 #endif
676                 virtual void pup(PUP::er &p){
677                         p((char *)&buf,sizeof(void *));  //supposed to work only with isomalloc
678                         p(count);
679                         p(type);
680                         p(src);p(tag);p(comm);
681 #if CMK_BLUEGENE_CHARM
682                 //needed for bigsim out-of-core emulation
683                 //as the "log" is not moved from memory, this pointer is safe
684                 //to be reused
685                         p((char *)&event, sizeof(void *));
686 #endif
687                 }
688         friend class ATAReq;
689         };
690         Request *myreqs;
691         int elmcount;
692         int idx;
693 public:
694         ATAReq(int c_):elmcount(c_),idx(0) { myreqs = new Request [c_]; isvalid=true; }
695         ATAReq(){};
696         ~ATAReq(void) { if(myreqs) delete [] myreqs; }
697         int addReq(void *buf_, int count_, int type_, int src_, int tag_, MPI_Comm comm_){
698                 myreqs[idx].buf=buf_;   myreqs[idx].count=count_;
699                 myreqs[idx].type=type_; myreqs[idx].src=src_;
700                 myreqs[idx].tag=tag_;   myreqs[idx].comm=comm_;
701                 return (++idx);
702         }
703         CmiBool test(MPI_Status *sts);
704         void complete(MPI_Status *sts);
705         int wait(MPI_Status *sts);
706         void receive(ampi *ptr, AmpiMsg *msg) {}
707         inline int getCount(void){ return elmcount; }
708         inline int getType(void){ return 3; }
709 //      inline void free(void){ isvalid=false; delete [] myreqs; }
710         virtual void pup(PUP::er &p){
711                 AmpiRequest::pup(p);
712                 p(elmcount);
713                 p(idx);
714                 if(p.isUnpacking()){
715                         myreqs = new Request[elmcount];
716                 }
717                 for(int i=0;i<idx;i++){
718                         myreqs[i].pup(p);
719                 }
720                 if(p.isDeleting()){
721                         delete []myreqs;
722                 }
723         }
724         //added due to BIGSIM_OOC DEBUGGING
725         virtual void print();
726 };
727
728 class SReq : public AmpiRequest {
729 public:
730         bool statusIreq;
731         SReq(MPI_Comm comm_): statusIreq(false) {
732                 comm = comm_; isvalid=true;
733         }
734         SReq(): statusIreq(false) {}
735         ~SReq(){ }
736         CmiBool test(MPI_Status *sts);
737         void complete(MPI_Status *sts);
738         int wait(MPI_Status *sts);
739         void receive(ampi *ptr, AmpiMsg *msg) {}
740         inline int getType(void){ return 4; }
741         virtual void pup(PUP::er &p){
742                 AmpiRequest::pup(p);
743                 p|statusIreq;
744         }
745         //added due to BIGSIM_OOC DEBUGGING
746         virtual void print();
747 };
748
749 class GPUReq : public AmpiRequest {
750     bool isComplete;
751
752 public:
753     GPUReq();
754     int getType() { return 5; }
755     CmiBool test(MPI_Status *sts);
756     void complete(MPI_Status *sts);
757     int wait(MPI_Status *sts);
758     void receive(ampi *ptr, AmpiMsg *msg);
759     void setComplete();
760 };
761
762 /// Special CkVec<AmpiRequest*> for AMPI. Most code copied from cklist.h
763 class AmpiRequestList : private CkSTLHelper<AmpiRequest *> {
764     AmpiRequest** block; //Elements of vector
765     int blklen; //Allocated size of block
766     int len; //Number of used elements in block
767     void makeBlock(int blklen_,int len_) {
768        block=new AmpiRequest* [blklen_];
769        blklen=blklen_; len=len_;
770     }
771     void freeBlock(void) {
772        len=0; blklen=0;
773        delete[] block; block=NULL;
774     }
775     void copyFrom(const AmpiRequestList &src) {
776        makeBlock(src.blklen, src.len);
777        elementCopy(block,src.block,blklen);
778     }
779   public:
780     AmpiRequestList() {block=NULL;blklen=len=0;}
781     ~AmpiRequestList() { freeBlock(); }
782     AmpiRequestList(const AmpiRequestList &src) {copyFrom(src);}
783     AmpiRequestList(int size) { makeBlock(size,size); }
784     AmpiRequestList &operator=(const AmpiRequestList &src) {
785       freeBlock();
786       copyFrom(src);
787       return *this;
788     }
789
790     AmpiRequest* operator[](size_t n) { return block[n]; }
791
792     int size(void) const {return len;}
793     void setSize(int blklen_) {
794       AmpiRequest **oldBlock=block;
795       makeBlock(blklen_,len);
796       elementCopy(block,oldBlock,len);
797       delete[] oldBlock; //WARNING: leaks if element copy throws exception
798     }
799     //Grow to contain at least this position:
800     void growAtLeast(int pos) {
801       if (pos>=blklen) setSize(pos*2+16);
802     }
803     void insertAt(int pos, AmpiRequest* elt) {
804       if (pos>=len) {
805         growAtLeast(pos);
806         len=pos+1;
807       }
808       block[pos] = elt;
809     }
810     void free(int pos) {
811       if (pos<0 || pos>=len) return;
812       block[pos]->free();
813       delete block[pos];
814       block[pos]=NULL;
815     }
816     void push_back(AmpiRequest* elt) {insertAt(len,elt);}
817     int insert(AmpiRequest* elt){
818       //search for invalidated slot
819       // disabled to make requests monotonously ascending 
820       // for multiple completion calls like MPI_Waitany
821        for(int i=0;i<len;i++){
822         if(block[i]==NULL){
823           block[i] = elt;
824           return i;
825         }
826       } 
827       push_back(elt);
828       return len-1;
829     }
830
831     inline void checkRequest(MPI_Request idx){
832       if(!(idx==-1 || (idx < this->len && (block[idx])->isValid())))
833         CkAbort("Invalid MPI_Request\n");
834     }
835
836     //find an AmpiRequest by its pointer value
837     //return -1 if not found!
838     int findRequestIndex(AmpiRequest *req){
839         for(int i=0; i<len; i++){
840             if(block[i]==req) return i;
841         }
842         return -1;
843     }
844
845     void pup(PUP::er &p);
846     
847     //BIGSIM_OOC DEBUGGING
848     void print(){
849         for(int i=0; i<len; i++){
850             if(block[i]==NULL) continue;
851             CmiPrintf("AmpiRequestList Element %d [%p]: \n", i+1, block[i]);
852             block[i]->print();
853         }
854     }
855 };
856
857 //A simple memory buffer
858 class memBuf {
859         CkVec<char> buf;
860 public:
861         memBuf() { }
862         memBuf(int size) :buf(size) {}
863         void setSize(int s) {buf.resize(s);}
864         int getSize(void) const {return buf.size();}
865         const void *getData(void) const {return (const void *)&buf[0];}
866         void *getData(void) {return (void *)&buf[0];}
867 };
868
869 template <class T>
870 inline void pupIntoBuf(memBuf &b,T &t) {
871         PUP::sizer ps;ps|t;
872         b.setSize(ps.size());
873         PUP::toMem pm(b.getData()); pm|t;
874 }
875
876 template <class T>
877 inline void pupFromBuf(const void *data,T &t) {
878         PUP::fromMem p(data); p|t;
879 }
880
881 class AmpiMsg : public CMessage_AmpiMsg {
882  public:
883   int seq; //Sequence number (for message ordering)
884   int tag; //MPI tag
885   int srcIdx; //Array index of source
886   int srcRank; //Communicator rank for source
887   MPI_Comm comm; //Communicator for source
888   int length; //Number of bytes in this message 
889 #if CMK_BLUEGENE_CHARM
890   void *event;
891   int  eventPe;  // the PE that the event is located
892 #endif
893   char *data;
894
895   AmpiMsg(void) { data = NULL; }
896   AmpiMsg(int _s, int t, int sIdx,int sRank, int l, int c) :
897     seq(_s), tag(t),srcIdx(sIdx), srcRank(sRank), comm(c), length(l) {
898   }  
899   static AmpiMsg* pup(PUP::er &p, AmpiMsg *m)
900   {
901     int seq, length, tag, srcIdx, srcRank, comm;
902     if(p.isPacking() || p.isSizing()) {
903       seq = m->seq;
904       tag = m->tag;
905       srcIdx = m->srcIdx;
906       srcRank = m->srcRank;
907       comm = m->comm;
908       length = m->length;
909     }
910     p(seq); p(tag); p(srcIdx); p(srcRank); p(comm); p(length); 
911     if(p.isUnpacking()) {
912       m = new (length, 0) AmpiMsg(seq, tag, srcIdx, srcRank, length, comm);
913     }
914     p(m->data, length);
915     if(p.isDeleting()) {
916       delete m;
917       m = 0;
918     }
919     return m;
920   }
921 };
922
923
924 #if 0
925 // Moved to cklists.h
926 /**
927   An array that is broken up into "pages"
928   which are separately allocated.  Useful for saving memory
929   when storing a long array with few used entries.
930   Should be almost as fast as a regular vector.
931   
932   When nothing is used, this reduces storage space by a factor
933   of pageSize*sizeof(T)/sizeof(T*).  For example, an 
934   array of 64-bit ints of length 1M takes 8MB; with paging
935   it takes just 32KB, plus 2KB for each used block.
936   
937   The class T must have a pup routine, a default constructor,
938   and a copy constructor.
939 */
940 template <class T>
941 class CkPagedVector : private CkNoncopyable {
942         /// This is the number of array entries per page.
943         ///  For any kind of efficiency, it must be a power of 2,
944         ///  which lets the compiler turn divides and mods into
945         ///  bitwise operations.
946         enum {pageSize=256u};
947         T **pages;
948         int nPages;
949         void allocatePages(int nEntries) {
950                 nPages=(nEntries+pageSize-1)/pageSize; // round up
951                 pages=new T*[nPages];
952                 for (int pg=0;pg<nPages;pg++) pages[pg]=0;
953         }
954         void allocatePage(int pg) {
955                 T *np=new T[pageSize];
956                 // Initialize the new elements of the page to 0:
957                 for (int of=0;of<pageSize;of++) np[of]=T();
958                 pages[pg]=np;
959         }
960 public:
961         CkPagedVector() :pages(0), nPages(0) {}
962         CkPagedVector(int nEntries) {init(nEntries);}
963         ~CkPagedVector() {
964                 for (int pg=0;pg<nPages;pg++) 
965                         if (pages[pg]) 
966                                 delete[] pages[pg];
967                 delete[] pages;
968         }
969         void init(int nEntries) {
970                 allocatePages(nEntries);
971         }
972         
973         inline T &operator[] (unsigned int i) {
974                 // This divide and mod should be fast, 
975                 //  as long as pageSize is a power of 2:
976                 int pg=i/pageSize; 
977                 int of=i%pageSize;
978                 if (pages[pg]==NULL) allocatePage(pg);
979                 return pages[pg][of];
980         }
981         
982         void pupPage(PUP::er &p,int pg) {
983                 if (p.isUnpacking()) allocatePage(pg);
984                 for (int of=0;of<pageSize;of++) 
985                         p|pages[pg][of];
986         }
987         void pup(PUP::er &p) {
988                 int pg, o;
989                 unsigned int i;
990                 p|nPages;
991                 if (!p.isUnpacking()) 
992                 { // Packing phase: save each used page
993                         for (pg=0;pg<nPages;pg++) 
994                         if (pages[pg]) {
995                                 p|pg;
996                                 pupPage(p,pg);
997                         }
998                         pg=-1;
999                         p|pg;
1000                 } else 
1001                 { // Unpacking phase: allocate and restore
1002                         allocatePages(nPages*pageSize);
1003                         p|pg;
1004                         while (pg!=-1) {
1005                                 pupPage(p,pg);
1006                                 p|pg;
1007                         }
1008                 }
1009         }
1010 };
1011 #endif
1012
1013
1014 /**
1015   Our local representation of another AMPI
1016  array element.  Used to keep track of incoming
1017  and outgoing message sequence numbers, and 
1018  the out-of-order message list.
1019 */
1020 class AmpiOtherElement {
1021 public:
1022 /// Next incoming and outgoing message sequence number
1023   int seqIncoming, seqOutgoing;
1024   
1025 /// Number of elements in out-of-order queue. (normally 0)
1026   int nOut;
1027
1028   AmpiOtherElement(void) {
1029     seqIncoming=0; seqOutgoing=0;
1030     nOut=0;
1031   }
1032   
1033   void pup(PUP::er &p) {
1034     p|seqIncoming; p|seqOutgoing;
1035     p|nOut;
1036   }
1037 };
1038
1039 class AmpiSeqQ : private CkNoncopyable {
1040   CkMsgQ<AmpiMsg> out;        // all out of order messages
1041   CkPagedVector<AmpiOtherElement>  elements;   // element info
1042
1043   void putOutOfOrder(int srcIdx, AmpiMsg *msg);
1044   
1045 public:
1046   AmpiSeqQ() {}
1047   void init(int numP);
1048   ~AmpiSeqQ ();
1049   void pup(PUP::er &p);
1050   
1051   /// Insert this message in the table.  Returns the number 
1052   ///  of messages now available for the element.
1053   ///   If 0, the message was out-of-order and is buffered.
1054   ///   If 1, this message can be immediately processed.
1055   ///   If >1, this message can be immediately processed,
1056   ///    and you should call "getOutOfOrder" repeatedly.
1057   inline int put(int srcIdx, AmpiMsg *msg) {
1058      AmpiOtherElement &el=elements[srcIdx];
1059      if (msg->seq==el.seqIncoming) { // In order:
1060        el.seqIncoming++;
1061        return 1+el.nOut;
1062      }
1063      else { // Out of order: stash message
1064        putOutOfOrder(srcIdx, msg);
1065        return 0;
1066      }
1067   }
1068   
1069   /// Get an out-of-order message from the table.
1070   ///  (in-order messages never go into the table)
1071   AmpiMsg *getOutOfOrder(int p);
1072   
1073   /// Return the next outgoing sequence number, and increment it.
1074   int nextOutgoing(int p) {
1075      return elements[p].seqOutgoing++;
1076   }
1077 };
1078 PUPmarshall(AmpiSeqQ)
1079
1080
1081 inline CProxy_ampi ampiCommStruct::getProxy(void) const {return ampiID;}
1082 const ampiCommStruct &universeComm2CommStruct(MPI_Comm universeNo);
1083
1084 /* KeyValue class for caching */
1085 class KeyvalNode {
1086 public:
1087         MPI_Copy_function *copy_fn;
1088         MPI_Delete_function *delete_fn;
1089         void *extra_state;
1090         /* value is associated with getKeyvals of communicator */
1091
1092         KeyvalNode(void): copy_fn(NULL), delete_fn(NULL), extra_state(NULL)
1093         { }
1094         KeyvalNode(MPI_Copy_function *cf, MPI_Delete_function *df, void* es):
1095                 copy_fn(cf), delete_fn(df), extra_state(es)
1096         { }
1097         // KeyvalNode is not supposed to be pup'ed
1098         void pup(PUP::er& p){ /* empty */ }
1099 };
1100
1101 /*
1102 An ampiParent holds all the communicators and the TCharm thread
1103 for its children, which are bound to it.
1104 */
1105 class ampiParent : public CBase_ampiParent {
1106     CProxy_TCharm threads;
1107     TCharm *thread;
1108     void prepareCtv(void);
1109
1110     MPI_Comm worldNo; //My MPI_COMM_WORLD
1111     ampi *worldPtr; //AMPI element corresponding to MPI_COMM_WORLD
1112     ampiCommStruct worldStruct;
1113     ampiCommStruct selfStruct;
1114
1115     CkPupPtrVec<ampiCommStruct> splitComm; //Communicators from MPI_Comm_split
1116     CkPupPtrVec<ampiCommStruct> groupComm; //Communicators from MPI_Comm_group
1117     CkPupPtrVec<ampiCommStruct> cartComm;  //Communicators from MPI_Cart_create
1118     CkPupPtrVec<ampiCommStruct> graphComm; //Communicators from MPI_Graph_create
1119     CkPupPtrVec<ampiCommStruct> interComm; //Communicators from MPI_Intercomm_create
1120     CkPupPtrVec<ampiCommStruct> intraComm; //Communicators from MPI_Intercomm_merge
1121
1122     CkPupPtrVec<groupStruct> groups; // "Wild" groups that don't have a communicator
1123     CkPupPtrVec<WinStruct> winStructList;   //List of windows for one-sided communication
1124     CkPupPtrVec<InfoStruct> infos; // list of all MPI_Infos
1125
1126     inline int isSplit(MPI_Comm comm) const {
1127       return (comm>=MPI_COMM_FIRST_SPLIT && comm<MPI_COMM_FIRST_GROUP);
1128     }
1129     const ampiCommStruct &getSplit(MPI_Comm comm) {
1130       int idx=comm-MPI_COMM_FIRST_SPLIT;
1131       if (idx>=splitComm.size()) CkAbort("Bad split communicator used");
1132       return *splitComm[idx];
1133     }
1134     void splitChildRegister(const ampiCommStruct &s);
1135
1136     inline int isGroup(MPI_Comm comm) const {
1137       return (comm>=MPI_COMM_FIRST_GROUP && comm<MPI_COMM_FIRST_CART);
1138     }
1139     const ampiCommStruct &getGroup(MPI_Comm comm) {
1140       int idx=comm-MPI_COMM_FIRST_GROUP;
1141       if (idx>=groupComm.size()) CkAbort("Bad group communicator used");
1142       return *groupComm[idx];
1143     }
1144     void groupChildRegister(const ampiCommStruct &s);
1145     inline int isInGroups(MPI_Group group) const {
1146       return (group>=0 && group<groups.size());
1147     }
1148
1149     void cartChildRegister(const ampiCommStruct &s);
1150     void graphChildRegister(const ampiCommStruct &s);
1151     void interChildRegister(const ampiCommStruct &s);
1152
1153     inline int isIntra(MPI_Comm comm) const {
1154       return (comm>=MPI_COMM_FIRST_INTRA && comm<MPI_COMM_FIRST_RESVD);
1155     }
1156     const ampiCommStruct &getIntra(MPI_Comm comm) {
1157       int idx=comm-MPI_COMM_FIRST_INTRA;
1158       if (idx>=intraComm.size()) CkAbort("Bad intra-communicator used");
1159       return *intraComm[idx];
1160     }
1161     void intraChildRegister(const ampiCommStruct &s);
1162
1163     // MPI MPI_Attr_get C binding returns a *pointer* to an integer,
1164     //  so there needs to be some storage somewhere to point to.
1165     int* kv_builtin_storage;
1166     int kv_is_builtin(int keyval);
1167     CkPupPtrVec<KeyvalNode> kvlist;
1168
1169     int RProxyCnt;
1170     CProxy_ampi tmpRProxy;
1171
1172 public:
1173     int ampiInitCallDone;
1174
1175 public:
1176     ampiParent(MPI_Comm worldNo_,CProxy_TCharm threads_);
1177     ampiParent(CkMigrateMessage *msg);
1178     void ckJustMigrated(void);
1179     void ckJustRestored(void);
1180     ~ampiParent();
1181
1182     ampi *lookupComm(MPI_Comm comm) {
1183        if (comm!=worldStruct.getComm())
1184           CkAbort("ampiParent::lookupComm> Bad communicator!");
1185        return worldPtr;
1186     }
1187
1188     //Children call this when they are first created, or just migrated
1189     TCharm *registerAmpi(ampi *ptr,ampiCommStruct s,bool forMigration);
1190
1191     // exchange proxy info between two ampi proxies
1192     void ExchangeProxy(CProxy_ampi rproxy){
1193       if(RProxyCnt==0){ tmpRProxy=rproxy; RProxyCnt=1; }
1194       else if(RProxyCnt==1) { tmpRProxy.setRemoteProxy(rproxy); rproxy.setRemoteProxy(tmpRProxy); RProxyCnt=0; }
1195       else CkAbort("ExchangeProxy: RProxyCnt>1");
1196     }
1197
1198     //Grab the next available split/group communicator
1199     MPI_Comm getNextSplit(void) const {return MPI_COMM_FIRST_SPLIT+splitComm.size();}
1200     MPI_Comm getNextGroup(void) const {return MPI_COMM_FIRST_GROUP+groupComm.size();}
1201     MPI_Comm getNextCart(void) const {return MPI_COMM_FIRST_CART+cartComm.size();}
1202     MPI_Comm getNextGraph(void) const {return MPI_COMM_FIRST_GRAPH+graphComm.size();}
1203     MPI_Comm getNextInter(void) const {return MPI_COMM_FIRST_INTER+interComm.size();}
1204     MPI_Comm getNextIntra(void) const {return MPI_COMM_FIRST_INTRA+intraComm.size();}
1205
1206     inline int isCart(MPI_Comm comm) const {
1207       return (comm>=MPI_COMM_FIRST_CART && comm<MPI_COMM_FIRST_GRAPH);
1208     }
1209     ampiCommStruct &getCart(MPI_Comm comm) {
1210       int idx=comm-MPI_COMM_FIRST_CART;
1211       if (idx>=cartComm.size()) CkAbort("Bad cartesian communicator used");
1212       return *cartComm[idx];
1213     }
1214     inline int isGraph(MPI_Comm comm) const {
1215       return (comm>=MPI_COMM_FIRST_GRAPH && comm<MPI_COMM_FIRST_INTER);
1216     }
1217     ampiCommStruct &getGraph(MPI_Comm comm) {
1218       int idx=comm-MPI_COMM_FIRST_GRAPH;
1219       if (idx>=graphComm.size()) CkAbort("Bad graph communicator used");
1220       return *graphComm[idx];
1221     }
1222     inline int isInter(MPI_Comm comm) const {
1223       return (comm>=MPI_COMM_FIRST_INTER && comm<MPI_COMM_FIRST_INTRA);
1224     }
1225     const ampiCommStruct &getInter(MPI_Comm comm) {
1226       int idx=comm-MPI_COMM_FIRST_INTER;
1227       if (idx>=interComm.size()) CkAbort("Bad inter-communicator used");
1228       return *interComm[idx];
1229     }
1230
1231     void pup(PUP::er &p);
1232
1233     inline void start_measure() {
1234       usesAutoMeasure = CmiFalse;
1235     }
1236     inline void stop_measure() {
1237       usesAutoMeasure = CmiTrue;
1238     }
1239     virtual void UserSetLBLoad(void) {
1240       // empty
1241     }
1242
1243     void startCheckpoint(const char* dname);
1244     void Checkpoint(int len, const char* dname);
1245     void ResumeThread(void);
1246     TCharm* getTCharmThread() {return thread;}
1247
1248     inline const ampiCommStruct &comm2CommStruct(MPI_Comm comm) {
1249       if (comm==MPI_COMM_WORLD) return worldStruct;
1250       if (comm==MPI_COMM_SELF) return selfStruct;
1251       if (comm==worldNo) return worldStruct;
1252       if (isSplit(comm)) return getSplit(comm);
1253       if (isGroup(comm)) return getGroup(comm);
1254       if (isCart(comm)) return getCart(comm);
1255       if (isGraph(comm)) return getGraph(comm);
1256       if (isInter(comm)) return getInter(comm);
1257       if (isIntra(comm)) return getIntra(comm);
1258       return universeComm2CommStruct(comm);
1259     }
1260     
1261      //ampi *comm2ampi(MPI_Comm comm);
1262     inline ampi *comm2ampi(MPI_Comm comm) {
1263       if (comm==MPI_COMM_WORLD) return worldPtr;
1264       if (comm==MPI_COMM_SELF) return worldPtr;
1265       if (comm==worldNo) return worldPtr;
1266       if (isSplit(comm)) {
1267          const ampiCommStruct &st=getSplit(comm);
1268          return st.getProxy()[thisIndex].ckLocal();
1269       }
1270       if (isGroup(comm)) {
1271          const ampiCommStruct &st=getGroup(comm);
1272          return st.getProxy()[thisIndex].ckLocal();
1273       }
1274       if (isCart(comm)) {
1275         const ampiCommStruct &st = getCart(comm);
1276         return st.getProxy()[thisIndex].ckLocal();
1277       }
1278       if (isGraph(comm)) {
1279         const ampiCommStruct &st = getGraph(comm);
1280         return st.getProxy()[thisIndex].ckLocal();
1281       }
1282       if (isInter(comm)) {
1283          const ampiCommStruct &st=getInter(comm);
1284          return st.getProxy()[thisIndex].ckLocal();
1285       }
1286       if (isIntra(comm)) {
1287          const ampiCommStruct &st=getIntra(comm);
1288          return st.getProxy()[thisIndex].ckLocal();
1289       }
1290       if (comm>MPI_COMM_WORLD) return worldPtr; //Use MPI_WORLD ampi for cross-world messages:
1291       CkAbort("Invalid communicator used!");
1292       return NULL;
1293     }
1294
1295     inline int hasComm(const MPI_Group group){
1296       MPI_Comm comm = (MPI_Comm)group;
1297       return ( comm==MPI_COMM_WORLD || comm==worldNo || isSplit(comm) || isGroup(comm) || isCart(comm) || isGraph(comm) || isIntra(comm) ); //isInter omitted because its comm number != its group number
1298     }
1299     inline const groupStruct group2vec(MPI_Group group){
1300       if(hasComm(group))
1301         return comm2CommStruct((MPI_Comm)group).getIndices();
1302       if(isInGroups(group))
1303         return *groups[group];
1304       CkAbort("ampiParent::group2vec: Invalid group id!");
1305       return *groups[0]; //meaningless return
1306     }
1307     inline MPI_Group saveGroupStruct(groupStruct vec){
1308       int idx = groups.size();
1309       groups.resize(idx+1);
1310       groups[idx]=new groupStruct(vec);
1311       return (MPI_Group)idx;
1312     }
1313     inline int getRank(const MPI_Group group){
1314       groupStruct vec = group2vec(group);
1315       return getPosOp(thisIndex,vec);
1316     }
1317     
1318     inline int getMyPe(void){
1319       return CkMyPe();
1320     }
1321     inline int hasWorld(void) const {
1322       return worldPtr!=NULL;
1323     }
1324     
1325     inline void checkComm(MPI_Comm comm){
1326       if ((comm > MPI_COMM_FIRST_RESVD && comm != MPI_COMM_SELF && comm != MPI_COMM_WORLD) 
1327        || (isSplit(comm) && comm-MPI_COMM_FIRST_SPLIT >= splitComm.size())
1328        || (isGroup(comm) && comm-MPI_COMM_FIRST_GROUP >= groupComm.size())
1329        || (isCart(comm)  && comm-MPI_COMM_FIRST_CART  >=  cartComm.size())
1330        || (isGraph(comm) && comm-MPI_COMM_FIRST_GRAPH >= graphComm.size())
1331        || (isInter(comm) && comm-MPI_COMM_FIRST_INTER >= interComm.size())
1332        || (isIntra(comm) && comm-MPI_COMM_FIRST_INTRA >= intraComm.size()) )
1333           CkAbort("Invalide MPI_Comm\n");
1334     }
1335
1336     /// if intra-communicator, return comm, otherwise return null group
1337     inline MPI_Group comm2group(const MPI_Comm comm){
1338       if(isInter(comm)) return MPI_GROUP_NULL;   // we don't support inter-communicator in such functions
1339       ampiCommStruct s = comm2CommStruct(comm);
1340       if(comm!=MPI_COMM_WORLD && comm!=s.getComm()) CkAbort("Error in ampiParent::comm2group()");
1341       return (MPI_Group)(s.getComm());
1342     }
1343
1344     inline int getRemoteSize(const MPI_Comm comm){
1345       if(isInter(comm)) return getInter(comm).getRemoteIndices().size();
1346       else return -1;
1347     }
1348     inline MPI_Group getRemoteGroup(const MPI_Comm comm){
1349       if(isInter(comm)) return saveGroupStruct(getInter(comm).getRemoteIndices());
1350       else return MPI_GROUP_NULL;
1351     }
1352
1353     int createKeyval(MPI_Copy_function *copy_fn, MPI_Delete_function *delete_fn,
1354                      int *keyval, void* extra_state);
1355     int freeKeyval(int *keyval);
1356     int putAttr(MPI_Comm comm, int keyval, void* attribute_val);
1357     int getAttr(MPI_Comm comm, int keyval, void *attribute_val, int *flag);
1358     int deleteAttr(MPI_Comm comm, int keyval);
1359
1360     CkDDT myDDTsto;
1361     CkDDT *myDDT;
1362     AmpiRequestList ampiReqs;
1363
1364     //added to make sure post_ireqs in ampi class share the same pointers
1365     //with those in ampiReqs after pupping routines.
1366     //AmpiRequestList oldAmpiReqs;
1367     
1368     int addWinStruct(WinStruct* win);
1369     WinStruct getWinStruct(MPI_Win win);
1370     void removeWinStruct(WinStruct win);
1371
1372 #if AMPI_COUNTER
1373 public:
1374     AmpiCounters counters;
1375 #endif    
1376     
1377 public:
1378     MPI_Info createInfo(void);
1379     MPI_Info dupInfo(MPI_Info info);
1380     void setInfo(MPI_Info info, char *key, char *value);
1381     int deleteInfo(MPI_Info info, char *key);    
1382     int getInfo(MPI_Info info, char *key, int valuelen, char *value);
1383     int getInfoValuelen(MPI_Info info, char *key, int *valuelen);
1384     int getInfoNkeys(MPI_Info info);
1385     int getInfoNthkey(MPI_Info info, int n, char *key);
1386     void freeInfo(MPI_Info info);
1387
1388  public:
1389 #ifdef AMPIMSGLOG
1390     /* message logging */
1391     int pupBytes;
1392 #if CMK_PROJECTIONS_USE_ZLIB
1393     gzFile fMsgLog;
1394     PUP::tozDisk *toPUPer;
1395     PUP::fromzDisk *fromPUPer;
1396 #else
1397     FILE* fMsgLog;
1398     PUP::tozDisk *toPUPer;
1399     PUP::fromzDisk *fromPUPer;
1400 #endif
1401 #endif
1402     void init();
1403     void finalize();
1404 };
1405
1406 /*
1407 An ampi manages the communication of one thread over
1408 one MPI communicator.
1409 */
1410 class ampi : public CBase_ampi {
1411 friend class IReq;
1412 friend class SReq;
1413     CProxy_ampiParent parentProxy;
1414     void findParent(bool forMigration);
1415     ampiParent *parent;
1416     TCharm *thread;
1417     bool resumeOnRecv;
1418
1419     ampiCommStruct myComm;
1420     int myRank;
1421     groupStruct tmpVec; // stores temp group info
1422     CProxy_ampi remoteProxy; // valid only for intercommunicator
1423
1424     /// A proxy used when delegating message sends to comlib
1425     CProxy_ampi comlibProxy;
1426     
1427     /// References to the comlib instance handles(currently just integers)
1428     ComlibInstanceHandle ciStreaming;
1429     ComlibInstanceHandle ciBcast;
1430     ComlibInstanceHandle ciAllgather;
1431     ComlibInstanceHandle ciAlltoall;
1432
1433     
1434     int seqEntries; //Number of elements in below arrays
1435     AmpiSeqQ oorder;
1436     void inorder(AmpiMsg *msg);
1437
1438     void init(void);
1439
1440   public: // entry methods
1441
1442     ampi();
1443     ampi(CkArrayID parent_,const ampiCommStruct &s);
1444     ampi(CkArrayID parent_,const ampiCommStruct &s,ComlibInstanceHandle ciStreaming_,
1445         ComlibInstanceHandle ciBcast_,ComlibInstanceHandle ciAllgather_,ComlibInstanceHandle ciAlltoall_);
1446     ampi(CkMigrateMessage *msg);
1447     void ckJustMigrated(void);
1448     void ckJustRestored(void);
1449     ~ampi();
1450
1451     virtual void pup(PUP::er &p);
1452
1453     void allInitDone(CkReductionMsg *m);
1454     void setInitDoneFlag();
1455   
1456     void block(void);
1457     void unblock(void);
1458     void yield(void);
1459     void generic(AmpiMsg *);
1460     void ssend_ack(int sreq);
1461     void reduceResult(CkReductionMsg *m);
1462     void splitPhase1(CkReductionMsg *msg);
1463     void commCreatePhase1(CkReductionMsg *msg);
1464     void cartCreatePhase1(CkReductionMsg *m);
1465     void graphCreatePhase1(CkReductionMsg *m);
1466     void intercommCreatePhase1(CkReductionMsg *m);
1467     void intercommMergePhase1(CkReductionMsg *msg);
1468
1469   public: // to be used by MPI_* functions
1470
1471     inline const ampiCommStruct &comm2CommStruct(MPI_Comm comm) {
1472       return parent->comm2CommStruct(comm);
1473     }
1474
1475     AmpiMsg *makeAmpiMsg(int destIdx,int t,int sRank,const void *buf,int count,
1476                          int type,MPI_Comm destcomm, int sync=0);
1477
1478     inline void comlibsend(int t, int s, const void* buf, int count, int type, int rank, MPI_Comm destcomm);
1479     inline void send(int t, int s, const void* buf, int count, int type, int rank, MPI_Comm destcomm, int sync=0);
1480     static void sendraw(int t, int s, void* buf, int len, CkArrayID aid,
1481                         int idx);
1482     void delesend(int t, int s, const void* buf, int count, int type,  
1483                   int rank, MPI_Comm destcomm, CProxy_ampi arrproxy, int sync=0);
1484     inline int processMessage(AmpiMsg *msg, int t, int s, void* buf, int count, int type);
1485     inline AmpiMsg * getMessage(int t, int s, int comm, int *sts);
1486     int recv(int t,int s,void* buf,int count,int type,int comm,int *sts=0);
1487     void probe(int t,int s,int comm,int *sts);
1488     int iprobe(int t,int s,int comm,int *sts);
1489     void bcast(int root, void* buf, int count, int type,MPI_Comm comm);
1490     static void bcastraw(void* buf, int len, CkArrayID aid);
1491     void split(int color,int key,MPI_Comm *dest, int type);
1492     void commCreate(const groupStruct vec,MPI_Comm *newcomm);
1493     void cartCreate(const groupStruct vec, MPI_Comm *newcomm);
1494     void graphCreate(const groupStruct vec, MPI_Comm *newcomm);
1495     void intercommCreate(const groupStruct rvec, int root, MPI_Comm *ncomm);
1496
1497     inline int isInter(void) { return myComm.isinter(); }
1498     void intercommMerge(int first, MPI_Comm *ncomm);
1499
1500     inline int getWorldRank(void) const {return parent->thisIndex;}
1501     /// Return our rank in this communicator
1502     inline int getRank(MPI_Comm comm) const {
1503         if (comm==MPI_COMM_SELF) return 0;
1504         else return myRank;
1505     }
1506     inline int getSize(MPI_Comm comm) const {
1507         if (comm==MPI_COMM_SELF) return 1;
1508         else return myComm.getSize();
1509     }
1510     inline MPI_Comm getComm(void) const {return myComm.getComm();}
1511     inline CkVec<int> getIndices(void) const { return myComm.getindices(); }
1512     inline const CProxy_ampi &getProxy(void) const {return thisProxy;}
1513     inline const CProxy_ampi &getComlibProxy(void) const { return comlibProxy; }
1514     inline const CProxy_ampi &getRemoteProxy(void) const {return remoteProxy;}
1515     inline void setRemoteProxy(CProxy_ampi rproxy) { remoteProxy = rproxy; thread->resume(); }
1516     inline int getIndexForRank(int r) const {return myComm.getIndexForRank(r);}
1517     inline int getIndexForRemoteRank(int r) const {return myComm.getIndexForRemoteRank(r);}
1518     inline ComlibInstanceHandle getStreaming(void) { return ciStreaming; }
1519     inline ComlibInstanceHandle getBcast(void) { return ciBcast; }
1520     inline ComlibInstanceHandle getAllgather(void) { return ciAllgather; }
1521     inline ComlibInstanceHandle getAlltoall(void) { return ciAlltoall; }
1522
1523 #if AMPI_COMLIB
1524     inline Strategy* getStreamingStrategy(void) { return CkpvAccess(conv_com_object).getStrategy(ciStreaming); }
1525     inline Strategy* getBcastStrategy(void) { return CkpvAccess(conv_com_object).getStrategy(ciBcast); }
1526     inline Strategy* getAllgatherStrategy(void) { return CkpvAccess(conv_com_object).getStrategy(ciAllgather); }
1527     inline Strategy* getAlltoallStrategy(void) { return CkpvAccess(conv_com_object).getStrategy(ciAlltoall); }
1528 #endif
1529     
1530     CkDDT *getDDT(void) {return parent->myDDT;}
1531     CthThread getThread() { return thread->getThread(); }
1532 #if CMK_LBDB_ON
1533     void setMigratable(int mig) { 
1534       if(mig) thread->setMigratable(CmiTrue); 
1535       else thread->setMigratable(CmiFalse);
1536     }
1537 #endif
1538   public:
1539     //These are directly used by API routines, which is hideous
1540     /*
1541     FIXME: CmmTable is only indexed by the tag, sender, and communicator.
1542     It should also be indexed by the source data type and length (if any).
1543     */
1544     CmmTable msgs;
1545     CmmTable posted_ireqs;         // posted irecv req
1546     int nbcasts;
1547     //------------------------ Added by YAN ---------------------
1548  private:
1549     CkPupPtrVec<win_obj> winObjects;
1550  public:
1551     MPI_Win createWinInstance(void *base, MPI_Aint size, int disp_unit, MPI_Info info); 
1552     int deleteWinInstance(MPI_Win win);
1553     int winGetGroup(WinStruct win, MPI_Group *group); 
1554     int winPut(void *orgaddr, int orgcnt, MPI_Datatype orgtype, int rank, 
1555                MPI_Aint targdisp, int targcnt, MPI_Datatype targtype, WinStruct win);
1556     int winGet(void *orgaddr, int orgcnt, MPI_Datatype orgtype, int rank, 
1557                MPI_Aint targdisp, int targcnt, MPI_Datatype targtype, WinStruct win);
1558     int winIGet(MPI_Aint orgdisp, int orgcnt, MPI_Datatype orgtype, int rank,
1559                MPI_Aint targdisp, int targcnt, MPI_Datatype targtype, WinStruct win, 
1560                MPI_Request *req);
1561     int winIGetWait(MPI_Request *request, MPI_Status *status);
1562     int winIGetFree(MPI_Request *request, MPI_Status *status);
1563     void winRemotePut(int orgtotalsize, char* orgaddr, int orgcnt, MPI_Datatype orgtype,
1564                       MPI_Aint targdisp, int targcnt, MPI_Datatype targtype, 
1565                       int winIndex, CkFutureID ftHandle, int pe_src);
1566     void winRemoteGet(int orgcnt, MPI_Datatype orgtype, MPI_Aint targdisp, 
1567                       int targcnt, MPI_Datatype targtype, 
1568                       int winIndex, CkFutureID ftHandle, int pe_src);
1569     AmpiMsg* winRemoteIGet(int orgdisp, int orgcnt, MPI_Datatype orgtype, MPI_Aint targdisp,
1570                            int targcnt, MPI_Datatype targtype, int winIndex);
1571     int winLock(int lock_type, int rank, WinStruct win);
1572     int winUnlock(int rank, WinStruct win);
1573     void winRemoteLock(int lock_type, int winIndex, CkFutureID ftHandle, int pe_src, int requestRank);
1574     void winRemoteUnlock(int winIndex, CkFutureID ftHandle, int pe_src, int requestRank);
1575     int winAccumulate(void *orgaddr, int orgcnt, MPI_Datatype orgtype, int rank,
1576                       MPI_Aint targdisp, int targcnt, MPI_Datatype targtype, 
1577                       MPI_Op op, WinStruct win);
1578     void winRemoteAccumulate(int orgtotalsize, char* orgaddr, int orgcnt, MPI_Datatype orgtype, MPI_Aint targdisp, 
1579                             int targcnt, MPI_Datatype targtype, 
1580                             MPI_Op op, int winIndex, CkFutureID ftHandle, 
1581                             int pe_src);
1582     void winSetName(WinStruct win, char *name);
1583     void winGetName(WinStruct win, char *name, int *length);
1584     win_obj* getWinObjInstance(WinStruct win); 
1585     int getNewSemaId(); 
1586
1587     AmpiMsg* Alltoall_RemoteIGet(int disp, int targcnt, MPI_Datatype targtype, int tag);
1588 private:
1589     int AlltoallGetFlag;
1590     void *Alltoallbuff;
1591 public:
1592     void setA2AIGetFlag(void* ptr) {AlltoallGetFlag=1;Alltoallbuff=ptr;}
1593     void resetA2AIGetFlag() {AlltoallGetFlag=0;Alltoallbuff=NULL;} 
1594     //------------------------ End of code by YAN ---------------------
1595 };
1596
1597 ampiParent *getAmpiParent(void);
1598 ampi *getAmpiInstance(MPI_Comm comm);
1599 void checkComm(MPI_Comm comm);
1600 void checkRequest(MPI_Request req);
1601
1602 //Use this to mark the start of AMPI interface routines:
1603 #define AMPIAPI(routineName) TCHARM_API_TRACE(routineName,"ampi")
1604
1605 #endif
1606