AMPI: Fix AmpiMsg::pup
[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 #include <algorithm>
6 #include <numeric>
7 #include <forward_list>
8 #include <bitset>
9
10 #include "ampi.h"
11 #include "ddt.h"
12 #include "charm++.h"
13
14 using std::vector;
15
16 //Uncomment for debug print statements
17 #define AMPI_DEBUG(...) //CkPrintf(__VA_ARGS__)
18
19 /*
20  * All MPI_* routines must be defined using the AMPI_API_IMPL macro.
21  * All calls inside AMPI to MPI_* routines must use MPI_* as the name.
22  * There are two reasons for this:
23  *
24  * 1. AMPI supports the PMPI interface only on Linux.
25  *
26  * 2. When AMPI is built on top of MPI, we rename the user's MPI_* calls as AMPI_*.
27  */
28 #define STRINGIFY(a) #a
29
30 #if defined(__linux__)
31   #if CMK_CONVERSE_MPI
32     #define AMPI_API_IMPL(ret, name, ...) \
33       _Pragma(STRINGIFY(weak A##name)) \
34       _Pragma(STRINGIFY(weak AP##name = A##name)) \
35       CLINKAGE \
36       ret A##name(__VA_ARGS__)
37   #else
38     #define AMPI_API_IMPL(ret, name, ...) \
39       _Pragma(STRINGIFY(weak name)) \
40       _Pragma(STRINGIFY(weak P##name = name)) \
41       CLINKAGE \
42       ret name(__VA_ARGS__)
43   #endif
44 #else // not Linux (no PMPI support):
45   #if CMK_CONVERSE_MPI
46     #define AMPI_API_IMPL(ret, name, ...) \
47       CLINKAGE \
48       ret A##name(__VA_ARGS__)
49   #else
50     #define AMPI_API_IMPL(ret, name, ...) \
51       CLINKAGE \
52       ret name(__VA_ARGS__)
53   #endif
54 #endif
55
56 #if AMPIMSGLOG
57 #include "ckliststring.h"
58 static CkListString msgLogRanks;
59 static int msgLogWrite;
60 static int msgLogRead;
61 static char *msgLogFilename;
62
63 #if CMK_USE_ZLIB && 0
64 #include <zlib.h>
65 namespace PUP{
66 class zdisk : public er {
67  protected:
68   gzFile F;//Disk file to read from/write to
69   zdisk(unsigned int type,gzFile f):er(type),F(f) {}
70   zdisk(const zdisk &p);                        //You don't want to copy
71   void operator=(const zdisk &p);       // You don't want to copy
72
73   //For seeking (pack/unpack in different orders)
74   virtual void impl_startSeek(seekBlock &s); /*Begin a seeking block*/
75   virtual int impl_tell(seekBlock &s); /*Give the current offset*/
76   virtual void impl_seek(seekBlock &s,int off); /*Seek to the given offset*/
77 };
78
79 //For packing to a disk file
80 class tozDisk : public zdisk {
81  protected:
82   //Generic bottleneck: pack n items of size itemSize from p.
83   virtual void bytes(void *p,int n,size_t itemSize,dataType t);
84  public:
85   //Write data to the given file pointer
86   // (must be opened for binary write)
87   // You must close the file yourself when done.
88   tozDisk(gzFile f):zdisk(IS_PACKING,f) {}
89 };
90
91 //For unpacking from a disk file
92 class fromzDisk : public zdisk {
93  protected:
94   //Generic bottleneck: unpack n items of size itemSize from p.
95   virtual void bytes(void *p,int n,size_t itemSize,dataType t);
96  public:
97   //Write data to the given file pointer
98   // (must be opened for binary read)
99   // You must close the file yourself when done.
100   fromzDisk(gzFile f):zdisk(IS_UNPACKING,f) {}
101 };
102 }; // namespace PUP
103 #endif
104 #endif // AMPIMSGLOG
105
106 /* AMPI sends messages inline to PE-local destination VPs if: BigSim is not being used and
107  * if tracing is not being used (see bug #1640 for more details on the latter). */
108 #ifndef AMPI_LOCAL_IMPL
109 #define AMPI_LOCAL_IMPL ( !CMK_BIGSIM_CHARM && !CMK_TRACE_ENABLED )
110 #endif
111
112 /* AMPI uses RDMA sends if BigSim is not being used and the underlying comm
113  * layer supports it (except for GNI, which has experimental RDMA support). */
114 #ifndef AMPI_RDMA_IMPL
115 #define AMPI_RDMA_IMPL ( !CMK_BIGSIM_CHARM && CMK_ONESIDED_IMPL && !CMK_CONVERSE_UGNI )
116 #endif
117
118 /* contiguous messages larger than or equal to this threshold are sent via RDMA */
119 #ifndef AMPI_RDMA_THRESHOLD_DEFAULT
120 #if CMK_USE_IBVERBS || CMK_OFI || CMK_CONVERSE_UGNI
121 #define AMPI_RDMA_THRESHOLD_DEFAULT 65536
122 #else
123 #define AMPI_RDMA_THRESHOLD_DEFAULT 32768
124 #endif
125 #endif
126
127 /* contiguous messages larger than or equal to this threshold that are being sent
128  * within a process are sent via RDMA. */
129 #ifndef AMPI_SMP_RDMA_THRESHOLD_DEFAULT
130 #define AMPI_SMP_RDMA_THRESHOLD_DEFAULT 16384
131 #endif
132
133 extern int AMPI_RDMA_THRESHOLD;
134 extern int AMPI_SMP_RDMA_THRESHOLD;
135
136 #define AMPI_ALLTOALL_THROTTLE   64
137 #define AMPI_ALLTOALL_SHORT_MSG  256
138 #if CMK_BIGSIM_CHARM
139 #define AMPI_ALLTOALL_LONG_MSG   4194304
140 #else
141 #define AMPI_ALLTOALL_LONG_MSG   32768
142 #endif
143
144 typedef void (*MPI_MigrateFn)(void);
145
146 /*
147  * AMPI Message Matching (Amm) Interface:
148  * messages are matched on 2 ints: [tag, src]
149  */
150 #define AMM_TAG   0
151 #define AMM_SRC   1
152 #define AMM_NTAGS 2
153
154 // Number of AmmEntry<T>'s in AmmEntryPool for pt2pt msgs:
155 #ifndef AMPI_AMM_PT2PT_POOL_SIZE
156 #define AMPI_AMM_PT2PT_POOL_SIZE 32
157 #endif
158
159 // Number of AmmEntry<T>'s in AmmEntryPool for coll msgs:
160 #ifndef AMPI_AMM_COLL_POOL_SIZE
161 #define AMPI_AMM_COLL_POOL_SIZE 4
162 #endif
163
164 class AmpiRequestList;
165
166 typedef void (*AmmPupMessageFn)(PUP::er& p, void **msg);
167
168 template <class T>
169 class AmmEntry {
170  public:
171   int tags[AMM_NTAGS]; // [tag, src]
172   AmmEntry<T>* next;
173   T msg; // T is either an AmpiRequest* or an AmpiMsg*
174   AmmEntry(T m) noexcept { tags[AMM_TAG] = m->getTag(); tags[AMM_SRC] = m->getSrcRank(); next = NULL; msg = m; }
175   AmmEntry(int tag, int src, T m) noexcept { tags[AMM_TAG] = tag; tags[AMM_SRC] = src; next = NULL; msg = m; }
176   AmmEntry() = default;
177   ~AmmEntry() = default;
178 };
179
180 template <class T, size_t N>
181 class Amm {
182  public:
183   AmmEntry<T>* first;
184   AmmEntry<T>** lasth;
185
186  private:
187   int startIdx;
188   std::bitset<N> validEntries;
189   std::array<AmmEntry<T>, N> entryPool;
190
191  public:
192   Amm() noexcept : first(NULL), lasth(&first), startIdx(0) { validEntries.reset();  }
193   ~Amm() = default;
194   inline AmmEntry<T>* newEntry(int tag, int src, T msg) noexcept {
195     if (validEntries.all()) {
196       return new AmmEntry<T>(tag, src, msg);
197     } else {
198       for (int i=startIdx; i<validEntries.size(); i++) {
199         if (!validEntries[i]) {
200           validEntries[i] = 1;
201           AmmEntry<T>* ent = new (&entryPool[i]) AmmEntry<T>(tag, src, msg);
202           startIdx = i+1;
203           return ent;
204         }
205       }
206       CkAbort("AMPI> failed to find a free entry in pool!");
207       return NULL;
208     }
209   }
210   inline AmmEntry<T>* newEntry(T msg) noexcept {
211     if (validEntries.all()) {
212       return new AmmEntry<T>(msg);
213     } else {
214       for (int i=startIdx; i<validEntries.size(); i++) {
215         if (!validEntries[i]) {
216           validEntries[i] = 1;
217           AmmEntry<T>* ent = new (&entryPool[i]) AmmEntry<T>(msg);
218           startIdx = i+1;
219           return ent;
220         }
221       }
222       CkAbort("AMPI> failed to find a free entry in pool!");
223       return NULL;
224     }
225   }
226   inline void deleteEntry(AmmEntry<T> *ent) noexcept {
227     if (ent >= &entryPool.front() && ent <= &entryPool.back()) {
228       int idx = (int)((intptr_t)ent - (intptr_t)&entryPool.front()) / sizeof(AmmEntry<T>);
229       validEntries[idx] = 0;
230       startIdx = std::min(idx, startIdx);
231     } else {
232       delete ent;
233     }
234   }
235   void freeAll() noexcept;
236   void flushMsgs() noexcept;
237   inline bool match(const int tags1[AMM_NTAGS], const int tags2[AMM_NTAGS]) const noexcept;
238   inline void put(T msg) noexcept;
239   inline void put(int tag, int src, T msg) noexcept;
240   inline T get(int tag, int src, int* rtags=NULL) noexcept;
241   inline T probe(int tag, int src, int* rtags) noexcept;
242   inline int size() const noexcept;
243   void pup(PUP::er& p, AmmPupMessageFn msgpup) noexcept;
244 };
245
246 PUPfunctionpointer(MPI_User_function*)
247
248 /*
249  * OpStruct's are used to lookup an MPI_User_function* and check its commutativity.
250  * They are also used to create AmpiOpHeader's, which are transmitted in reductions
251  * that are user-defined or else lack an equivalent Charm++ reducer type.
252  */
253 class OpStruct {
254  public:
255   MPI_User_function* func;
256   bool isCommutative;
257  private:
258   bool isValid;
259
260  public:
261   OpStruct() = default;
262   OpStruct(MPI_User_function* f) noexcept : func(f), isCommutative(true), isValid(true) {}
263   OpStruct(MPI_User_function* f, bool c) noexcept : func(f), isCommutative(c), isValid(true) {}
264   void init(MPI_User_function* f, bool c) noexcept {
265     func = f;
266     isCommutative = c;
267     isValid = true;
268   }
269   bool isFree() const noexcept { return !isValid; }
270   void free() noexcept { isValid = false; }
271   void pup(PUP::er &p) {
272     p|func;  p|isCommutative;  p|isValid;
273   }
274 };
275
276 class AmpiOpHeader {
277  public:
278   MPI_User_function* func;
279   MPI_Datatype dtype;
280   int len;
281   int szdata;
282   AmpiOpHeader(MPI_User_function* f,MPI_Datatype d,int l,int szd) noexcept :
283     func(f),dtype(d),len(l),szdata(szd) { }
284 };
285
286 //------------------- added by YAN for one-sided communication -----------
287 /* the index is unique within a communicator */
288 class WinStruct{
289  public:
290   MPI_Comm comm;
291   int index;
292
293 private:
294   bool areRecvsPosted;
295   bool inEpoch;
296   vector<int> exposureRankList;
297   vector<int> accessRankList;
298   vector<MPI_Request> requestList;
299
300 public:
301   WinStruct() noexcept : comm(MPI_COMM_NULL), index(-1), areRecvsPosted(false), inEpoch(false) {
302     exposureRankList.clear(); accessRankList.clear(); requestList.clear();
303   }
304   WinStruct(MPI_Comm comm_, int index_) noexcept : comm(comm_), index(index_), areRecvsPosted(false), inEpoch(false) {
305     exposureRankList.clear(); accessRankList.clear(); requestList.clear();
306   }
307   void pup(PUP::er &p) noexcept {
308     p|comm; p|index; p|areRecvsPosted; p|inEpoch; p|exposureRankList; p|accessRankList; p|requestList;
309   }
310   void clearEpochAccess() noexcept {
311     accessRankList.clear(); inEpoch = false;
312   }
313   void clearEpochExposure() noexcept {
314     exposureRankList.clear(); areRecvsPosted = false; requestList.clear(); inEpoch=false;
315   }
316   vector<int>& getExposureRankList() noexcept {return exposureRankList;}
317   vector<int>& getAccessRankList() noexcept {return accessRankList;}
318   void setExposureRankList(vector<int> &tmpExposureRankList) noexcept {exposureRankList = tmpExposureRankList;}
319   void setAccessRankList(vector<int> &tmpAccessRankList) noexcept {accessRankList = tmpAccessRankList;}
320   vector<int>& getRequestList() noexcept {return requestList;}
321   bool AreRecvsPosted() const noexcept {return areRecvsPosted;}
322   void setAreRecvsPosted(bool setR) noexcept {areRecvsPosted = setR;}
323   bool isInEpoch() const noexcept {return inEpoch;}
324   void setInEpoch(bool arg) noexcept {inEpoch = arg;}
325 };
326
327 class lockQueueEntry {
328  public:
329   int requestRank;
330   int lock_type;
331   lockQueueEntry (int _requestRank, int _lock_type) noexcept
332     : requestRank(_requestRank), lock_type(_lock_type) {}
333   lockQueueEntry() = default;
334 };
335
336 typedef CkQ<lockQueueEntry *> LockQueue;
337
338 class ampiParent;
339
340 class win_obj {
341  public:
342   void *baseAddr;
343   MPI_Aint winSize;
344   int disp_unit;
345   MPI_Comm comm;
346
347   int owner; // Rank of owner of the lock, -1 if not locked
348   LockQueue lockQueue; // queue of waiting processors for the lock
349                        // top of queue is the one holding the lock
350                        // queue is empty if lock is not applied
351   std::string winName;
352   bool initflag;
353
354   vector<int> keyvals; // list of keyval attributes
355
356   void setName(const char *src) noexcept;
357   void getName(char *src,int *len) noexcept;
358
359  public:
360   void pup(PUP::er &p) noexcept;
361
362   win_obj() noexcept;
363   win_obj(const char *name, void *base, MPI_Aint size, int disp_unit, MPI_Comm comm) noexcept;
364   ~win_obj() noexcept;
365
366   int create(const char *name, void *base, MPI_Aint size, int disp_unit,
367              MPI_Comm comm) noexcept;
368   int free() noexcept;
369
370   vector<int>& getKeyvals() { return keyvals; }
371
372   int put(void *orgaddr, int orgcnt, int orgunit,
373           MPI_Aint targdisp, int targcnt, int targunit) noexcept;
374
375   int get(void *orgaddr, int orgcnt, int orgunit,
376           MPI_Aint targdisp, int targcnt, int targunit) noexcept;
377   int accumulate(void *orgaddr, int count, MPI_Aint targdisp, MPI_Datatype targtype,
378                  MPI_Op op, ampiParent* pptr) noexcept;
379
380   int iget(int orgcnt, MPI_Datatype orgtype,
381           MPI_Aint targdisp, int targcnt, MPI_Datatype targtype) noexcept;
382   int igetWait(MPI_Request *req, MPI_Status *status) noexcept;
383   int igetFree(MPI_Request *req, MPI_Status *status) noexcept;
384
385   int fence() noexcept;
386
387   int lock(int requestRank, int lock_type) noexcept;
388   int unlock(int requestRank) noexcept;
389
390   int wait() noexcept;
391   int post() noexcept;
392   int start() noexcept;
393   int complete() noexcept;
394
395   void lockTopQueue() noexcept;
396   void enqueue(int requestRank, int lock_type) noexcept;
397   void dequeue() noexcept;
398   bool emptyQueue() noexcept;
399 };
400 //-----------------------End of code by YAN ----------------------
401
402 class KeyvalPair{
403  protected:
404   std::string key;
405   std::string val;
406  public:
407   KeyvalPair() = default;
408   KeyvalPair(const char* k, const char* v) noexcept;
409   ~KeyvalPair() = default;
410   void pup(PUP::er& p) noexcept {
411     p|key;
412     p|val;
413   }
414   friend class InfoStruct;
415 };
416
417 class InfoStruct{
418   CkPupPtrVec<KeyvalPair> nodes;
419   bool valid;
420  public:
421   InfoStruct() noexcept : valid(true) { }
422   void setvalid(bool valid_) noexcept { valid = valid_; }
423   bool getvalid() const noexcept { return valid; }
424   int set(const char* k, const char* v) noexcept;
425   int dup(InfoStruct& src) noexcept;
426   int get(const char* k, int vl, char*& v, int *flag) const noexcept;
427   int deletek(const char* k) noexcept;
428   int get_valuelen(const char* k, int* vl, int *flag) const noexcept;
429   int get_nkeys(int *nkeys) const noexcept;
430   int get_nthkey(int n,char* k) const noexcept;
431   void myfree() noexcept;
432   void pup(PUP::er& p) noexcept;
433 };
434
435 class CProxy_ampi;
436 class CProxyElement_ampi;
437
438 //Virtual class describing a virtual topology: Cart, Graph, DistGraph
439 class ampiTopology {
440  private:
441   vector<int> v; // dummy variable for const& returns from virtual functions
442
443  public:
444   virtual ~ampiTopology() noexcept {};
445   virtual void pup(PUP::er &p) noexcept =0;
446   virtual int getType() const noexcept =0;
447   virtual void dup(ampiTopology* topo) noexcept =0;
448   virtual const vector<int> &getnbors() const noexcept =0;
449   virtual void setnbors(const vector<int> &nbors_) noexcept =0;
450
451   virtual const vector<int> &getdims() const noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class."); return v;}
452   virtual const vector<int> &getperiods() const noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class."); return v;}
453   virtual int getndims() const noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class."); return -1;}
454   virtual void setdims(const vector<int> &dims_) noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class.");}
455   virtual void setperiods(const vector<int> &periods_) noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class.");}
456   virtual void setndims(int ndims_) noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class.");}
457
458   virtual int getnvertices() const noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class."); return -1;}
459   virtual const vector<int> &getindex() const noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class."); return v;}
460   virtual const vector<int> &getedges() const noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class."); return v;}
461   virtual void setnvertices(int nvertices_) noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class.");}
462   virtual void setindex(const vector<int> &index_) noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class.");}
463   virtual void setedges(const vector<int> &edges_) noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class.");}
464
465   virtual int getInDegree() const noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class."); return -1;}
466   virtual const vector<int> &getSources() const noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class."); return v;}
467   virtual const vector<int> &getSourceWeights() const noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class."); return v;}
468   virtual int getOutDegree() const noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class."); return -1;}
469   virtual const vector<int> &getDestinations() const noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class."); return v;}
470   virtual const vector<int> &getDestWeights() const noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class."); return v;}
471   virtual bool areSourcesWeighted() const noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class."); return false;}
472   virtual bool areDestsWeighted() const noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class."); return false;}
473   virtual void setAreSourcesWeighted(bool val) noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class.");}
474   virtual void setAreDestsWeighted(bool val) noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class.");}
475   virtual void setInDegree(int degree) noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class.");}
476   virtual void setSources(const vector<int> &sources) noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class.");}
477   virtual void setSourceWeights(const vector<int> &sourceWeights) noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class.");}
478   virtual void setOutDegree(int degree) noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class.");}
479   virtual void setDestinations(const vector<int> &destinations) noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class.");}
480   virtual void setDestWeights(const vector<int> &destWeights) noexcept {CkAbort("AMPI: instance of invalid Virtual Topology class.");}
481 };
482
483 class ampiCartTopology final : public ampiTopology {
484  private:
485   int ndims;
486   vector<int> dims, periods, nbors;
487
488  public:
489   ampiCartTopology() noexcept : ndims(-1) {}
490
491   void pup(PUP::er &p) noexcept {
492     p|ndims;
493     p|dims;
494     p|periods;
495     p|nbors;
496   }
497
498   inline int getType() const noexcept {return MPI_CART;}
499   inline void dup(ampiTopology* topo) noexcept {
500     CkAssert(topo->getType() == MPI_CART);
501     setndims(topo->getndims());
502     setdims(topo->getdims());
503     setperiods(topo->getperiods());
504     setnbors(topo->getnbors());
505   }
506
507   inline const vector<int> &getdims() const noexcept {return dims;}
508   inline const vector<int> &getperiods() const noexcept {return periods;}
509   inline int getndims() const noexcept {return ndims;}
510   inline const vector<int> &getnbors() const noexcept {return nbors;}
511
512   inline void setdims(const vector<int> &d) noexcept {dims = d; dims.shrink_to_fit();}
513   inline void setperiods(const vector<int> &p) noexcept {periods = p; periods.shrink_to_fit();}
514   inline void setndims(int nd) noexcept {ndims = nd;}
515   inline void setnbors(const vector<int> &n) noexcept {nbors = n; nbors.shrink_to_fit();}
516 };
517
518 class ampiGraphTopology final : public ampiTopology {
519  private:
520   int nvertices;
521   vector<int> index, edges, nbors;
522
523  public:
524   ampiGraphTopology() noexcept : nvertices(-1) {}
525
526   void pup(PUP::er &p) noexcept {
527     p|nvertices;
528     p|index;
529     p|edges;
530     p|nbors;
531   }
532
533   inline int getType() const noexcept {return MPI_GRAPH;}
534   inline void dup(ampiTopology* topo) noexcept {
535     CkAssert(topo->getType() == MPI_GRAPH);
536     setnvertices(topo->getnvertices());
537     setindex(topo->getindex());
538     setedges(topo->getedges());
539     setnbors(topo->getnbors());
540   }
541
542   inline int getnvertices() const noexcept {return nvertices;}
543   inline const vector<int> &getindex() const noexcept {return index;}
544   inline const vector<int> &getedges() const noexcept {return edges;}
545   inline const vector<int> &getnbors() const noexcept {return nbors;}
546
547   inline void setnvertices(int nv) noexcept {nvertices = nv;}
548   inline void setindex(const vector<int> &i) noexcept {index = i; index.shrink_to_fit();}
549   inline void setedges(const vector<int> &e) noexcept {edges = e; edges.shrink_to_fit();}
550   inline void setnbors(const vector<int> &n) noexcept {nbors = n; nbors.shrink_to_fit();}
551 };
552
553 class ampiDistGraphTopology final : public ampiTopology {
554  private:
555   int inDegree, outDegree;
556   bool sourcesWeighted, destsWeighted;
557   vector<int> sources, sourceWeights, destinations, destWeights, nbors;
558
559  public:
560   ampiDistGraphTopology() noexcept : inDegree(-1), outDegree(-1), sourcesWeighted(false), destsWeighted(false) {}
561
562   void pup(PUP::er &p) noexcept {
563     p|inDegree;
564     p|outDegree;
565     p|sourcesWeighted;
566     p|destsWeighted;
567     p|sources;
568     p|sourceWeights;
569     p|destinations;
570     p|destWeights;
571     p|nbors;
572   }
573
574   inline int getType() const noexcept {return MPI_DIST_GRAPH;}
575   inline void dup(ampiTopology* topo) noexcept {
576     CkAssert(topo->getType() == MPI_DIST_GRAPH);
577     setAreSourcesWeighted(topo->areSourcesWeighted());
578     setAreDestsWeighted(topo->areDestsWeighted());
579     setInDegree(topo->getInDegree());
580     setSources(topo->getSources());
581     setSourceWeights(topo->getSourceWeights());
582     setOutDegree(topo->getOutDegree());
583     setDestinations(topo->getDestinations());
584     setDestWeights(topo->getDestWeights());
585     setnbors(topo->getnbors());
586   }
587
588   inline int getInDegree() const noexcept {return inDegree;}
589   inline const vector<int> &getSources() const noexcept {return sources;}
590   inline const vector<int> &getSourceWeights() const noexcept {return sourceWeights;}
591   inline int getOutDegree() const noexcept {return outDegree;}
592   inline const vector<int> &getDestinations() const noexcept {return destinations;}
593   inline const vector<int> &getDestWeights() const noexcept {return destWeights;}
594   inline bool areSourcesWeighted() const noexcept {return sourcesWeighted;}
595   inline bool areDestsWeighted() const noexcept {return destsWeighted;}
596   inline const vector<int> &getnbors() const noexcept {return nbors;}
597
598   inline void setAreSourcesWeighted(bool v) noexcept {sourcesWeighted = v ? 1 : 0;}
599   inline void setAreDestsWeighted(bool v) noexcept {destsWeighted = v ? 1 : 0;}
600   inline void setInDegree(int d) noexcept {inDegree = d;}
601   inline void setSources(const vector<int> &s) noexcept {sources = s; sources.shrink_to_fit();}
602   inline void setSourceWeights(const vector<int> &sw) noexcept {sourceWeights = sw; sourceWeights.shrink_to_fit();}
603   inline void setOutDegree(int d) noexcept {outDegree = d;}
604   inline void setDestinations(const vector<int> &d) noexcept {destinations = d; destinations.shrink_to_fit();}
605   inline void setDestWeights(const vector<int> &dw) noexcept {destWeights = dw; destWeights.shrink_to_fit();}
606   inline void setnbors(const vector<int> &nbors_) noexcept {nbors = nbors_; nbors.shrink_to_fit();}
607 };
608
609 /* KeyValue class for attribute caching */
610 class KeyvalNode {
611  public:
612   void *val;
613   MPI_Copy_function *copy_fn;
614   MPI_Delete_function *delete_fn;
615   void *extra_state;
616   int refCount;
617   bool isValSet;
618
619   KeyvalNode() : val(NULL), copy_fn(NULL), delete_fn(NULL), extra_state(NULL), refCount(1), isValSet(false) { }
620   KeyvalNode(MPI_Copy_function *cf, MPI_Delete_function *df, void* es) :
621              val(NULL), copy_fn(cf), delete_fn(df), extra_state(es), refCount(1), isValSet(false) { }
622   bool hasVal() const { return isValSet; }
623   void clearVal() { isValSet = false; }
624   void setVal(void *v) { val = v; isValSet = true; }
625   void* getVal() const { return val; }
626   void incRefCount() { refCount++; }
627   int decRefCount() { CkAssert(refCount > 0); refCount--; return refCount; }
628   void pup(PUP::er& p) {
629     p((char *)val, sizeof(void *));
630     p((char *)copy_fn, sizeof(void *));
631     p((char *)delete_fn, sizeof(void *));
632     p((char *)extra_state, sizeof(void *));
633     p|refCount;
634     p|isValSet;
635   }
636 };
637
638 // Only store Group ranks explicitly when they can't be
639 // lazily and transiently created via std::iota()
640 class groupStruct {
641  private:
642   int sz; // -1 if ranks is valid, otherwise the size to pass to std::iota()
643   vector<int> ranks;
644
645  private:
646   bool ranksIsIota() const noexcept {
647     for (int i=0; i<ranks.size(); i++)
648       if (ranks[i] != i)
649         return false;
650     return true;
651   }
652
653  public:
654   groupStruct() noexcept : sz(0) {}
655   groupStruct(int s) noexcept : sz(s) {}
656   groupStruct(vector<int> r) noexcept : sz(-1), ranks(std::move(r)) {
657     if (ranksIsIota()) {
658       sz = ranks.size();
659       ranks.clear();
660     }
661     ranks.shrink_to_fit();
662   }
663   groupStruct &operator=(const groupStruct &obj) noexcept {
664     sz = obj.sz;
665     ranks = obj.ranks;
666     return *this;
667   }
668   ~groupStruct() = default;
669   void pup(PUP::er& p) noexcept {
670     p|sz;
671     p|ranks;
672   }
673   bool isIota() const noexcept {return (sz != -1);}
674   int operator[](int i) const noexcept {return (isIota()) ? i : ranks[i];}
675   int size() const noexcept {return (isIota()) ? sz : ranks.size();}
676   vector<int> getRanks() const noexcept {
677     if (isIota()) {
678       // Lazily create ranks:
679       vector<int> tmpRanks(sz);
680       std::iota(tmpRanks.begin(), tmpRanks.end(), 0);
681       tmpRanks.shrink_to_fit();
682       return tmpRanks;
683     }
684     else {
685       return ranks;
686     }
687   }
688 };
689
690 enum AmpiCommType : uint8_t {
691    WORLD = 0
692   ,INTRA = 1
693   ,INTER = 2
694 };
695
696 //Describes an AMPI communicator
697 class ampiCommStruct {
698  private:
699   MPI_Comm comm; //Communicator
700   CkArrayID ampiID; //ID of corresponding ampi array
701   int size; //Number of processes in communicator
702   AmpiCommType commType; //COMM_WORLD, intracomm, intercomm?
703   groupStruct indices;  //indices[r] gives the array index for rank r
704   groupStruct remoteIndices;  // remote group for inter-communicator
705
706   ampiTopology *ampiTopo; // Virtual topology
707   int topoType; // Type of virtual topology: MPI_CART, MPI_GRAPH, MPI_DIST_GRAPH, or MPI_UNDEFINED
708
709   // For communicator attributes (MPI_*_get_attr): indexed by keyval
710   vector<int> keyvals;
711
712   // For communicator names
713   std::string commName;
714
715  public:
716   ampiCommStruct(int ignored=0) noexcept
717     : size(-1), commType(INTRA), ampiTopo(NULL), topoType(MPI_UNDEFINED)
718   {}
719   ampiCommStruct(MPI_Comm comm_,const CkArrayID &id_,int size_) noexcept
720     : comm(comm_), ampiID(id_),size(size_), commType(WORLD), indices(size_),
721       ampiTopo(NULL), topoType(MPI_UNDEFINED)
722   {}
723   ampiCommStruct(MPI_Comm comm_,const CkArrayID &id_, const vector<int> &indices_) noexcept
724     : comm(comm_), ampiID(id_), size(indices_.size()), commType(INTRA), indices(indices_),
725       ampiTopo(NULL), topoType(MPI_UNDEFINED)
726   {}
727   ampiCommStruct(MPI_Comm comm_, const CkArrayID &id_, const vector<int> &indices_,
728                  const vector<int> &remoteIndices_) noexcept
729     : comm(comm_), ampiID(id_), size(indices_.size()), commType(INTER), indices(indices_),
730       remoteIndices(remoteIndices_), ampiTopo(NULL), topoType(MPI_UNDEFINED)
731   {}
732
733   ~ampiCommStruct() noexcept {
734     if (ampiTopo != NULL)
735       delete ampiTopo;
736   }
737
738   // Overloaded copy constructor. Used when creating virtual topologies.
739   ampiCommStruct(const ampiCommStruct &obj, int topoNumber=MPI_UNDEFINED) noexcept {
740     switch (topoNumber) {
741       case MPI_CART:
742         ampiTopo = new ampiCartTopology();
743         break;
744       case MPI_GRAPH:
745         ampiTopo = new ampiGraphTopology();
746         break;
747       case MPI_DIST_GRAPH:
748         ampiTopo = new ampiDistGraphTopology();
749         break;
750       default:
751         ampiTopo = NULL;
752         break;
753     }
754     topoType       = topoNumber;
755     comm           = obj.comm;
756     ampiID         = obj.ampiID;
757     size           = obj.size;
758     commType       = obj.commType;
759     indices        = obj.indices;
760     remoteIndices  = obj.remoteIndices;
761     keyvals        = obj.keyvals;
762     commName       = obj.commName;
763   }
764
765   ampiCommStruct &operator=(const ampiCommStruct &obj) noexcept {
766     if (this == &obj) {
767       return *this;
768     }
769     switch (obj.topoType) {
770       case MPI_CART:
771         ampiTopo = new ampiCartTopology(*(static_cast<ampiCartTopology*>(obj.ampiTopo)));
772         break;
773       case MPI_GRAPH:
774         ampiTopo = new ampiGraphTopology(*(static_cast<ampiGraphTopology*>(obj.ampiTopo)));
775         break;
776       case MPI_DIST_GRAPH:
777         ampiTopo = new ampiDistGraphTopology(*(static_cast<ampiDistGraphTopology*>(obj.ampiTopo)));
778         break;
779       default:
780         ampiTopo = NULL;
781         break;
782     }
783     topoType       = obj.topoType;
784     comm           = obj.comm;
785     ampiID         = obj.ampiID;
786     size           = obj.size;
787     commType       = obj.commType;
788     indices        = obj.indices;
789     remoteIndices  = obj.remoteIndices;
790     keyvals        = obj.keyvals;
791     commName       = obj.commName;
792     return *this;
793   }
794
795   const ampiTopology* getTopologyforNeighbors() const noexcept {
796     return ampiTopo;
797   }
798
799   ampiTopology* getTopology() noexcept {
800     return ampiTopo;
801   }
802
803   inline bool isinter() const noexcept {return commType==INTER;}
804   void setArrayID(const CkArrayID &nID) noexcept {ampiID=nID;}
805
806   MPI_Comm getComm() const noexcept {return comm;}
807   inline vector<int> getIndices() const noexcept {return indices.getRanks();}
808   inline vector<int> getRemoteIndices() const noexcept {return remoteIndices.getRanks();}
809   vector<int> &getKeyvals() noexcept {return keyvals;}
810
811   void setName(const char *src) noexcept {
812     CkDDT_SetName(commName, src);
813   }
814
815   void getName(char *name, int *len) const noexcept {
816     int length = *len = commName.size();
817     memcpy(name, commName.data(), length);
818     name[length] = '\0';
819   }
820
821   //Get the proxy for the entire array
822   CProxy_ampi getProxy() const noexcept;
823
824   //Get the array index for rank r in this communicator
825   int getIndexForRank(int r) const noexcept {
826 #if CMK_ERROR_CHECKING
827     if (r>=size) CkAbort("AMPI> You passed in an out-of-bounds process rank!");
828 #endif
829     return indices[r];
830   }
831   int getIndexForRemoteRank(int r) const noexcept {
832 #if CMK_ERROR_CHECKING
833     if (r>=remoteIndices.size()) CkAbort("AMPI> You passed in an out-of-bounds intercomm remote process rank!");
834 #endif
835     return remoteIndices[r];
836   }
837   //Get the rank for this array index (Warning: linear time)
838   int getRankForIndex(int i) const noexcept {
839     if (indices.isIota()) return i;
840     else {
841       const vector<int>& ind = indices.getRanks();
842       for (int r=0;r<ind.size();r++)
843         if (ind[r]==i) return r;
844       return -1; /*That index isn't in this communicator*/
845     }
846   }
847
848   int getSize() const noexcept {return size;}
849
850   void pup(PUP::er &p) noexcept {
851     p|comm;
852     p|ampiID;
853     p|size;
854     p|commType;
855     p|indices;
856     p|remoteIndices;
857     p|keyvals;
858     p|commName;
859     p|topoType;
860     if (topoType != MPI_UNDEFINED) {
861       if (p.isUnpacking()) {
862         switch (topoType) {
863           case MPI_CART:
864             ampiTopo = new ampiCartTopology();
865             break;
866           case MPI_GRAPH:
867             ampiTopo = new ampiGraphTopology();
868             break;
869           case MPI_DIST_GRAPH:
870             ampiTopo = new ampiDistGraphTopology();
871             break;
872           default:
873             CkAbort("AMPI> Communicator has an invalid topology!");
874             break;
875         }
876       }
877       ampiTopo->pup(p);
878     } else {
879       ampiTopo = NULL;
880     }
881     if (p.isDeleting()) {
882       delete ampiTopo; ampiTopo = NULL;
883     }
884   }
885 };
886 PUPmarshall(ampiCommStruct)
887
888 class mpi_comm_worlds{
889   ampiCommStruct comms[MPI_MAX_COMM_WORLDS];
890  public:
891   ampiCommStruct &operator[](int i) noexcept {return comms[i];}
892   void pup(PUP::er &p) noexcept {
893     for (int i=0;i<MPI_MAX_COMM_WORLDS;i++)
894       comms[i].pup(p);
895   }
896 };
897
898 // group operations
899 inline void outputOp(const vector<int>& vec) noexcept {
900   if (vec.size() > 50) {
901     CkPrintf("vector too large to output!\n");
902     return;
903   }
904   CkPrintf("output vector: size=%d {",vec.size());
905   for (int i=0; i<vec.size(); i++) {
906     CkPrintf(" %d ", vec[i]);
907   }
908   CkPrintf("}\n");
909 }
910
911 inline int getPosOp(int idx, const vector<int>& vec) noexcept {
912   for (int r=0; r<vec.size(); r++) {
913     if (vec[r] == idx) {
914       return r;
915     }
916   }
917   return MPI_UNDEFINED;
918 }
919
920 inline vector<int> unionOp(const vector<int>& vec1, const vector<int>& vec2) noexcept {
921   vector<int> newvec(vec1);
922   for (int i=0; i<vec2.size(); i++) {
923     if (getPosOp(vec2[i], vec1) == MPI_UNDEFINED) {
924       newvec.push_back(vec2[i]);
925     }
926   }
927   return newvec;
928 }
929
930 inline vector<int> intersectOp(const vector<int>& vec1, const vector<int>& vec2) noexcept {
931   vector<int> newvec;
932   for (int i=0; i<vec1.size(); i++) {
933     if (getPosOp(vec1[i], vec2) != MPI_UNDEFINED) {
934       newvec.push_back(vec1[i]);
935     }
936   }
937   return newvec;
938 }
939
940 inline vector<int> diffOp(const vector<int>& vec1, const vector<int>& vec2) noexcept {
941   vector<int> newvec;
942   for (int i=0; i<vec1.size(); i++) {
943     if (getPosOp(vec1[i], vec2) == MPI_UNDEFINED) {
944       newvec.push_back(vec1[i]);
945     }
946   }
947   return newvec;
948 }
949
950 inline int* translateRanksOp(int n, const vector<int>& vec1, const int* ranks1,
951                              const vector<int>& vec2, int *ret) noexcept {
952   for (int i=0; i<n; i++) {
953     ret[i] = (ranks1[i] == MPI_PROC_NULL) ? MPI_PROC_NULL : getPosOp(vec1[ranks1[i]], vec2);
954   }
955   return ret;
956 }
957
958 inline int compareVecOp(const vector<int>& vec1, const vector<int>& vec2) noexcept {
959   int pos, ret = MPI_IDENT;
960   if (vec1.size() != vec2.size()) {
961     return MPI_UNEQUAL;
962   }
963   for (int i=0; i<vec1.size(); i++) {
964     pos = getPosOp(vec1[i], vec2);
965     if (pos == MPI_UNDEFINED) {
966       return MPI_UNEQUAL;
967     }
968     else if (pos != i) {
969       ret = MPI_SIMILAR;
970     }
971   }
972   return ret;
973 }
974
975 inline vector<int> inclOp(int n, const int* ranks, const vector<int>& vec) noexcept {
976   vector<int> retvec(n);
977   for (int i=0; i<n; i++) {
978     retvec[i] = vec[ranks[i]];
979   }
980   return retvec;
981 }
982
983 inline vector<int> exclOp(int n, const int* ranks, const vector<int>& vec) noexcept {
984   vector<int> retvec;
985   bool add = true;
986   for (int j=0; j<vec.size(); j++) {
987     for (int i=0; i<n; i++) {
988       if (j == ranks[i]) {
989         add = false;
990         break;
991       }
992     }
993     if (add) {
994       retvec.push_back(vec[j]);
995     }
996     else {
997       add = true;
998     }
999   }
1000   return retvec;
1001 }
1002
1003 inline vector<int> rangeInclOp(int n, int ranges[][3], const vector<int>& vec,
1004                                int *flag) noexcept {
1005   vector<int> retvec;
1006   int first, last, stride;
1007   for (int i=0; i<n; i++) {
1008     first  = ranges[i][0];
1009     last   = ranges[i][1];
1010     stride = ranges[i][2];
1011     if (stride != 0) {
1012       for (int j=0; j<=(last-first)/stride; j++) {
1013         retvec.push_back(vec[first+stride*j]);
1014       }
1015     }
1016     else {
1017       *flag = MPI_ERR_ARG;
1018       return vector<int>();
1019     }
1020   }
1021   *flag = MPI_SUCCESS;
1022   return retvec;
1023 }
1024
1025 inline vector<int> rangeExclOp(int n, int ranges[][3], const vector<int>& vec,
1026                                int *flag) noexcept {
1027   vector<int> ranks;
1028   int first, last, stride;
1029   for (int i=0; i<n; i++) {
1030     first  = ranges[i][0];
1031     last   = ranges[i][1];
1032     stride = ranges[i][2];
1033     if (stride != 0) {
1034       for (int j=0; j<=(last-first)/stride; j++) {
1035         ranks.push_back(first+stride*j);
1036       }
1037     }
1038     else {
1039       *flag = MPI_ERR_ARG;
1040       return vector<int>();
1041     }
1042   }
1043   *flag = MPI_SUCCESS;
1044   return exclOp(ranks.size(), &ranks[0], vec);
1045 }
1046
1047 #include "tcharm.h"
1048 #include "tcharmc.h"
1049
1050 #include "ampi.decl.h"
1051 #include "charm-api.h"
1052 #include <sys/stat.h> // for mkdir
1053
1054 extern int _mpi_nworlds;
1055
1056 //MPI_ANY_TAG is defined in ampi.h to MPI_TAG_UB_VALUE+1
1057 #define MPI_ATA_SEQ_TAG     MPI_TAG_UB_VALUE+2
1058 #define MPI_BCAST_TAG       MPI_TAG_UB_VALUE+3
1059 #define MPI_REDN_TAG        MPI_TAG_UB_VALUE+4
1060 #define MPI_SCATTER_TAG     MPI_TAG_UB_VALUE+5
1061 #define MPI_SCAN_TAG        MPI_TAG_UB_VALUE+6
1062 #define MPI_EXSCAN_TAG      MPI_TAG_UB_VALUE+7
1063 #define MPI_ATA_TAG         MPI_TAG_UB_VALUE+8
1064 #define MPI_NBOR_TAG        MPI_TAG_UB_VALUE+9
1065 #define MPI_RMA_TAG         MPI_TAG_UB_VALUE+10
1066 #define MPI_EPOCH_START_TAG MPI_TAG_UB_VALUE+11
1067 #define MPI_EPOCH_END_TAG   MPI_TAG_UB_VALUE+12
1068
1069 #define AMPI_COLL_SOURCE 0
1070 #define AMPI_COLL_COMM   MPI_COMM_WORLD
1071
1072 enum AmpiReqType : uint8_t {
1073   AMPI_INVALID_REQ = 0,
1074   AMPI_I_REQ       = 1,
1075   AMPI_ATA_REQ     = 2,
1076   AMPI_SEND_REQ    = 3,
1077   AMPI_SSEND_REQ   = 4,
1078   AMPI_REDN_REQ    = 5,
1079   AMPI_GATHER_REQ  = 6,
1080   AMPI_GATHERV_REQ = 7,
1081   AMPI_G_REQ       = 8,
1082 #if CMK_CUDA
1083   AMPI_GPU_REQ     = 9
1084 #endif
1085 };
1086
1087 inline void operator|(PUP::er &p, AmpiReqType &r) {
1088   pup_bytes(&p, (void *)&r, sizeof(AmpiReqType));
1089 }
1090
1091 enum AmpiReqSts : char {
1092   AMPI_REQ_PENDING   = 0,
1093   AMPI_REQ_BLOCKED   = 1,
1094   AMPI_REQ_COMPLETED = 2
1095 };
1096
1097 enum AmpiSendType : bool {
1098   BLOCKING_SEND = false,
1099   I_SEND = true
1100 };
1101
1102 #define MyAlign8(x) (((x)+7)&(~7))
1103
1104 /**
1105 Represents an MPI request that has been initiated
1106 using Isend, Irecv, Ialltoall, Send_init, etc.
1107 */
1108 class AmpiRequest {
1109  public:
1110   void *buf          = nullptr;
1111   int count          = 0;
1112   MPI_Datatype type  = MPI_DATATYPE_NULL;
1113   int tag            = MPI_ANY_TAG; // the order must match MPI_Status
1114   int src            = MPI_ANY_SOURCE;
1115   MPI_Comm comm      = MPI_COMM_NULL;
1116   MPI_Request reqIdx = MPI_REQUEST_NULL;
1117   bool complete      = false;
1118   bool blocked       = false; // this req is currently blocked on
1119
1120 #if CMK_BIGSIM_CHARM
1121  public:
1122   void *event        = nullptr; // the event point that corresponds to this message
1123   int eventPe        = -1; // the PE that the event is located on
1124 #endif
1125
1126  public:
1127   AmpiRequest() =default;
1128   virtual ~AmpiRequest() =default;
1129
1130   /// Activate this persistent request.
1131   ///  Only meaningful for persistent Ireq, SendReq, and SsendReq requests.
1132   virtual void start(MPI_Request reqIdx) noexcept {}
1133
1134   /// Used by AmmEntry's constructor
1135   virtual int getTag() const noexcept { return tag; }
1136   virtual int getSrcRank() const noexcept { return src; }
1137
1138   /// Return true if this request is finished (progress):
1139   virtual bool test(MPI_Status *sts=MPI_STATUS_IGNORE) noexcept =0;
1140
1141   /// Block until this request is finished,
1142   ///  returning a valid MPI error code.
1143   virtual int wait(MPI_Status *sts) noexcept =0;
1144
1145   /// Mark this request for cancellation.
1146   /// Supported only for IReq requests
1147   virtual void cancel() noexcept {}
1148
1149   /// Mark this request persistent.
1150   /// Supported only for IReq, SendReq, and SsendReq requests
1151   virtual void setPersistent(bool p) noexcept {}
1152   virtual bool isPersistent() const noexcept { return false; }
1153
1154   /// Receive an AmpiMsg
1155   virtual void receive(ampi *ptr, AmpiMsg *msg, bool deleteMsg=true) noexcept =0;
1156
1157   /// Receive a CkReductionMsg
1158   virtual void receive(ampi *ptr, CkReductionMsg *msg) noexcept =0;
1159
1160   /// Receive an Rdma message
1161   virtual void receiveRdma(ampi *ptr, char *sbuf, int slength, int ssendReq,
1162                            int srcRank, MPI_Comm scomm) noexcept { }
1163
1164   /// Set the request's index into AmpiRequestList
1165   void setReqIdx(MPI_Request idx) noexcept { reqIdx = idx; }
1166   MPI_Request getReqIdx() const noexcept { return reqIdx; }
1167
1168   /// Free the request's datatype
1169   void free(CkDDT* ddt) noexcept {
1170     if (type != MPI_DATATYPE_NULL) ddt->freeType(type);
1171   }
1172
1173   /// Set whether the request is currently blocked on
1174   void setBlocked(bool b) noexcept { blocked = b; }
1175   bool isBlocked() const noexcept { return blocked; }
1176
1177   /// Returns the type of request:
1178   ///  AMPI_I_REQ, AMPI_ATA_REQ, AMPI_SEND_REQ, AMPI_SSEND_REQ,
1179   ///  AMPI_REDN_REQ, AMPI_GATHER_REQ, AMPI_GATHERV_REQ, AMPI_G_REQ
1180   virtual AmpiReqType getType() const noexcept =0;
1181
1182   /// Returns whether this request will need to be matched.
1183   /// It is used to determine whether this request should be inserted into postedReqs.
1184   /// AMPI_SEND_REQ, AMPI_SSEND_REQ, and AMPI_ATA_REQ should not be posted.
1185   virtual bool isUnmatched() const noexcept =0;
1186
1187   /// Returns whether this type is pooled or not:
1188   /// Only AMPI_I_REQ, AMPI_SEND_REQ, and AMPI_SSEND_REQs are pooled.
1189   virtual bool isPooledType() const noexcept { return false; }
1190
1191   /// Return the actual number of bytes that were received.
1192   virtual int getNumReceivedBytes(CkDDT *ddt) const noexcept {
1193     // by default, return number of bytes requested
1194     return count * ddt->getSize(type);
1195   }
1196
1197   virtual void pup(PUP::er &p) noexcept {
1198     p((char *)&buf, sizeof(void *)); //supposed to work only with Isomalloc
1199     p(count);
1200     p(type);
1201     p(tag);
1202     p(src);
1203     p(comm);
1204     p(reqIdx);
1205     p(complete);
1206     p(blocked);
1207 #if CMK_BIGSIM_CHARM
1208     //needed for bigsim out-of-core emulation
1209     //as the "log" is not moved from memory, this pointer is safe
1210     //to be reused
1211     p((char *)&event, sizeof(void *));
1212     p(eventPe);
1213 #endif
1214   }
1215
1216   virtual void print() const noexcept =0;
1217 };
1218
1219 // This is used in the constructors of the AmpiRequest types below,
1220 // assuming arguments: (MPI_Datatype type_, CkDDT* ddt_, AmpiReqSts sts_)
1221 #define AMPI_REQUEST_COMMON_INIT           \
1222 {                                          \
1223   complete = (sts_ == AMPI_REQ_COMPLETED); \
1224   blocked  = (sts_ == AMPI_REQ_BLOCKED);   \
1225   if (type_ != MPI_DATATYPE_NULL) {        \
1226     ddt_->getType(type_)->incRefCount();   \
1227   }                                        \
1228 }
1229
1230 class IReq final : public AmpiRequest {
1231  public:
1232   bool cancelled  = false; // track if request is cancelled
1233   bool persistent = false; // Is this a persistent recv request?
1234   int length      = 0; // recv'ed length in bytes
1235
1236   IReq(void *buf_, int count_, MPI_Datatype type_, int src_, int tag_,
1237        MPI_Comm comm_, CkDDT *ddt_, AmpiReqSts sts_=AMPI_REQ_PENDING) noexcept
1238   {
1239     buf   = buf_;
1240     count = count_;
1241     type  = type_;
1242     src   = src_;
1243     tag   = tag_;
1244     comm  = comm_;
1245     AMPI_REQUEST_COMMON_INIT
1246   }
1247   IReq() =default;
1248   ~IReq() =default;
1249   bool test(MPI_Status *sts=MPI_STATUS_IGNORE) noexcept override;
1250   int wait(MPI_Status *sts) noexcept override;
1251   void cancel() noexcept override { if (!complete) cancelled = true; }
1252   AmpiReqType getType() const noexcept override { return AMPI_I_REQ; }
1253   bool isUnmatched() const noexcept override { return !complete; }
1254   bool isPooledType() const noexcept override { return true; }
1255   void setPersistent(bool p) noexcept override { persistent = p; }
1256   bool isPersistent() const noexcept override { return persistent; }
1257   void start(MPI_Request reqIdx) noexcept override;
1258   void receive(ampi *ptr, AmpiMsg *msg, bool deleteMsg=true) noexcept override;
1259   void receive(ampi *ptr, CkReductionMsg *msg) noexcept override {}
1260   void receiveRdma(ampi *ptr, char *sbuf, int slength, int ssendReq, int srcRank, MPI_Comm scomm) noexcept override;
1261   int getNumReceivedBytes(CkDDT *ptr) const noexcept override {
1262     return length;
1263   }
1264   void pup(PUP::er &p) noexcept override {
1265     AmpiRequest::pup(p);
1266     p|cancelled;
1267     p|persistent;
1268     p|length;
1269   }
1270   void print() const noexcept override;
1271 };
1272
1273 class RednReq final : public AmpiRequest {
1274  public:
1275   MPI_Op op = MPI_OP_NULL;
1276
1277   RednReq(void *buf_, int count_, MPI_Datatype type_, MPI_Comm comm_,
1278           MPI_Op op_, CkDDT* ddt_, AmpiReqSts sts_=AMPI_REQ_PENDING) noexcept
1279   {
1280     buf   = buf_;
1281     count = count_;
1282     type  = type_;
1283     src   = AMPI_COLL_SOURCE;
1284     tag   = MPI_REDN_TAG;
1285     comm  = comm_;
1286     op    = op_;
1287     AMPI_REQUEST_COMMON_INIT
1288   }
1289   RednReq() =default;
1290   ~RednReq() =default;
1291   bool test(MPI_Status *sts=MPI_STATUS_IGNORE) noexcept override;
1292   int wait(MPI_Status *sts) noexcept override;
1293   void cancel() noexcept override {}
1294   AmpiReqType getType() const noexcept override { return AMPI_REDN_REQ; }
1295   bool isUnmatched() const noexcept override { return !complete; }
1296   void receive(ampi *ptr, AmpiMsg *msg, bool deleteMsg=true) noexcept override {}
1297   void receive(ampi *ptr, CkReductionMsg *msg) noexcept override;
1298   void pup(PUP::er &p) noexcept override {
1299     AmpiRequest::pup(p);
1300     p|op;
1301   }
1302   void print() const noexcept override;
1303 };
1304
1305 class GatherReq final : public AmpiRequest {
1306  public:
1307   GatherReq(void *buf_, int count_, MPI_Datatype type_, MPI_Comm comm_,
1308             CkDDT *ddt_, AmpiReqSts sts_=AMPI_REQ_PENDING) noexcept
1309   {
1310     buf   = buf_;
1311     count = count_;
1312     type  = type_;
1313     src   = AMPI_COLL_SOURCE;
1314     tag   = MPI_REDN_TAG;
1315     comm  = comm_;
1316     AMPI_REQUEST_COMMON_INIT
1317   }
1318   GatherReq() =default;
1319   ~GatherReq() =default;
1320   bool test(MPI_Status *sts=MPI_STATUS_IGNORE) noexcept override;
1321   int wait(MPI_Status *sts) noexcept override;
1322   void cancel() noexcept override {}
1323   AmpiReqType getType() const noexcept override { return AMPI_GATHER_REQ; }
1324   bool isUnmatched() const noexcept override { return !complete; }
1325   void receive(ampi *ptr, AmpiMsg *msg, bool deleteMsg=true) noexcept override {}
1326   void receive(ampi *ptr, CkReductionMsg *msg) noexcept override;
1327   void pup(PUP::er &p) noexcept override {
1328     AmpiRequest::pup(p);
1329   }
1330   void print() const noexcept override;
1331 };
1332
1333 class GathervReq final : public AmpiRequest {
1334  public:
1335   vector<int> recvCounts;
1336   vector<int> displs;
1337
1338   GathervReq(void *buf_, int count_, MPI_Datatype type_, MPI_Comm comm_, const int *rc,
1339              const int *d, CkDDT* ddt_, AmpiReqSts sts_=AMPI_REQ_PENDING) noexcept
1340   {
1341     buf   = buf_;
1342     count = count_;
1343     type  = type_;
1344     src   = AMPI_COLL_SOURCE;
1345     tag   = MPI_REDN_TAG;
1346     comm  = comm_;
1347     recvCounts.assign(rc, rc+count);
1348     displs.assign(d, d+count);
1349     AMPI_REQUEST_COMMON_INIT
1350   }
1351   GathervReq() =default;
1352   ~GathervReq() =default;
1353   bool test(MPI_Status *sts=MPI_STATUS_IGNORE) noexcept override;
1354   int wait(MPI_Status *sts) noexcept override;
1355   AmpiReqType getType() const noexcept override { return AMPI_GATHERV_REQ; }
1356   bool isUnmatched() const noexcept override { return !complete; }
1357   void receive(ampi *ptr, AmpiMsg *msg, bool deleteMsg=true) noexcept override {}
1358   void receive(ampi *ptr, CkReductionMsg *msg) noexcept override;
1359   void pup(PUP::er &p) noexcept override {
1360     AmpiRequest::pup(p);
1361     p|recvCounts;
1362     p|displs;
1363   }
1364   void print() const noexcept override;
1365 };
1366
1367 class SendReq final : public AmpiRequest {
1368   bool persistent = false; // is this a persistent send request?
1369
1370  public:
1371   SendReq(MPI_Datatype type_, MPI_Comm comm_, CkDDT* ddt_, AmpiReqSts sts_=AMPI_REQ_PENDING) noexcept
1372   {
1373     type = type_;
1374     comm = comm_;
1375     AMPI_REQUEST_COMMON_INIT
1376   }
1377   SendReq(void* buf_, int count_, MPI_Datatype type_, int dest_, int tag_,
1378           MPI_Comm comm_, CkDDT* ddt_, AmpiReqSts sts_=AMPI_REQ_PENDING) noexcept
1379   {
1380     buf   = buf_;
1381     count = count_;
1382     type  = type_;
1383     src   = dest_;
1384     tag   = tag_;
1385     comm  = comm_;
1386     AMPI_REQUEST_COMMON_INIT
1387   }
1388   SendReq() noexcept {}
1389   ~SendReq() noexcept {}
1390   bool test(MPI_Status *sts=MPI_STATUS_IGNORE) noexcept override;
1391   int wait(MPI_Status *sts) noexcept override;
1392   void setPersistent(bool p) noexcept override { persistent = p; }
1393   bool isPersistent() const noexcept override { return persistent; }
1394   void start(MPI_Request reqIdx) noexcept override;
1395   void receive(ampi *ptr, AmpiMsg *msg, bool deleteMsg=true) noexcept override {}
1396   void receive(ampi *ptr, CkReductionMsg *msg) noexcept override {}
1397   AmpiReqType getType() const noexcept override { return AMPI_SEND_REQ; }
1398   bool isUnmatched() const noexcept override { return false; }
1399   bool isPooledType() const noexcept override { return true; }
1400   void pup(PUP::er &p) noexcept override {
1401     AmpiRequest::pup(p);
1402     p|persistent;
1403   }
1404   void print() const noexcept override;
1405 };
1406
1407 class SsendReq final : public AmpiRequest {
1408  private:
1409   bool persistent = false; // is this a persistent Ssend request?
1410
1411  public:
1412   SsendReq(MPI_Datatype type_, MPI_Comm comm_, CkDDT* ddt_, AmpiReqSts sts_=AMPI_REQ_PENDING) noexcept
1413   {
1414     type = type_;
1415     comm = comm_;
1416     AMPI_REQUEST_COMMON_INIT
1417   }
1418   SsendReq(void* buf_, int count_, MPI_Datatype type_, int dest_, int tag_, MPI_Comm comm_,
1419            CkDDT* ddt_, AmpiReqSts sts_=AMPI_REQ_PENDING) noexcept
1420   {
1421     buf   = buf_;
1422     count = count_;
1423     type  = type_;
1424     src   = dest_;
1425     tag   = tag_;
1426     comm  = comm_;
1427     AMPI_REQUEST_COMMON_INIT
1428   }
1429   SsendReq(void* buf_, int count_, MPI_Datatype type_, int dest_, int tag_, MPI_Comm comm_,
1430            int src_, CkDDT* ddt_, AmpiReqSts sts_=AMPI_REQ_PENDING) noexcept
1431   {
1432     buf   = buf_;
1433     count = count_;
1434     type  = type_;
1435     src   = dest_;
1436     tag   = tag_;
1437     comm  = comm_;
1438     AMPI_REQUEST_COMMON_INIT
1439   }
1440   SsendReq() =default;
1441   ~SsendReq() =default;
1442   bool test(MPI_Status *sts=MPI_STATUS_IGNORE) noexcept override;
1443   int wait(MPI_Status *sts) noexcept override;
1444   void setPersistent(bool p) noexcept override { persistent = p; }
1445   bool isPersistent() const noexcept override { return persistent; }
1446   void start(MPI_Request reqIdx) noexcept override;
1447   void receive(ampi *ptr, AmpiMsg *msg, bool deleteMsg=true) noexcept override {}
1448   void receive(ampi *ptr, CkReductionMsg *msg) noexcept override {}
1449   AmpiReqType getType() const noexcept override { return AMPI_SSEND_REQ; }
1450   bool isUnmatched() const noexcept override { return false; }
1451   bool isPooledType() const noexcept override { return true; }
1452   void pup(PUP::er &p) noexcept override {
1453     AmpiRequest::pup(p);
1454     p|persistent;
1455   }
1456   void print() const noexcept override;
1457 };
1458
1459 #if CMK_CUDA
1460 class GPUReq : public AmpiRequest {
1461  public:
1462   GPUReq() noexcept;
1463   ~GPUReq() =default;
1464   bool test(MPI_Status *sts=MPI_STATUS_IGNORE) noexcept override;
1465   int wait(MPI_Status *sts) noexcept override;
1466   void receive(ampi *ptr, AmpiMsg *msg, bool deleteMsg=true) noexcept override;
1467   void receive(ampi *ptr, CkReductionMsg *msg) noexcept override;
1468   AmpiReqType getType() const noexcept override { return AMPI_GPU_REQ; }
1469   bool isUnmatched() const noexcept override { return false; }
1470   void setComplete() noexcept;
1471   void print() const noexcept override;
1472 };
1473 #endif
1474
1475 class ATAReq final : public AmpiRequest {
1476  public:
1477   vector<MPI_Request> reqs;
1478
1479   ATAReq(int numReqs_) noexcept : reqs(numReqs_) {}
1480   ATAReq() =default;
1481   ~ATAReq() =default;
1482   bool test(MPI_Status *sts=MPI_STATUS_IGNORE) noexcept override;
1483   int wait(MPI_Status *sts) noexcept override;
1484   void receive(ampi *ptr, AmpiMsg *msg, bool deleteMsg=true) noexcept override {}
1485   void receive(ampi *ptr, CkReductionMsg *msg) noexcept override {}
1486   int getCount() const noexcept { return reqs.size(); }
1487   AmpiReqType getType() const noexcept override { return AMPI_ATA_REQ; }
1488   bool isUnmatched() const noexcept override { return false; }
1489   void pup(PUP::er &p) noexcept override {
1490     AmpiRequest::pup(p);
1491     p|reqs;
1492   }
1493   void print() const noexcept override;
1494 };
1495
1496 class GReq final : public AmpiRequest {
1497  private:
1498   MPI_Grequest_query_function* queryFn;
1499   MPI_Grequest_free_function* freeFn;
1500   MPI_Grequest_cancel_function* cancelFn;
1501   MPIX_Grequest_poll_function* pollFn;
1502   MPIX_Grequest_wait_function* waitFn;
1503   void* extraState;
1504
1505  public:
1506   GReq(MPI_Grequest_query_function* q, MPI_Grequest_free_function* f, MPI_Grequest_cancel_function* c, void* es) noexcept
1507     : queryFn(q), freeFn(f), cancelFn(c), pollFn(nullptr), waitFn(nullptr), extraState(es) {}
1508   GReq(MPI_Grequest_query_function *q, MPI_Grequest_free_function* f, MPI_Grequest_cancel_function* c, MPIX_Grequest_poll_function* p, void* es) noexcept
1509     : queryFn(q), freeFn(f), cancelFn(c), pollFn(p), waitFn(nullptr), extraState(es) {}
1510   GReq(MPI_Grequest_query_function *q, MPI_Grequest_free_function* f, MPI_Grequest_cancel_function* c, MPIX_Grequest_poll_function* p, MPIX_Grequest_wait_function* w, void* es) noexcept
1511     : queryFn(q), freeFn(f), cancelFn(c), pollFn(p), waitFn(w), extraState(es) {}
1512   GReq() =default;
1513   ~GReq() noexcept { (*freeFn)(extraState); }
1514   bool test(MPI_Status *sts=MPI_STATUS_IGNORE) noexcept override;
1515   int wait(MPI_Status *sts) noexcept override;
1516   void receive(ampi *ptr, AmpiMsg *msg, bool deleteMsg=true) noexcept override {}
1517   void receive(ampi *ptr, CkReductionMsg *msg) noexcept override {}
1518   void cancel() noexcept override { (*cancelFn)(extraState, complete); }
1519   AmpiReqType getType() const noexcept override { return AMPI_G_REQ; }
1520   bool isUnmatched() const noexcept override { return false; }
1521   void pup(PUP::er &p) noexcept override {
1522     AmpiRequest::pup(p);
1523     p((char *)queryFn, sizeof(void *));
1524     p((char *)freeFn, sizeof(void *));
1525     p((char *)cancelFn, sizeof(void *));
1526     p((char *)pollFn, sizeof(void *));
1527     p((char *)waitFn, sizeof(void *));
1528     p((char *)extraState, sizeof(void *));
1529   }
1530   void print() const noexcept override;
1531 };
1532
1533 class AmpiRequestPool;
1534
1535 class AmpiRequestList {
1536  private:
1537   vector<AmpiRequest*> reqs; // indexed by MPI_Request
1538   int startIdx; // start next search from this index
1539   AmpiRequestPool* reqPool;
1540  public:
1541   AmpiRequestList() noexcept : startIdx(0) {}
1542   AmpiRequestList(int size, AmpiRequestPool* reqPoolPtr) noexcept
1543     : reqs(size), startIdx(0), reqPool(reqPoolPtr) {}
1544   ~AmpiRequestList() noexcept {}
1545
1546   inline AmpiRequest* operator[](int n) noexcept {
1547 #if CMK_ERROR_CHECKING
1548     return reqs.at(n);
1549 #else
1550     return reqs[n];
1551 #endif
1552   }
1553   void free(AmpiRequestPool& reqPool, int idx, CkDDT *ddt) noexcept;
1554   void freeNonPersReq(int &idx) noexcept;
1555   inline int insert(AmpiRequest* req) noexcept {
1556     for (int i=startIdx; i<reqs.size(); i++) {
1557       if (reqs[i] == NULL) {
1558         req->setReqIdx(i);
1559         reqs[i] = req;
1560         startIdx = i+1;
1561         return i;
1562       }
1563     }
1564     reqs.push_back(req);
1565     int idx = reqs.size()-1;
1566     req->setReqIdx(idx);
1567     startIdx = idx+1;
1568     return idx;
1569   }
1570
1571   inline void checkRequest(MPI_Request idx) const noexcept {
1572     if (idx != MPI_REQUEST_NULL && (idx < 0 || idx >= reqs.size()))
1573       CkAbort("Invalid MPI_Request\n");
1574   }
1575
1576   inline void unblockReqs(MPI_Request *requests, int numReqs) noexcept {
1577     for (int i=0; i<numReqs; i++) {
1578       if (requests[i] != MPI_REQUEST_NULL) {
1579         reqs[requests[i]]->setBlocked(false);
1580       }
1581     }
1582   }
1583
1584   void pup(PUP::er &p, AmpiRequestPool* reqPool) noexcept;
1585
1586   void print() const noexcept {
1587     for (int i=0; i<reqs.size(); i++) {
1588       if (reqs[i] == NULL) continue;
1589       CkPrintf("AmpiRequestList Element %d [%p]: \n", i+1, reqs[i]);
1590       reqs[i]->print();
1591     }
1592   }
1593 };
1594
1595 //A simple memory buffer
1596 class memBuf {
1597   CkVec<char> buf;
1598  public:
1599   memBuf() =default;
1600   memBuf(int size) noexcept : buf(size) {}
1601   void setSize(int s) noexcept {buf.resize(s);}
1602   int getSize() const noexcept {return buf.size();}
1603   const void *getData() const noexcept {return (const void *)&buf[0];}
1604   void *getData() noexcept {return (void *)&buf[0];}
1605 };
1606
1607 template <class T>
1608 inline void pupIntoBuf(memBuf &b,T &t) noexcept {
1609   PUP::sizer ps;ps|t;
1610   b.setSize(ps.size());
1611   PUP::toMem pm(b.getData()); pm|t;
1612 }
1613
1614 template <class T>
1615 inline void pupFromBuf(const void *data,T &t) noexcept {
1616   PUP::fromMem p(data); p|t;
1617 }
1618
1619 #define COLL_SEQ_IDX      -1
1620
1621 class AmpiMsg final : public CMessage_AmpiMsg {
1622  private:
1623   int ssendReq; //Index to the sender's request
1624   int tag; //MPI tag
1625   int srcRank; //Communicator rank for source
1626   int length; //Number of bytes in this message
1627   int origLength; // true size of allocation
1628   MPI_Comm comm; // Communicator
1629  public:
1630   char *data; //Payload
1631 #if CMK_BIGSIM_CHARM
1632  public:
1633   void *event;
1634   int  eventPe; // the PE that the event is located
1635 #endif
1636
1637  public:
1638   AmpiMsg() noexcept { data = NULL; }
1639   AmpiMsg(int sreq, int t, int sRank, int l) noexcept :
1640     ssendReq(sreq), tag(t), srcRank(sRank), length(l), origLength(l)
1641   { /* only called from AmpiMsg::pup() since the refnum (seq) will get pup'ed by the runtime */ }
1642   AmpiMsg(CMK_REFNUM_TYPE seq, int sreq, int t, int sRank, int l) noexcept :
1643     ssendReq(sreq), tag(t), srcRank(sRank), length(l), origLength(l)
1644   { CkSetRefNum(this, seq); }
1645   inline void setSsendReq(int s) noexcept { CkAssert(s >= 0); ssendReq = s; }
1646   inline void setSeq(CMK_REFNUM_TYPE s) noexcept { CkAssert(s >= 0); UsrToEnv(this)->setRef(s); }
1647   inline void setSrcRank(int sr) noexcept { srcRank = sr; }
1648   inline void setLength(int l) noexcept { length = l; }
1649   inline void setTag(int t) noexcept { tag = t; }
1650   inline void setComm(MPI_Comm c) noexcept { comm = c; }
1651   inline CMK_REFNUM_TYPE getSeq() const noexcept { return UsrToEnv(this)->getRef(); }
1652   inline int getSsendReq() const noexcept { return ssendReq; }
1653   inline int getSeqIdx() const noexcept {
1654     // seqIdx is srcRank, unless this message was part of a collective
1655     if (tag >= MPI_BCAST_TAG && tag <= MPI_ATA_TAG) {
1656       return COLL_SEQ_IDX;
1657     }
1658     else {
1659       return srcRank;
1660     }
1661   }
1662   inline int getSrcRank() const noexcept { return srcRank; }
1663   inline int getLength() const noexcept { return length; }
1664   inline char* getData() const noexcept { return data; }
1665   inline int getTag() const noexcept { return tag; }
1666   inline MPI_Comm getComm() const noexcept { return comm; }
1667   static AmpiMsg* pup(PUP::er &p, AmpiMsg *m) noexcept
1668   {
1669     int ref, ssendReq, tag, srcRank, length, origLength;
1670     MPI_Comm comm;
1671     if(p.isPacking() || p.isSizing()) {
1672       ref = CkGetRefNum(m);
1673       ssendReq = m->ssendReq;
1674       tag = m->tag;
1675       srcRank = m->srcRank;
1676       length = m->length;
1677       origLength = m->origLength;
1678       comm = m->comm;
1679     }
1680     p(ref); p(ssendReq); p(tag); p(srcRank); p(length); p(origLength); p(comm);
1681     if(p.isUnpacking()) {
1682       m = new (origLength, 0) AmpiMsg(ref, ssendReq, tag, srcRank, origLength);
1683       m->setLength(length);
1684       m->setComm(comm);
1685     }
1686     p(m->data, length);
1687     if(p.isDeleting()) {
1688       delete m;
1689       m = 0;
1690     }
1691     return m;
1692   }
1693 };
1694
1695 #define AMPI_MSG_POOL_SIZE   32 // Max # of AmpiMsg's allowed in the pool
1696 #define AMPI_POOLED_MSG_SIZE 64 // Max # of Bytes in pooled msgs' payload
1697
1698 class AmpiMsgPool {
1699  private:
1700   std::forward_list<AmpiMsg *> msgs; // list of free msgs
1701   int msgLength; // AmpiMsg::length of messages in the pool
1702   int msgUsersize; // usersize of message envelopes in the pool
1703   int maxMsgs; // max # of msgs in the pool
1704   int currMsgs; // current # of msgs in the pool
1705
1706  public:
1707   AmpiMsgPool() noexcept : msgLength(0), msgUsersize(0), maxMsgs(0), currMsgs(0) {}
1708   AmpiMsgPool(int _numMsgs, int _msgLength) noexcept {
1709     msgLength = _msgLength;
1710     maxMsgs = _numMsgs;
1711     if (maxMsgs > 0 && msgLength > 0) {
1712       /* Construct an AmpiMsg to find the usersize (and add it to the pool while it's here).
1713        * The rest of the pool can be filled lazily. */
1714       AmpiMsg* msg = new (msgLength, 0) AmpiMsg(0, 0, 0, 0, msgLength);
1715       msgs.push_front(msg);
1716       currMsgs = 1;
1717       /* Usersize is the true size of the message envelope, not the length member
1718        * of the AmpiMsg. AmpiMsg::length is used by Ssend msgs to convey the real
1719        * msg payload's length, and is not the length of the Ssend msg itself, so
1720        * it cannot be trusted when returning msgs to the pool. */
1721       msgUsersize = UsrToEnv(msgs.front())->getUsersize();
1722     }
1723     else {
1724       currMsgs = 0;
1725       msgUsersize = 0;
1726     }
1727   }
1728   ~AmpiMsgPool() =default;
1729   inline void clear() noexcept {
1730     while (!msgs.empty()) {
1731       delete msgs.front();
1732       msgs.pop_front();
1733     }
1734     currMsgs = 0;
1735   }
1736   inline AmpiMsg* newAmpiMsg(CMK_REFNUM_TYPE seq, int ssendReq, int tag, int srcRank, int len) noexcept {
1737     if (len > msgLength || msgs.empty()) {
1738       return new (len, 0) AmpiMsg(seq, ssendReq, tag, srcRank, len);
1739     } else {
1740       AmpiMsg* msg = msgs.front();
1741       CkAssert(msg != NULL);
1742       msgs.pop_front();
1743       currMsgs--;
1744       msg->setSeq(seq);
1745       msg->setSsendReq(ssendReq);
1746       msg->setTag(tag);
1747       msg->setSrcRank(srcRank);
1748       msg->setLength(len);
1749       return msg;
1750     }
1751   }
1752   inline void deleteAmpiMsg(AmpiMsg* msg) noexcept {
1753     if (currMsgs != maxMsgs && UsrToEnv(msg)->getUsersize() == msgUsersize) {
1754       CkAssert(msg != NULL);
1755       msgs.push_front(msg);
1756       currMsgs++;
1757     } else {
1758       delete msg;
1759     }
1760   }
1761   void pup(PUP::er& p) {
1762     p|msgLength;
1763     p|msgUsersize;
1764     p|maxMsgs;
1765     // Don't PUP the msgs in the free list or currMsgs, let the pool fill lazily
1766   }
1767 };
1768
1769 // Number of requests in the pool
1770 #ifndef AMPI_REQ_POOL_SIZE
1771 #define AMPI_REQ_POOL_SIZE 64
1772 #endif
1773
1774 // Helper macro for pool size and alignment calculations
1775 #define DefinePooledReqX(name, func) \
1776 static const size_t ireq##name = func(IReq); \
1777 static const size_t sreq##name = func(SendReq); \
1778 static const size_t ssreq##name = func(SsendReq); \
1779 static const size_t pooledReq##name = (ireq##name >= sreq##name && ireq##name >= ssreq##name) ? ireq##name : \
1780                                       (sreq##name >= ireq##name && sreq##name >= ssreq##name) ? sreq##name : \
1781                                       (ssreq##name);
1782
1783 // This defines 'static const size_t pooledReqSize = ... ;'
1784 DefinePooledReqX(Size, sizeof)
1785
1786 // This defines 'static const size_t pooledReqAlign = ... ;'
1787 DefinePooledReqX(Align, alignof)
1788
1789 // Pool of IReq, SendReq, and SsendReq objects:
1790 // These are different sizes, but we use a single pool for them so
1791 // that iteration over these objects is fast, as in AMPI_Waitall.
1792 // We also try to always allocate new requests from the start to the end
1793 // of the pool, so that forward iteration over requests is fast.
1794 class AmpiRequestPool {
1795  private:
1796   std::bitset<AMPI_REQ_POOL_SIZE> validReqs; // reqs in the pool are either valid (being used by a real req) or invalid
1797   int startIdx = 0; // start next search from this index
1798   alignas(pooledReqAlign) std::array<char, AMPI_REQ_POOL_SIZE*pooledReqSize> reqs; // pool of memory for requests
1799
1800  public:
1801   AmpiRequestPool() =default;
1802   ~AmpiRequestPool() =default;
1803   template <typename T, typename... Args>
1804   inline T* newReq(Args&&... args) noexcept {
1805     if (validReqs.all()) {
1806       return new T(std::forward<Args>(args)...);
1807     } else {
1808       for (int i=startIdx; i<validReqs.size(); i++) {
1809         if (!validReqs[i]) {
1810           validReqs[i] = 1;
1811           startIdx = i+1;
1812           T* req = new (&reqs[i*pooledReqSize]) T(std::forward<Args>(args)...);
1813           return req;
1814         }
1815       }
1816       CkAbort("AMPI> failed to find a free request in pool!");
1817       return NULL;
1818     }
1819   }
1820   inline void deleteReq(AmpiRequest* req) noexcept {
1821     if (req->isPooledType() &&
1822         ((char*)req >= &reqs.front() && (char*)req <= &reqs.back()))
1823     {
1824       int idx = (int)((intptr_t)req - (intptr_t)&reqs[0]) / pooledReqSize;
1825       validReqs[idx] = 0;
1826       startIdx = std::min(idx, startIdx);
1827     } else {
1828       delete req;
1829     }
1830   }
1831   void pup(PUP::er& p) noexcept {
1832     // Nothing to do here, because AmpiRequestList::pup will be the
1833     // one to actually PUP the AmpiRequest objects to/from the pool
1834   }
1835 };
1836
1837 /**
1838   Our local representation of another AMPI
1839  array element.  Used to keep track of incoming
1840  and outgoing message sequence numbers, and
1841  the out-of-order message list.
1842 */
1843 class AmpiOtherElement {
1844 private:
1845   /// Next incoming and outgoing message sequence number
1846   CMK_REFNUM_TYPE seqIncoming, seqOutgoing;
1847
1848   /// Number of messages in out-of-order queue (normally 0)
1849   uint16_t numOutOfOrder;
1850
1851 public:
1852   /// seqIncoming starts from 1, b/c 0 means unsequenced
1853   /// seqOutgoing starts from 0, b/c this will be incremented for the first real seq #
1854   AmpiOtherElement() noexcept : seqIncoming(1), seqOutgoing(0), numOutOfOrder(0) {}
1855
1856   /// Handle wrap around of unsigned type CMK_REFNUM_TYPE
1857   inline void incSeqIncoming() noexcept { seqIncoming++; if (seqIncoming==0) seqIncoming=1; }
1858   inline CMK_REFNUM_TYPE getSeqIncoming() const noexcept { return seqIncoming; }
1859
1860   inline void incSeqOutgoing() noexcept { seqOutgoing++; if (seqOutgoing==0) seqOutgoing=1; }
1861   inline CMK_REFNUM_TYPE getSeqOutgoing() const noexcept { return seqOutgoing; }
1862
1863   inline void incNumOutOfOrder() noexcept { numOutOfOrder++; }
1864   inline void decNumOutOfOrder() noexcept { numOutOfOrder--; }
1865   inline uint16_t getNumOutOfOrder() const noexcept { return numOutOfOrder; }
1866 };
1867 PUPbytes(AmpiOtherElement)
1868
1869 class AmpiSeqQ : private CkNoncopyable {
1870   CkMsgQ<AmpiMsg> out; // all out of order messages
1871   std::unordered_map<int, AmpiOtherElement> elements; // element info: indexed by seqIdx (comm rank)
1872
1873 public:
1874   AmpiSeqQ() =default;
1875   AmpiSeqQ(int commSize) noexcept {
1876     elements.reserve(std::min(commSize, 64));
1877   }
1878   ~AmpiSeqQ() =default;
1879   void pup(PUP::er &p) noexcept;
1880
1881   /// Insert this message in the table.  Returns the number
1882   /// of messages now available for the element.
1883   ///   If 0, the message was out-of-order and is buffered.
1884   ///   If 1, this message can be immediately processed.
1885   ///   If >1, this message can be immediately processed,
1886   ///     and you should call "getOutOfOrder" repeatedly.
1887   inline int put(int seqIdx, AmpiMsg *msg) noexcept {
1888     AmpiOtherElement &el = elements[seqIdx];
1889     if (msg->getSeq() == el.getSeqIncoming()) { // In order:
1890       el.incSeqIncoming();
1891       return 1+el.getNumOutOfOrder();
1892     }
1893     else { // Out of order: stash message
1894       putOutOfOrder(seqIdx, msg);
1895       return 0;
1896     }
1897   }
1898
1899   /// Is this message in order (return >0) or not (return 0)?
1900   /// Same as put() except we don't call putOutOfOrder() here,
1901   /// so the caller should do that separately
1902   inline int isInOrder(int srcRank, CMK_REFNUM_TYPE seq) noexcept {
1903     AmpiOtherElement &el = elements[srcRank];
1904     if (seq == el.getSeqIncoming()) { // In order:
1905       el.incSeqIncoming();
1906       return 1+el.getNumOutOfOrder();
1907     }
1908     else { // Out of order: caller should stash message
1909       return 0;
1910     }
1911   }
1912
1913   /// Get an out-of-order message from the table.
1914   /// (in-order messages never go into the table)
1915   AmpiMsg *getOutOfOrder(int seqIdx) noexcept;
1916
1917   /// Stash an out-of-order message
1918   void putOutOfOrder(int seqIdx, AmpiMsg *msg) noexcept;
1919
1920   /// Increment the outgoing sequence number.
1921   inline void incCollSeqOutgoing() noexcept {
1922     elements[COLL_SEQ_IDX].incSeqOutgoing();
1923   }
1924
1925   /// Return the next outgoing sequence number, and increment it.
1926   inline CMK_REFNUM_TYPE nextOutgoing(int destRank) noexcept {
1927     AmpiOtherElement &el = elements[destRank];
1928     el.incSeqOutgoing();
1929     return el.getSeqOutgoing();
1930   }
1931 };
1932 PUPmarshall(AmpiSeqQ)
1933
1934
1935 inline CProxy_ampi ampiCommStruct::getProxy() const noexcept {return ampiID;}
1936 const ampiCommStruct &universeComm2CommStruct(MPI_Comm universeNo) noexcept;
1937
1938 // Max value of a predefined MPI_Op (values defined in ampi.h)
1939 #define AMPI_MAX_PREDEFINED_OP 13
1940
1941 /*
1942 An ampiParent holds all the communicators and the TCharm thread
1943 for its children, which are bound to it.
1944 */
1945 class ampiParent final : public CBase_ampiParent {
1946  private:
1947   TCharm *thread;
1948   CProxy_TCharm threads;
1949
1950  public: // Communication state:
1951   int numBlockedReqs; // number of requests currently blocked on
1952   bool resumeOnRecv, resumeOnColl;
1953   AmpiRequestList ampiReqs;
1954   AmpiRequestPool reqPool;
1955   AmpiRequest *blockingReq;
1956   CkDDT myDDT;
1957
1958  private:
1959   MPI_Comm worldNo; //My MPI_COMM_WORLD
1960   ampi *worldPtr; //AMPI element corresponding to MPI_COMM_WORLD
1961
1962   CkPupPtrVec<ampiCommStruct> splitComm;     //Communicators from MPI_Comm_split
1963   CkPupPtrVec<ampiCommStruct> groupComm;     //Communicators from MPI_Comm_group
1964   CkPupPtrVec<ampiCommStruct> cartComm;      //Communicators from MPI_Cart_create
1965   CkPupPtrVec<ampiCommStruct> graphComm;     //Communicators from MPI_Graph_create
1966   CkPupPtrVec<ampiCommStruct> distGraphComm; //Communicators from MPI_Dist_graph_create
1967   CkPupPtrVec<ampiCommStruct> interComm;     //Communicators from MPI_Intercomm_create
1968   CkPupPtrVec<ampiCommStruct> intraComm;     //Communicators from MPI_Intercomm_merge
1969
1970   CkPupPtrVec<groupStruct> groups; // "Wild" groups that don't have a communicator
1971   CkPupPtrVec<WinStruct> winStructList; //List of windows for one-sided communication
1972   CkPupPtrVec<InfoStruct> infos; // list of all MPI_Infos
1973   const std::array<MPI_User_function*, AMPI_MAX_PREDEFINED_OP+1>& predefinedOps; // owned by ampiNodeMgr
1974   vector<OpStruct> userOps; // list of any user-defined MPI_Ops
1975   vector<AmpiMsg *> matchedMsgs; // for use with MPI_Mprobe and MPI_Mrecv
1976
1977   /* MPI_*_get_attr C binding returns a *pointer* to an integer,
1978    *  so there needs to be some storage somewhere to point to.
1979    * All builtin keyvals are ints, except for MPI_WIN_BASE, which
1980    *  is a pointer, and MPI_WIN_SIZE, which is an MPI_Aint. */
1981   int* kv_builtin_storage;
1982   MPI_Aint* win_size_storage;
1983   void** win_base_storage;
1984   CkPupPtrVec<KeyvalNode> kvlist;
1985   void* bsendBuffer;   // NOTE: we don't actually use this for buffering of MPI_Bsend's,
1986   int bsendBufferSize; //       we only keep track of it to return it from MPI_Buffer_detach
1987
1988   // Intercommunicator creation:
1989   bool isTmpRProxySet;
1990   CProxy_ampi tmpRProxy;
1991
1992   MPI_MigrateFn userAboutToMigrateFn, userJustMigratedFn;
1993
1994  public:
1995   bool ampiInitCallDone;
1996
1997  private:
1998   bool kv_set_builtin(int keyval, void* attribute_val) noexcept;
1999   bool kv_get_builtin(int keyval) noexcept;
2000
2001  public:
2002   void prepareCtv() noexcept;
2003
2004   MPI_Message putMatchedMsg(AmpiMsg* msg) noexcept {
2005     // Search thru matchedMsgs for any NULL ones first:
2006     for (int i=0; i<matchedMsgs.size(); i++) {
2007       if (matchedMsgs[i] == NULL) {
2008         matchedMsgs[i] = msg;
2009         return i;
2010       }
2011     }
2012     // No NULL entries, so create a new one:
2013     matchedMsgs.push_back(msg);
2014     return matchedMsgs.size() - 1;
2015   }
2016   AmpiMsg* getMatchedMsg(MPI_Message message) noexcept {
2017     if (message == MPI_MESSAGE_NO_PROC || message == MPI_MESSAGE_NULL) {
2018       return NULL;
2019     }
2020     CkAssert(message >= 0 && message < matchedMsgs.size());
2021     AmpiMsg* msg = matchedMsgs[message];
2022     // Mark this matchedMsg index NULL and free from back of vector:
2023     matchedMsgs[message] = NULL;
2024     while (matchedMsgs.back() == NULL) {
2025       matchedMsgs.pop_back();
2026     }
2027     return msg;
2028   }
2029
2030   inline void attachBuffer(void *buffer, int size) noexcept {
2031     bsendBuffer = buffer;
2032     bsendBufferSize = size;
2033   }
2034   inline void detachBuffer(void *buffer, int *size) noexcept {
2035     *(void **)buffer = bsendBuffer;
2036     *size = bsendBufferSize;
2037   }
2038   inline bool isSplit(MPI_Comm comm) const noexcept {
2039     return (comm>=MPI_COMM_FIRST_SPLIT && comm<MPI_COMM_FIRST_GROUP);
2040   }
2041   const ampiCommStruct &getSplit(MPI_Comm comm) const noexcept {
2042     int idx=comm-MPI_COMM_FIRST_SPLIT;
2043     if (idx>=splitComm.size()) CkAbort("Bad split communicator used");
2044     return *splitComm[idx];
2045   }
2046   void splitChildRegister(const ampiCommStruct &s) noexcept;
2047
2048   inline bool isGroup(MPI_Comm comm) const noexcept {
2049     return (comm>=MPI_COMM_FIRST_GROUP && comm<MPI_COMM_FIRST_CART);
2050   }
2051   const ampiCommStruct &getGroup(MPI_Comm comm) const noexcept {
2052     int idx=comm-MPI_COMM_FIRST_GROUP;
2053     if (idx>=groupComm.size()) CkAbort("Bad group communicator used");
2054     return *groupComm[idx];
2055   }
2056   void groupChildRegister(const ampiCommStruct &s) noexcept;
2057   inline bool isInGroups(MPI_Group group) const noexcept {
2058     return (group>=0 && group<groups.size());
2059   }
2060
2061   void cartChildRegister(const ampiCommStruct &s) noexcept;
2062   void graphChildRegister(const ampiCommStruct &s) noexcept;
2063   void distGraphChildRegister(const ampiCommStruct &s) noexcept;
2064   void interChildRegister(const ampiCommStruct &s) noexcept;
2065   void intraChildRegister(const ampiCommStruct &s) noexcept;
2066
2067  public:
2068   ampiParent(MPI_Comm worldNo_,CProxy_TCharm threads_,int nRanks_) noexcept;
2069   ampiParent(CkMigrateMessage *msg) noexcept;
2070   void ckAboutToMigrate() noexcept;
2071   void ckJustMigrated() noexcept;
2072   void ckJustRestored() noexcept;
2073   void setUserAboutToMigrateFn(MPI_MigrateFn f) noexcept;
2074   void setUserJustMigratedFn(MPI_MigrateFn f) noexcept;
2075   ~ampiParent() noexcept;
2076
2077   //Children call this when they are first created, or just migrated
2078   TCharm *registerAmpi(ampi *ptr,ampiCommStruct s,bool forMigration) noexcept;
2079
2080   // exchange proxy info between two ampi proxies
2081   void ExchangeProxy(CProxy_ampi rproxy) noexcept {
2082     if(!isTmpRProxySet){ tmpRProxy=rproxy; isTmpRProxySet=true; }
2083     else{ tmpRProxy.setRemoteProxy(rproxy); rproxy.setRemoteProxy(tmpRProxy); isTmpRProxySet=false; }
2084   }
2085
2086   //Grab the next available split/group communicator
2087   MPI_Comm getNextSplit() const noexcept {return MPI_COMM_FIRST_SPLIT+splitComm.size();}
2088   MPI_Comm getNextGroup() const noexcept {return MPI_COMM_FIRST_GROUP+groupComm.size();}
2089   MPI_Comm getNextCart() const noexcept {return MPI_COMM_FIRST_CART+cartComm.size();}
2090   MPI_Comm getNextGraph() const noexcept {return MPI_COMM_FIRST_GRAPH+graphComm.size();}
2091   MPI_Comm getNextDistGraph() const noexcept {return MPI_COMM_FIRST_DIST_GRAPH+distGraphComm.size();}
2092   MPI_Comm getNextInter() const noexcept {return MPI_COMM_FIRST_INTER+interComm.size();}
2093   MPI_Comm getNextIntra() const noexcept {return MPI_COMM_FIRST_INTRA+intraComm.size();}
2094
2095   inline bool isCart(MPI_Comm comm) const noexcept {
2096     return (comm>=MPI_COMM_FIRST_CART && comm<MPI_COMM_FIRST_GRAPH);
2097   }
2098   ampiCommStruct &getCart(MPI_Comm comm) const noexcept {
2099     int idx=comm-MPI_COMM_FIRST_CART;
2100     if (idx>=cartComm.size()) CkAbort("AMPI> Bad cartesian communicator used!\n");
2101     return *cartComm[idx];
2102   }
2103   inline bool isGraph(MPI_Comm comm) const noexcept {
2104     return (comm>=MPI_COMM_FIRST_GRAPH && comm<MPI_COMM_FIRST_DIST_GRAPH);
2105   }
2106   ampiCommStruct &getGraph(MPI_Comm comm) const noexcept {
2107     int idx=comm-MPI_COMM_FIRST_GRAPH;
2108     if (idx>=graphComm.size()) CkAbort("AMPI> Bad graph communicator used!\n");
2109     return *graphComm[idx];
2110   }
2111   inline bool isDistGraph(MPI_Comm comm) const noexcept {
2112     return (comm >= MPI_COMM_FIRST_DIST_GRAPH && comm < MPI_COMM_FIRST_INTER);
2113   }
2114   ampiCommStruct &getDistGraph(MPI_Comm comm) const noexcept {
2115     int idx = comm-MPI_COMM_FIRST_DIST_GRAPH;
2116     if (idx>=distGraphComm.size()) CkAbort("Bad distributed graph communicator used");
2117     return *distGraphComm[idx];
2118   }
2119   inline bool isInter(MPI_Comm comm) const noexcept {
2120     return (comm>=MPI_COMM_FIRST_INTER && comm<MPI_COMM_FIRST_INTRA);
2121   }
2122   const ampiCommStruct &getInter(MPI_Comm comm) const noexcept {
2123     int idx=comm-MPI_COMM_FIRST_INTER;
2124     if (idx>=interComm.size()) CkAbort("AMPI> Bad inter-communicator used!\n");
2125     return *interComm[idx];
2126   }
2127   inline bool isIntra(MPI_Comm comm) const noexcept {
2128     return (comm>=MPI_COMM_FIRST_INTRA && comm<MPI_COMM_FIRST_RESVD);
2129   }
2130   const ampiCommStruct &getIntra(MPI_Comm comm) const noexcept {
2131     int idx=comm-MPI_COMM_FIRST_INTRA;
2132     if (idx>=intraComm.size()) CkAbort("Bad intra-communicator used");
2133     return *intraComm[idx];
2134   }
2135
2136   void pup(PUP::er &p) noexcept;
2137
2138   void startCheckpoint(const char* dname) noexcept;
2139   void Checkpoint(int len, const char* dname) noexcept;
2140   void ResumeThread() noexcept;
2141   TCharm* getTCharmThread() const noexcept {return thread;}
2142   inline ampiParent* blockOnRecv() noexcept;
2143   inline CkDDT* getDDT() noexcept { return &myDDT; }
2144
2145 #if CMK_LBDB_ON
2146   void setMigratable(bool mig) noexcept {
2147     thread->setMigratable(mig);
2148   }
2149 #endif
2150
2151   const ampiCommStruct &getWorldStruct() const noexcept;
2152
2153   inline const ampiCommStruct &comm2CommStruct(MPI_Comm comm) const noexcept {
2154     if (comm==MPI_COMM_WORLD) return getWorldStruct();
2155     if (comm==worldNo) return getWorldStruct();
2156     if (isSplit(comm)) return getSplit(comm);
2157     if (isGroup(comm)) return getGroup(comm);
2158     if (isCart(comm)) return getCart(comm);
2159     if (isGraph(comm)) return getGraph(comm);
2160     if (isDistGraph(comm)) return getDistGraph(comm);
2161     if (isInter(comm)) return getInter(comm);
2162     if (isIntra(comm)) return getIntra(comm);
2163     return universeComm2CommStruct(comm);
2164   }
2165
2166   inline vector<int>& getKeyvals(MPI_Comm comm) noexcept {
2167     ampiCommStruct &cs = *(ampiCommStruct *)&comm2CommStruct(comm);
2168     return cs.getKeyvals();
2169   }
2170
2171   inline ampi *comm2ampi(MPI_Comm comm) const noexcept {
2172     if (comm==MPI_COMM_WORLD) return worldPtr;
2173     if (comm==worldNo) return worldPtr;
2174     if (isSplit(comm)) {
2175       const ampiCommStruct &st=getSplit(comm);
2176       return st.getProxy()[thisIndex].ckLocal();
2177     }
2178     if (isGroup(comm)) {
2179       const ampiCommStruct &st=getGroup(comm);
2180       return st.getProxy()[thisIndex].ckLocal();
2181     }
2182     if (isCart(comm)) {
2183       const ampiCommStruct &st = getCart(comm);
2184       return st.getProxy()[thisIndex].ckLocal();
2185     }
2186     if (isGraph(comm)) {
2187       const ampiCommStruct &st = getGraph(comm);
2188       return st.getProxy()[thisIndex].ckLocal();
2189     }
2190     if (isDistGraph(comm)) {
2191       const ampiCommStruct &st = getDistGraph(comm);
2192       return st.getProxy()[thisIndex].ckLocal();
2193     }
2194     if (isInter(comm)) {
2195       const ampiCommStruct &st=getInter(comm);
2196       return st.getProxy()[thisIndex].ckLocal();
2197     }
2198     if (isIntra(comm)) {
2199       const ampiCommStruct &st=getIntra(comm);
2200       return st.getProxy()[thisIndex].ckLocal();
2201     }
2202     if (comm>MPI_COMM_WORLD) return worldPtr; //Use MPI_WORLD ampi for cross-world messages:
2203     CkAbort("Invalid communicator used!");
2204     return NULL;
2205   }
2206
2207   inline bool hasComm(const MPI_Group group) const noexcept {
2208     MPI_Comm comm = (MPI_Comm)group;
2209     return ( comm==MPI_COMM_WORLD || comm==worldNo || isSplit(comm) || isGroup(comm) ||
2210              isCart(comm) || isGraph(comm) || isDistGraph(comm) || isIntra(comm) );
2211     //isInter omitted because its comm number != its group number
2212   }
2213   inline vector<int> group2vec(MPI_Group group) const noexcept {
2214     if (group == MPI_GROUP_NULL || group == MPI_GROUP_EMPTY) {
2215       return vector<int>();
2216     }
2217     else if (hasComm(group)) {
2218       return comm2CommStruct((MPI_Comm)group).getIndices();
2219     }
2220     else {
2221       CkAssert(isInGroups(group));
2222       return groups[group]->getRanks();
2223     }
2224   }
2225   inline MPI_Group saveGroupStruct(const vector<int>& vec) noexcept {
2226     if (vec.empty()) return MPI_GROUP_EMPTY;
2227     int idx = groups.size();
2228     groups.resize(idx+1);
2229     groups[idx]=new groupStruct(vec);
2230     return (MPI_Group)idx;
2231   }
2232   inline int getRank(const MPI_Group group) const noexcept {
2233     vector<int> vec = group2vec(group);
2234     return getPosOp(thisIndex,vec);
2235   }
2236   inline AmpiRequestList &getReqs() noexcept { return ampiReqs; }
2237   inline int getMyPe() const noexcept {
2238     return CkMyPe();
2239   }
2240   inline bool hasWorld() const noexcept {
2241     return worldPtr!=NULL;
2242   }
2243
2244   inline void checkComm(MPI_Comm comm) const noexcept {
2245     if ((comm != MPI_COMM_SELF && comm != MPI_COMM_WORLD)
2246      || (isSplit(comm) && comm-MPI_COMM_FIRST_SPLIT >= splitComm.size())
2247      || (isGroup(comm) && comm-MPI_COMM_FIRST_GROUP >= groupComm.size())
2248      || (isCart(comm)  && comm-MPI_COMM_FIRST_CART  >=  cartComm.size())
2249      || (isGraph(comm) && comm-MPI_COMM_FIRST_GRAPH >= graphComm.size())
2250      || (isDistGraph(comm) && comm-MPI_COMM_FIRST_DIST_GRAPH >= distGraphComm.size())
2251      || (isInter(comm) && comm-MPI_COMM_FIRST_INTER >= interComm.size())
2252      || (isIntra(comm) && comm-MPI_COMM_FIRST_INTRA >= intraComm.size()) )
2253       CkAbort("Invalid MPI_Comm\n");
2254   }
2255
2256   /// if intra-communicator, return comm, otherwise return null group
2257   inline MPI_Group comm2group(const MPI_Comm comm) const noexcept {
2258     if(isInter(comm)) return MPI_GROUP_NULL;   // we don't support inter-communicator in such functions
2259     ampiCommStruct s = comm2CommStruct(comm);
2260     if(comm!=MPI_COMM_WORLD && comm!=s.getComm()) CkAbort("Error in ampiParent::comm2group()");
2261     return (MPI_Group)(s.getComm());
2262   }
2263
2264   inline int getRemoteSize(const MPI_Comm comm) const noexcept {
2265     if(isInter(comm)) return getInter(comm).getRemoteIndices().size();
2266     else return -1;
2267   }
2268   inline MPI_Group getRemoteGroup(const MPI_Comm comm) noexcept {
2269     if(isInter(comm)) return saveGroupStruct(getInter(comm).getRemoteIndices());
2270     else return MPI_GROUP_NULL;
2271   }
2272
2273   int createKeyval(MPI_Copy_function *copy_fn, MPI_Delete_function *delete_fn,
2274                   int *keyval, void* extra_state) noexcept;
2275   bool getBuiltinKeyval(int keyval, void *attribute_val) noexcept;
2276   int setUserKeyval(MPI_Comm comm, int keyval, void *attribute_val) noexcept;
2277   bool getUserKeyval(MPI_Comm comm, vector<int>& keyvals, int keyval, void *attribute_val, int *flag) noexcept;
2278   int dupUserKeyvals(MPI_Comm old_comm, MPI_Comm new_comm) noexcept;
2279   int freeUserKeyval(int context, vector<int>& keyvals, int *keyval) noexcept;
2280   int freeUserKeyvals(int context, vector<int>& keyvals) noexcept;
2281
2282   int setAttr(MPI_Comm comm, vector<int>& keyvals, int keyval, void *attribute_val) noexcept;
2283   int getAttr(MPI_Comm comm, vector<int>& keyvals, int keyval, void *attribute_val, int *flag) noexcept;
2284   int deleteAttr(MPI_Comm comm, vector<int>& keyvals, int keyval) noexcept;
2285
2286   int addWinStruct(WinStruct *win) noexcept;
2287   WinStruct *getWinStruct(MPI_Win win) const noexcept;
2288   void removeWinStruct(WinStruct *win) noexcept;
2289
2290   int createInfo(MPI_Info *newinfo) noexcept;
2291   int dupInfo(MPI_Info info, MPI_Info *newinfo) noexcept;
2292   int setInfo(MPI_Info info, const char *key, const char *value) noexcept;
2293   int deleteInfo(MPI_Info info, const char *key) noexcept;
2294   int getInfo(MPI_Info info, const char *key, int valuelen, char *value, int *flag) const noexcept;
2295   int getInfoValuelen(MPI_Info info, const char *key, int *valuelen, int *flag) const noexcept;
2296   int getInfoNkeys(MPI_Info info, int *nkeys) const noexcept;
2297   int getInfoNthkey(MPI_Info info, int n, char *key) const noexcept;
2298   int freeInfo(MPI_Info info) noexcept;
2299   void defineInfoEnv(int nRanks_) noexcept;
2300   void defineInfoMigration() noexcept;
2301
2302   // An 'MPI_Op' is an integer that indexes into either:
2303   //   A) an array of predefined ops owned by ampiNodeMgr, or
2304   //   B) a vector of user-defined ops owned by ampiParent
2305   // The MPI_Op is compared to AMPI_MAX_PREDEFINED_OP to disambiguate.
2306   inline int createOp(MPI_User_function *fn, bool isCommutative) noexcept {
2307     // Search thru non-predefined op's for any invalidated ones:
2308     for (int i=0; i<userOps.size(); i++) {
2309       if (userOps[i].isFree()) {
2310         userOps[i].init(fn, isCommutative);
2311         return AMPI_MAX_PREDEFINED_OP + 1 + i;
2312       }
2313     }
2314     // No invalid entries, so create a new one:
2315     userOps.emplace_back(fn, isCommutative);
2316     return AMPI_MAX_PREDEFINED_OP + userOps.size();
2317   }
2318   inline void freeOp(MPI_Op op) noexcept {
2319     // Don't free predefined op's:
2320     if (!opIsPredefined(op)) {
2321       // Invalidate op, then free all invalid op's from the back of the userOp's vector
2322       int opIdx = op - 1 - AMPI_MAX_PREDEFINED_OP;
2323       CkAssert(opIdx < userOps.size());
2324       userOps[opIdx].free();
2325       while (!userOps.empty() && userOps.back().isFree()) {
2326         userOps.pop_back();
2327       }
2328     }
2329   }
2330   inline bool opIsPredefined(MPI_Op op) const noexcept {
2331     return (op <= AMPI_MAX_PREDEFINED_OP);
2332   }
2333   inline bool opIsCommutative(MPI_Op op) const noexcept {
2334     if (opIsPredefined(op)) {
2335       return true; // all predefined ops are commutative
2336     }
2337     else {
2338       int opIdx = op - 1 - AMPI_MAX_PREDEFINED_OP;
2339       CkAssert(opIdx < userOps.size());
2340       return userOps[opIdx].isCommutative;
2341     }
2342   }
2343   inline MPI_User_function* op2User_function(MPI_Op op) const noexcept {
2344     if (opIsPredefined(op)) {
2345       return predefinedOps[op];
2346     }
2347     else {
2348       int opIdx = op - 1 - AMPI_MAX_PREDEFINED_OP;
2349       CkAssert(opIdx < userOps.size());
2350       return userOps[opIdx].func;
2351     }
2352   }
2353   inline AmpiOpHeader op2AmpiOpHeader(MPI_Op op, MPI_Datatype type, int count) const noexcept {
2354     if (opIsPredefined(op)) {
2355       int size = myDDT.getType(type)->getSize(count);
2356       return AmpiOpHeader(predefinedOps[op], type, count, size);
2357     }
2358     else {
2359       int opIdx = op - 1 - AMPI_MAX_PREDEFINED_OP;
2360       CkAssert(opIdx < userOps.size());
2361       int size = myDDT.getType(type)->getSize(count);
2362       return AmpiOpHeader(userOps[opIdx].func, type, count, size);
2363     }
2364   }
2365   inline void applyOp(MPI_Datatype datatype, MPI_Op op, int count, const void* invec, void* inoutvec) const noexcept {
2366     // inoutvec[i] = invec[i] op inoutvec[i]
2367     MPI_User_function *func = op2User_function(op);
2368     (func)((void*)invec, inoutvec, &count, &datatype);
2369   }
2370
2371   void init() noexcept;
2372   void finalize() noexcept;
2373   void block() noexcept;
2374   void yield() noexcept;
2375
2376 #if AMPI_PRINT_MSG_SIZES
2377 // Map of AMPI routine names to message sizes and number of messages:
2378 // ["AMPI_Routine"][ [msg_size][num_msgs] ]
2379   std::unordered_map<std::string, std::map<int, int> > msgSizes;
2380   inline bool isRankRecordingMsgSizes() noexcept;
2381   inline void recordMsgSize(const char* func, int msgSize) noexcept;
2382   void printMsgSizes() noexcept;
2383 #endif
2384
2385 #if AMPIMSGLOG
2386   /* message logging */
2387   int pupBytes;
2388 #if CMK_USE_ZLIB && 0
2389   gzFile fMsgLog;
2390   PUP::tozDisk *toPUPer;
2391   PUP::fromzDisk *fromPUPer;
2392 #else
2393   FILE* fMsgLog;
2394   PUP::toDisk *toPUPer;
2395   PUP::fromDisk *fromPUPer;
2396 #endif
2397 #endif
2398 };
2399
2400 // Store a generalized request class created by MPIX_Grequest_class_create
2401 class greq_class_desc {
2402 public:
2403   MPI_Grequest_query_function *query_fn;
2404   MPI_Grequest_free_function *free_fn;
2405   MPI_Grequest_cancel_function *cancel_fn;
2406   MPIX_Grequest_poll_function *poll_fn;
2407   MPIX_Grequest_wait_function *wait_fn;
2408
2409   void pup(PUP::er &p) noexcept {
2410     p((char *)query_fn, sizeof(void *));
2411     p((char *)free_fn, sizeof(void *));
2412     p((char *)cancel_fn, sizeof(void *));
2413     p((char *)poll_fn, sizeof(void *));
2414     p((char *)wait_fn, sizeof(void *));
2415   }
2416 };
2417
2418 /*
2419 An ampi manages the communication of one thread over
2420 one MPI communicator.
2421 */
2422 class ampi final : public CBase_ampi {
2423  private:
2424   friend class IReq; // for checking resumeOnRecv
2425   friend class SendReq;
2426   friend class SsendReq;
2427   friend class RednReq;
2428   friend class GatherReq;
2429   friend class GathervReq;
2430
2431   ampiParent *parent;
2432   CProxy_ampiParent parentProxy;
2433   TCharm *thread;
2434   int myRank;
2435   AmpiSeqQ oorder;
2436
2437  public:
2438   /*
2439    * AMPI Message Matching (Amm) queues are indexed by the tag and sender.
2440    * Since ampi objects are per-communicator, there are separate Amm's per communicator.
2441    */
2442   Amm<AmpiRequest *, AMPI_AMM_PT2PT_POOL_SIZE> postedReqs;
2443   Amm<AmpiMsg *, AMPI_AMM_PT2PT_POOL_SIZE> unexpectedMsgs;
2444
2445   // Bcast requests / msgs must be kept separate from pt2pt,
2446   // so we don't match them to wildcard recv's
2447   Amm<AmpiRequest *, AMPI_AMM_COLL_POOL_SIZE> postedBcastReqs;
2448   Amm<AmpiMsg *, AMPI_AMM_COLL_POOL_SIZE> unexpectedBcastMsgs;
2449
2450   // Store generalized request classes created by MPIX_Grequest_class_create
2451   vector<greq_class_desc> greq_classes;
2452
2453  private:
2454   ampiCommStruct myComm;
2455   vector<int> tmpVec; // stores temp group info
2456   CProxy_ampi remoteProxy; // valid only for intercommunicator
2457   CkPupPtrVec<win_obj> winObjects;
2458
2459  private:
2460   void inorder(AmpiMsg *msg) noexcept;
2461   void inorderBcast(AmpiMsg *msg, bool deleteMsg) noexcept;
2462   void inorderRdma(char* buf, int size, CMK_REFNUM_TYPE seq, int tag, int srcRank,
2463                    MPI_Comm comm, int ssendReq) noexcept;
2464
2465   void init() noexcept;
2466   void findParent(bool forMigration) noexcept;
2467
2468  public: // entry methods
2469   ampi() noexcept;
2470   ampi(CkArrayID parent_,const ampiCommStruct &s) noexcept;
2471   ampi(CkMigrateMessage *msg) noexcept;
2472   void ckJustMigrated() noexcept;
2473   void ckJustRestored() noexcept;
2474   ~ampi() noexcept;
2475
2476   void pup(PUP::er &p) noexcept;
2477
2478   void allInitDone() noexcept;
2479   void setInitDoneFlag() noexcept;
2480
2481   void unblock() noexcept;
2482   void injectMsg(int size, char* buf) noexcept;
2483   void generic(AmpiMsg *) noexcept;
2484   void genericRdma(char* buf, int size, CMK_REFNUM_TYPE seq, int tag, int srcRank,
2485                    MPI_Comm destcomm, int ssendReq) noexcept;
2486   void completedRdmaSend(CkDataMsg *msg) noexcept;
2487   void ssend_ack(int sreq) noexcept;
2488   void barrierResult() noexcept;
2489   void ibarrierResult() noexcept;
2490   void bcastResult(AmpiMsg *msg) noexcept;
2491   void rednResult(CkReductionMsg *msg) noexcept;
2492   void irednResult(CkReductionMsg *msg) noexcept;
2493
2494   void splitPhase1(CkReductionMsg *msg) noexcept;
2495   void splitPhaseInter(CkReductionMsg *msg) noexcept;
2496   void commCreatePhase1(MPI_Comm nextGroupComm) noexcept;
2497   void intercommCreatePhase1(MPI_Comm nextInterComm) noexcept;
2498   void intercommMergePhase1(MPI_Comm nextIntraComm) noexcept;
2499
2500  private: // Used by the above entry methods that create new MPI_Comm objects
2501   CProxy_ampi createNewChildAmpiSync() noexcept;
2502   void insertNewChildAmpiElements(MPI_Comm newComm, CProxy_ampi newAmpi) noexcept;
2503
2504   inline void handleBlockedReq(AmpiRequest* req) noexcept {
2505     if (req->isBlocked() && parent->numBlockedReqs != 0) {
2506       parent->numBlockedReqs--;
2507     }
2508   }
2509   inline void resumeThreadIfReady() noexcept {
2510     if (parent->resumeOnRecv && parent->numBlockedReqs == 0) {
2511       thread->resume();
2512     }
2513   }
2514
2515  public: // to be used by MPI_* functions
2516   inline const ampiCommStruct &comm2CommStruct(MPI_Comm comm) const noexcept {
2517     return parent->comm2CommStruct(comm);
2518   }
2519   inline const ampiCommStruct &getCommStruct() const noexcept { return myComm; }
2520
2521   inline ampi* blockOnRecv() noexcept;
2522   inline ampi* blockOnColl() noexcept;
2523   inline void setBlockingReq(AmpiRequest *req) noexcept;
2524   MPI_Request postReq(AmpiRequest* newreq) noexcept;
2525
2526   inline CMK_REFNUM_TYPE getSeqNo(int destRank, MPI_Comm destcomm, int tag) noexcept;
2527   AmpiMsg *makeBcastMsg(const void *buf,int count,MPI_Datatype type,int root,MPI_Comm destcomm) noexcept;
2528   AmpiMsg *makeAmpiMsg(int destRank,int t,int sRank,const void *buf,int count,
2529                        MPI_Datatype type,MPI_Comm destcomm, int ssendReq=0) noexcept;
2530
2531   MPI_Request send(int t, int s, const void* buf, int count, MPI_Datatype type, int rank,
2532                    MPI_Comm destcomm, int ssendReq=0, AmpiSendType sendType=BLOCKING_SEND) noexcept;
2533   static void sendraw(int t, int s, void* buf, int len, CkArrayID aid, int idx) noexcept;
2534   inline MPI_Request sendLocalMsg(int t, int sRank, const void* buf, int size, MPI_Datatype type, int destRank,
2535                                   MPI_Comm destcomm, ampi* destPtr, int ssendReq, AmpiSendType sendType) noexcept;
2536   inline MPI_Request sendRdmaMsg(int t, int sRank, const void* buf, int size, MPI_Datatype type, int destIdx,
2537                                  int destRank, MPI_Comm destcomm, CProxy_ampi arrProxy, int ssendReq) noexcept;
2538   inline bool destLikelyWithinProcess(CProxy_ampi arrProxy, int destIdx) const noexcept {
2539     CkArray* localBranch = arrProxy.ckLocalBranch();
2540     int destPe = localBranch->lastKnown(CkArrayIndex1D(destIdx));
2541     return (CkNodeOf(destPe) == CkMyNode());
2542   }
2543   MPI_Request delesend(int t, int s, const void* buf, int count, MPI_Datatype type, int rank,
2544                        MPI_Comm destcomm, CProxy_ampi arrproxy, int ssend, AmpiSendType sendType) noexcept;
2545   inline void processAmpiMsg(AmpiMsg *msg, void* buf, MPI_Datatype type, int count) noexcept;
2546   inline void processRdmaMsg(const void *sbuf, int slength, int ssendReq, int srank, void* rbuf,
2547                              int rcount, MPI_Datatype rtype, MPI_Comm comm) noexcept;
2548   inline void processRednMsg(CkReductionMsg *msg, void* buf, MPI_Datatype type, int count) noexcept;
2549   inline void processNoncommutativeRednMsg(CkReductionMsg *msg, void* buf, MPI_Datatype type, int count,
2550                                            MPI_User_function* func) noexcept;
2551   inline void processGatherMsg(CkReductionMsg *msg, void* buf, MPI_Datatype type, int recvCount) noexcept;
2552   inline void processGathervMsg(CkReductionMsg *msg, void* buf, MPI_Datatype type,
2553                                int* recvCounts, int* displs) noexcept;
2554   inline AmpiMsg * getMessage(int t, int s, MPI_Comm comm, int *sts) const noexcept;
2555   int recv(int t,int s,void* buf,int count,MPI_Datatype type,MPI_Comm comm,MPI_Status *sts=NULL) noexcept;
2556   void irecv(void *buf, int count, MPI_Datatype type, int src,
2557              int tag, MPI_Comm comm, MPI_Request *request) noexcept;
2558   void mrecv(int tag, int src, void* buf, int count, MPI_Datatype datatype, MPI_Comm comm,
2559              MPI_Status* status, MPI_Message* message) noexcept;
2560   void imrecv(void* buf, int count, MPI_Datatype datatype, int src, int tag, MPI_Comm comm,
2561               MPI_Request* request, MPI_Message* message) noexcept;
2562   void irecvBcast(void *buf, int count, MPI_Datatype type, int src,
2563                   MPI_Comm comm, MPI_Request *request) noexcept;
2564   void sendrecv(const void *sbuf, int scount, MPI_Datatype stype, int dest, int stag,
2565                 void *rbuf, int rcount, MPI_Datatype rtype, int src, int rtag,
2566                 MPI_Comm comm, MPI_Status *sts) noexcept;
2567   void sendrecv_replace(void* buf, int count, MPI_Datatype datatype,
2568                         int dest, int sendtag, int source, int recvtag,
2569                         MPI_Comm comm, MPI_Status *status) noexcept;
2570   void probe(int t,int s,MPI_Comm comm,MPI_Status *sts) noexcept;
2571   void mprobe(int t, int s, MPI_Comm comm, MPI_Status *sts, MPI_Message *message) noexcept;
2572   int iprobe(int t,int s,MPI_Comm comm,MPI_Status *sts) noexcept;
2573   int improbe(int t, int s, MPI_Comm comm, MPI_Status *sts, MPI_Message *message) noexcept;
2574   void barrier() noexcept;
2575   void ibarrier(MPI_Request *request) noexcept;
2576   void bcast(int root, void* buf, int count, MPI_Datatype type, MPI_Comm comm) noexcept;
2577   int intercomm_bcast(int root, void* buf, int count, MPI_Datatype type, MPI_Comm intercomm) noexcept;
2578   void ibcast(int root, void* buf, int count, MPI_Datatype type, MPI_Comm comm, MPI_Request* request) noexcept;
2579   int intercomm_ibcast(int root, void* buf, int count, MPI_Datatype type, MPI_Comm intercomm, MPI_Request *request) noexcept;
2580   static void bcastraw(void* buf, int len, CkArrayID aid) noexcept;
2581   void split(int color,int key,MPI_Comm *dest, int type) noexcept;
2582   void commCreate(const vector<int>& vec,MPI_Comm *newcomm) noexcept;
2583   MPI_Comm cartCreate0D() noexcept;
2584   MPI_Comm cartCreate(vector<int>& vec, int ndims, const int* dims) noexcept;
2585   void graphCreate(const vector<int>& vec, MPI_Comm *newcomm) noexcept;
2586   void distGraphCreate(const vector<int>& vec, MPI_Comm *newcomm) noexcept;
2587   void intercommCreate(const vector<int>& rvec, int root, MPI_Comm tcomm, MPI_Comm *ncomm) noexcept;
2588
2589   inline bool isInter() const noexcept { return myComm.isinter(); }
2590   void intercommMerge(int first, MPI_Comm *ncomm) noexcept;
2591
2592   inline int getWorldRank() const noexcept {return parent->thisIndex;}
2593   /// Return our rank in this communicator
2594   inline int getRank() const noexcept {return myRank;}
2595   inline int getSize() const noexcept {return myComm.getSize();}
2596   inline MPI_Comm getComm() const noexcept {return myComm.getComm();}
2597   inline void setCommName(const char *name) noexcept {myComm.setName(name);}
2598   inline void getCommName(char *name, int *len) const noexcept {myComm.getName(name,len);}
2599   inline vector<int> getIndices() const noexcept { return myComm.getIndices(); }
2600   inline vector<int> getRemoteIndices() const noexcept { return myComm.getRemoteIndices(); }
2601   inline const CProxy_ampi &getProxy() const noexcept {return thisProxy;}
2602   inline const CProxy_ampi &getRemoteProxy() const noexcept {return remoteProxy;}
2603   inline void setRemoteProxy(CProxy_ampi rproxy) noexcept { remoteProxy = rproxy; thread->resume(); }
2604   inline int getIndexForRank(int r) const noexcept {return myComm.getIndexForRank(r);}
2605   inline int getIndexForRemoteRank(int r) const noexcept {return myComm.getIndexForRemoteRank(r);}
2606   void findNeighbors(MPI_Comm comm, int rank, vector<int>& neighbors) const noexcept;
2607   inline const vector<int>& getNeighbors() const noexcept { return myComm.getTopologyforNeighbors()->getnbors(); }
2608   inline bool opIsCommutative(MPI_Op op) const noexcept { return parent->opIsCommutative(op); }
2609   inline MPI_User_function* op2User_function(MPI_Op op) const noexcept { return parent->op2User_function(op); }
2610   void topoDup(int topoType, int rank, MPI_Comm comm, MPI_Comm *newcomm) noexcept;
2611
2612   inline AmpiRequestList& getReqs() noexcept { return parent->ampiReqs; }
2613   CkDDT *getDDT() noexcept {return &parent->myDDT;}
2614   CthThread getThread() const noexcept { return thread->getThread(); }
2615
2616  public:
2617   MPI_Win createWinInstance(void *base, MPI_Aint size, int disp_unit, MPI_Info info) noexcept;
2618   int deleteWinInstance(MPI_Win win) noexcept;
2619   int winGetGroup(WinStruct *win, MPI_Group *group) const noexcept;
2620   int winPut(const void *orgaddr, int orgcnt, MPI_Datatype orgtype, int rank,
2621              MPI_Aint targdisp, int targcnt, MPI_Datatype targtype, WinStruct *win) noexcept;
2622   int winGet(void *orgaddr, int orgcnt, MPI_Datatype orgtype, int rank,
2623              MPI_Aint targdisp, int targcnt, MPI_Datatype targtype, WinStruct *win) noexcept;
2624   int winIget(MPI_Aint orgdisp, int orgcnt, MPI_Datatype orgtype, int rank,
2625               MPI_Aint targdisp, int targcnt, MPI_Datatype targtype, WinStruct *win,
2626               MPI_Request *req) noexcept;
2627   int winIgetWait(MPI_Request *request, MPI_Status *status) noexcept;
2628   int winIgetFree(MPI_Request *request, MPI_Status *status) noexcept;
2629   void winRemotePut(int orgtotalsize, char* orgaddr, int orgcnt, MPI_Datatype orgtype,
2630                     MPI_Aint targdisp, int targcnt, MPI_Datatype targtype, int winIndex) noexcept;
2631   char* winLocalGet(int orgcnt, MPI_Datatype orgtype, MPI_Aint targdisp, int targcnt,
2632                     MPI_Datatype targtype, int winIndex) noexcept;
2633   AmpiMsg* winRemoteGet(int orgcnt, MPI_Datatype orgtype, MPI_Aint targdisp,
2634                     int targcnt, MPI_Datatype targtype, int winIndex) noexcept;
2635   AmpiMsg* winRemoteIget(MPI_Aint orgdisp, int orgcnt, MPI_Datatype orgtype, MPI_Aint targdisp,
2636                          int targcnt, MPI_Datatype targtype, int winIndex) noexcept;
2637   int winLock(int lock_type, int rank, WinStruct *win) noexcept;
2638   int winUnlock(int rank, WinStruct *win) noexcept;
2639   void winRemoteLock(int lock_type, int winIndex, int requestRank) noexcept;
2640   void winRemoteUnlock(int winIndex, int requestRank) noexcept;
2641   int winAccumulate(const void *orgaddr, int orgcnt, MPI_Datatype orgtype, int rank,
2642                     MPI_Aint targdisp, int targcnt, MPI_Datatype targtype,
2643                     MPI_Op op, WinStruct *win) noexcept;
2644   void winRemoteAccumulate(int orgtotalsize, char* orgaddr, int orgcnt, MPI_Datatype orgtype,
2645                            MPI_Aint targdisp, int targcnt, MPI_Datatype targtype,
2646                            MPI_Op op, int winIndex) noexcept;
2647   int winGetAccumulate(const void *orgaddr, int orgcnt, MPI_Datatype orgtype, void *resaddr,
2648                        int rescnt, MPI_Datatype restype, int rank, MPI_Aint targdisp,
2649                        int targcnt, MPI_Datatype targtype, MPI_Op op, WinStruct *win) noexcept;
2650   void winLocalGetAccumulate(int orgtotalsize, char* sorgaddr, int orgcnt, MPI_Datatype orgtype,
2651                              MPI_Aint targdisp, int targcnt, MPI_Datatype targtype, MPI_Op op,
2652                              char *resaddr, int winIndex) noexcept;
2653   AmpiMsg* winRemoteGetAccumulate(int orgtotalsize, char* sorgaddr, int orgcnt, MPI_Datatype orgtype,
2654                                   MPI_Aint targdisp, int targcnt, MPI_Datatype targtype, MPI_Op op,
2655                                   int winIndex) noexcept;
2656   int winCompareAndSwap(const void *orgaddr, const void *compaddr, void *resaddr, MPI_Datatype type,
2657                         int rank, MPI_Aint targdisp, WinStruct *win) noexcept;
2658   char* winLocalCompareAndSwap(int size, char* sorgaddr, char* compaddr, MPI_Datatype type,
2659                                MPI_Aint targdisp, int winIndex) noexcept;
2660   AmpiMsg* winRemoteCompareAndSwap(int size, char *sorgaddr, char *compaddr, MPI_Datatype type,
2661                                    MPI_Aint targdisp, int winIndex) noexcept;
2662   void winSetName(WinStruct *win, const char *name) noexcept;
2663   void winGetName(WinStruct *win, char *name, int *length) const noexcept;
2664   win_obj* getWinObjInstance(WinStruct *win) const noexcept;
2665   int getNewSemaId() noexcept;
2666
2667   int intercomm_scatter(int root, const void *sendbuf, int sendcount, MPI_Datatype sendtype,
2668                         void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm intercomm) noexcept;
2669   int intercomm_iscatter(int root, const void *sendbuf, int sendcount, MPI_Datatype sendtype,
2670                          void *recvbuf, int recvcount, MPI_Datatype recvtype,
2671                          MPI_Comm intercomm, MPI_Request *request) noexcept;
2672   int intercomm_scatterv(int root, const void* sendbuf, const int* sendcounts, const int* displs,
2673                          MPI_Datatype sendtype, void* recvbuf, int recvcount,
2674                          MPI_Datatype recvtype, MPI_Comm intercomm) noexcept;
2675   int intercomm_iscatterv(int root, const void* sendbuf, const int* sendcounts, const int* displs,
2676                           MPI_Datatype sendtype, void* recvbuf, int recvcount,
2677                           MPI_Datatype recvtype, MPI_Comm intercomm, MPI_Request* request) noexcept;
2678 };
2679
2680 ampiParent *getAmpiParent() noexcept;
2681 bool isAmpiThread() noexcept;
2682 ampi *getAmpiInstance(MPI_Comm comm) noexcept;
2683 void checkComm(MPI_Comm comm) noexcept;
2684 void checkRequest(MPI_Request req) noexcept;
2685 void handle_MPI_BOTTOM(void* &buf, MPI_Datatype type) noexcept;
2686 void handle_MPI_BOTTOM(void* &buf1, MPI_Datatype type1, void* &buf2, MPI_Datatype type2) noexcept;
2687
2688 #if AMPI_ERROR_CHECKING
2689 int ampiErrhandler(const char* func, int errcode) noexcept;
2690 #else
2691 #define ampiErrhandler(func, errcode) (errcode)
2692 #endif
2693
2694
2695 #if CMK_TRACE_ENABLED
2696
2697 // List of AMPI functions to trace:
2698 static const char *funclist[] = {"AMPI_Abort", "AMPI_Add_error_class", "AMPI_Add_error_code", "AMPI_Add_error_string",
2699 "AMPI_Address", "AMPI_Allgather", "AMPI_Allgatherv", "AMPI_Allreduce", "AMPI_Alltoall",
2700 "AMPI_Alltoallv", "AMPI_Alltoallw", "AMPI_Attr_delete", "AMPI_Attr_get",
2701 "AMPI_Attr_put", "AMPI_Barrier", "AMPI_Bcast", "AMPI_Bsend", "AMPI_Cancel",
2702 "AMPI_Cart_coords", "AMPI_Cart_create", "AMPI_Cart_get", "AMPI_Cart_map",
2703 "AMPI_Cart_rank", "AMPI_Cart_shift", "AMPI_Cart_sub", "AMPI_Cartdim_get",
2704 "AMPI_Comm_call_errhandler", "AMPI_Comm_compare", "AMPI_Comm_create", "AMPI_Comm_create_group",
2705 "AMPI_Comm_create_errhandler", "AMPI_Comm_create_keyval", "AMPI_Comm_delete_attr",
2706 "AMPI_Comm_dup", "AMPI_Comm_dup_with_info", "AMPI_Comm_free",
2707 "AMPI_Comm_free_errhandler", "AMPI_Comm_free_keyval", "AMPI_Comm_get_attr",
2708 "AMPI_Comm_get_errhandler", "AMPI_Comm_get_info", "AMPI_Comm_get_name",
2709 "AMPI_Comm_group", "AMPI_Comm_rank", "AMPI_Comm_remote_group", "AMPI_Comm_remote_size",
2710 "AMPI_Comm_set_attr", "AMPI_Comm_set_errhandler", "AMPI_Comm_set_info", "AMPI_Comm_set_name",
2711 "AMPI_Comm_size", "AMPI_Comm_split", "AMPI_Comm_split_type", "AMPI_Comm_test_inter",
2712 "AMPI_Dims_create", "AMPI_Dist_graph_create", "AMPI_Dist_graph_create_adjacent",
2713 "AMPI_Dist_graph_neighbors", "AMPI_Dist_graph_neighbors_count",
2714 "AMPI_Errhandler_create", "AMPI_Errhandler_free", "AMPI_Errhandler_get",
2715 "AMPI_Errhandler_set", "AMPI_Error_class", "AMPI_Error_string", "AMPI_Exscan", "AMPI_Finalize",
2716 "AMPI_Finalized", "AMPI_Gather", "AMPI_Gatherv", "AMPI_Get_address", "AMPI_Get_count",
2717 "AMPI_Get_elements", "AMPI_Get_library_version", "AMPI_Get_processor_name", "AMPI_Get_version",
2718 "AMPI_Graph_create", "AMPI_Graph_get", "AMPI_Graph_map", "AMPI_Graph_neighbors",
2719 "AMPI_Graph_neighbors_count", "AMPI_Graphdims_get", "AMPI_Group_compare", "AMPI_Group_difference",
2720 "AMPI_Group_excl", "AMPI_Group_free", "AMPI_Group_incl", "AMPI_Group_intersection",
2721 "AMPI_Group_range_excl", "AMPI_Group_range_incl", "AMPI_Group_rank", "AMPI_Group_size",
2722 "AMPI_Group_translate_ranks", "AMPI_Group_union", "AMPI_Iallgather", "AMPI_Iallgatherv",
2723 "AMPI_Iallreduce", "AMPI_Ialltoall", "AMPI_Ialltoallv", "AMPI_Ialltoallw", "AMPI_Ibarrier",
2724 "AMPI_Ibcast", "AMPI_Iexscan", "AMPI_Igather", "AMPI_Igatherv", "AMPI_Ineighbor_allgather",
2725 "AMPI_Ineighbor_allgatherv", "AMPI_Ineighbor_alltoall", "AMPI_Ineighbor_alltoallv",
2726 "AMPI_Ineighbor_alltoallw", "AMPI_Init", "AMPI_Init_thread", "AMPI_Initialized", "AMPI_Intercomm_create",
2727 "AMPI_Intercomm_merge", "AMPI_Iprobe", "AMPI_Irecv", "AMPI_Ireduce", "AMPI_Ireduce_scatter",
2728 "AMPI_Ireduce_scatter_block", "AMPI_Is_thread_main", "AMPI_Iscan", "AMPI_Iscatter", "AMPI_Iscatterv",
2729 "AMPI_Isend", "AMPI_Issend", "AMPI_Keyval_create", "AMPI_Keyval_free", "AMPI_Neighbor_allgather",
2730 "AMPI_Neighbor_allgatherv", "AMPI_Neighbor_alltoall", "AMPI_Neighbor_alltoallv", "AMPI_Neighbor_alltoallw",
2731 "AMPI_Op_commutative", "AMPI_Op_create", "AMPI_Op_free", "AMPI_Pack", "AMPI_Pack_size",
2732 "AMPI_Pcontrol", "AMPI_Probe", "AMPI_Query_thread", "AMPI_Recv", "AMPI_Recv_init", "AMPI_Reduce",
2733 "AMPI_Reduce_local", "AMPI_Reduce_scatter", "AMPI_Reduce_scatter_block", "AMPI_Request_free",
2734 "AMPI_Request_get_status", "AMPI_Rsend", "AMPI_Scan", "AMPI_Scatter", "AMPI_Scatterv", "AMPI_Send",
2735 "AMPI_Send_init",  "AMPI_Sendrecv", "AMPI_Sendrecv_replace", "AMPI_Ssend", "AMPI_Ssend_init",
2736 "AMPI_Start", "AMPI_Startall", "AMPI_Status_set_cancelled", "AMPI_Status_set_elements", "AMPI_Test",
2737 "AMPI_Test_cancelled", "AMPI_Testall", "AMPI_Testany", "AMPI_Testsome", "AMPI_Topo_test",
2738 "AMPI_Type_commit", "AMPI_Type_contiguous", "AMPI_Type_create_hindexed",
2739 "AMPI_Type_create_hindexed_block", "AMPI_Type_create_hvector", "AMPI_Type_create_indexed_block",
2740 "AMPI_Type_create_keyval", "AMPI_Type_create_resized", "AMPI_Type_create_struct",
2741 "AMPI_Type_delete_attr", "AMPI_Type_dup", "AMPI_Type_extent", "AMPI_Type_free",
2742 "AMPI_Type_free_keyval", "AMPI_Type_get_attr", "AMPI_Type_get_contents", "AMPI_Type_get_envelope",
2743 "AMPI_Type_get_extent", "AMPI_Type_get_name", "AMPI_Type_get_true_extent", "AMPI_Type_hindexed",
2744 "AMPI_Type_hvector", "AMPI_Type_indexed", "AMPI_Type_lb", "AMPI_Type_set_attr",
2745 "AMPI_Type_set_name", "AMPI_Type_size", "AMPI_Type_struct", "AMPI_Type_ub", "AMPI_Type_vector",
2746 "AMPI_Unpack", "AMPI_Wait", "AMPI_Waitall", "AMPI_Waitany", "AMPI_Waitsome", "AMPI_Wtick", "AMPI_Wtime",
2747 "AMPI_Accumulate", "AMPI_Compare_and_swap", "AMPI_Fetch_and_op", "AMPI_Get", "AMPI_Get_accumulate",
2748 "AMPI_Info_create", "AMPI_Info_delete", "AMPI_Info_dup", "AMPI_Info_free", "AMPI_Info_get",
2749 "AMPI_Info_get_nkeys", "AMPI_Info_get_nthkey", "AMPI_Info_get_valuelen",
2750 "AMPI_Info_set", "AMPI_Put", "AMPI_Raccumulate", "AMPI_Rget", "AMPI_Rget_accumulate",
2751 "AMPI_Rput", "AMPI_Win_complete", "AMPI_Win_create", "AMPI_Win_create_errhandler",
2752 "AMPI_Win_create_keyval", "AMPI_Win_delete_attr", "AMPI_Win_fence", "AMPI_Win_free",
2753 "AMPI_Win_free_keyval", "AMPI_Win_get_attr", "AMPI_Win_get_errhandler",
2754 "AMPI_Win_get_group", "AMPI_Win_get_info", "AMPI_Win_get_name", "AMPI_Win_lock",
2755 "AMPI_Win_post", "AMPI_Win_set_attr", "AMPI_Win_set_errhandler", "AMPI_Win_set_info",
2756 "AMPI_Win_set_name", "AMPI_Win_start", "AMPI_Win_test", "AMPI_Win_unlock",
2757 "AMPI_Win_wait", "AMPI_Exit" /*AMPI extensions:*/, "AMPI_Migrate",
2758 "AMPI_Load_start_measure", "AMPI_Load_stop_measure",
2759 "AMPI_Load_set_value", "AMPI_Migrate_to_pe", "AMPI_Set_migratable",
2760 "AMPI_Register_pup", "AMPI_Get_pup_data", "AMPI_Register_main",
2761 "AMPI_Register_about_to_migrate", "AMPI_Register_just_migrated",
2762 "AMPI_Iget", "AMPI_Iget_wait", "AMPI_Iget_free", "AMPI_Iget_data",
2763 "AMPI_Type_is_contiguous", "AMPI_Yield", "AMPI_Suspend",
2764 "AMPI_Resume", "AMPI_Print", "AMPI_Alltoall_medium",
2765 "AMPI_Alltoall_long", "AMPI_System"};
2766
2767 // not traced: AMPI_Trace_begin, AMPI_Trace_end
2768
2769 #endif // CMK_TRACE_ENABLED
2770
2771 //Use this to mark the start of AMPI interface routines that can only be called on AMPI threads:
2772 #if CMK_ERROR_CHECKING
2773 #define AMPI_API(routineName) \
2774   if (!isAmpiThread()) { CkAbort("AMPI> cannot call MPI routines from non-AMPI threads!"); } \
2775   TCHARM_API_TRACE(routineName, "ampi");
2776 #else
2777 #define AMPI_API(routineName) TCHARM_API_TRACE(routineName, "ampi")
2778 #endif
2779
2780 //Use this for MPI_Init and routines than can be called before AMPI threads have been initialized:
2781 #define AMPI_API_INIT(routineName) TCHARM_API_TRACE(routineName, "ampi")
2782
2783 #endif // _AMPIIMPL_H