43723808d50f5a2673615ce7a052bd455ca1061f
[charm.git] / src / arch / cuda / hybridAPI / cuda-hybrid-api.cu
1 /* 
2  * cuda-hybrid-api.cu
3  *
4  * by Lukasz Wesolowski
5  * 04.01.2008
6  *
7  * an interface for execution on the GPU
8  *
9  * description: 
10  * -user enqueues one or more work requests to the work
11  * request queue (wrQueue) to be executed on the GPU
12  * - a converse function (gpuProgressFn) executes periodically to
13  * offload work requests to the GPU one at a time
14  *
15  */
16
17 #include "wrqueue.h"
18 #include "cuda-hybrid-api.h"
19 #include "stdio.h"
20 #include <cutil.h>
21
22 /* A function in ck.C which casts the void * to a CkCallback object
23  *  and executes the callback 
24  */ 
25 extern void CUDACallbackManager(void * fn); 
26 extern int CmiMyPe();
27
28 /* initial size of the user-addressed portion of host/device buffer
29  * arrays; the system-addressed portion of host/device buffer arrays
30  * (used when there is no need to share buffers between work requests)
31  * will be equivalant in size.  
32  */ 
33 #define NUM_BUFFERS 256
34 #define MAX_PINNED_REQ 64  
35 #define MAX_DELAYED_FREE_REQS 64  
36
37 /* a flag which tells the system to record the time for invocation and
38  *  completion of GPU events: memory allocation, transfer and
39  *  kernel execution
40  */  
41 //#define GPU_PROFILE
42 //#define GPU_DEBUG
43 //#define GPU_TRACE
44 //#define _DEBUG
45
46 /* work request queue */
47 workRequestQueue *wrQueue = NULL; 
48
49 /* pending page-locked memory allocation requests */
50 unsigned int pinnedMemQueueIndex = 0; 
51 pinnedMemReq pinnedMemQueue[MAX_PINNED_REQ];
52
53 unsigned int currentDfr = 0;
54 void *delayedFreeReqs[MAX_DELAYED_FREE_REQS];
55
56 #ifdef GPU_MEMPOOL
57 #define GPU_MEMPOOL_NUM_SLOTS 15
58
59 CkVec<BufferPool> memPoolFreeBufs;
60 CkVec<int> memPoolBoundaries;
61 //int memPoolBoundaries[GPU_MEMPOOL_NUM_SLOTS];
62 #endif
63
64 /* The runtime system keeps track of all allocated buffers on the GPU.
65  * The following arrays contain pointers to host (CPU) data and the
66  * corresponding data on the device (GPU). 
67  */ 
68
69 /* host buffers  */ 
70 void **hostBuffers = NULL; 
71
72 /* device buffers */
73 void **devBuffers = NULL; 
74
75 /* used to assign bufferIDs automatically by the system if the user 
76    specifies an invalid bufferID */
77 unsigned int nextBuffer; 
78
79 unsigned int timerHandle; 
80
81 #ifdef GPU_PROFILE
82
83 /* event types */
84 #define DATA_SETUP          1            
85 #define KERNEL_EXECUTION    2
86 #define DATA_CLEANUP        3
87
88 typedef struct gpuEventTimer {
89   float startTime; 
90   float endTime; 
91   int eventType;
92   int ID; 
93 #ifdef GPU_TRACE
94   int stage; 
95   double cmistartTime; 
96   double cmiendTime; 
97 #endif
98 } gpuEventTimer; 
99
100 gpuEventTimer gpuEvents[QUEUE_SIZE_INIT * 3]; 
101 unsigned int timeIndex = 0; 
102 unsigned int runningKernelIndex = 0; 
103 unsigned int dataSetupIndex = 0; 
104 unsigned int dataCleanupIndex = 0; 
105
106 #if defined GPU_TRACE || defined GPU_INSTRUMENT_WRS
107 extern "C" double CmiWallTimer(); 
108 #endif
109
110 #ifdef GPU_TRACE
111 extern "C" int traceRegisterUserEvent(const char*x, int e);
112 extern "C" void traceUserBracketEvent(int e, double beginT, double endT);
113
114 #define GPU_MEM_SETUP 8800
115 #define GPU_KERNEL_EXEC 8801
116 #define GPU_MEM_CLEANUP 8802
117
118 #endif
119
120 #endif
121
122 #ifdef GPU_INSTRUMENT_WRS
123 CkVec<CkVec<CkVec<RequestTimeInfo> > > avgTimes;
124 bool initialized_instrument;
125 bool initializedInstrument();
126 #endif
127
128 /* There are separate CUDA streams for kernel execution, data transfer
129  *  into the device, and data transfer out. This allows prefetching of
130  *  data for a subsequent kernel while the previous kernel is
131  *  executing and transferring data out of the device. 
132  */
133 cudaStream_t kernel_stream; 
134 cudaStream_t data_in_stream;
135 cudaStream_t data_out_stream; 
136
137 /* pinnedMallocHost
138  *
139  * schedules a pinned memory allocation so that it does not impede
140  * concurrent asynchronous execution 
141  *
142  */
143 void pinnedMallocHost(pinnedMemReq *reqs) {
144
145   if ( (cudaStreamQuery(kernel_stream) == cudaSuccess) &&
146        (cudaStreamQuery(data_in_stream) == cudaSuccess) &&
147        (cudaStreamQuery(data_out_stream) == cudaSuccess) ) {    
148
149
150
151     for (int i=0; i<reqs->nBuffers; i++) {
152       CUDA_SAFE_CALL_NO_SYNC(cudaMallocHost((void **) reqs->hostPtrs[i], 
153                                             reqs->sizes[i])); 
154     }
155
156     free(reqs->hostPtrs);
157     free(reqs->sizes);
158
159     CUDACallbackManager(reqs->callbackFn);
160
161   }
162   else {
163     pinnedMemQueue[pinnedMemQueueIndex].hostPtrs = reqs->hostPtrs;
164     pinnedMemQueue[pinnedMemQueueIndex].sizes = reqs->sizes; 
165     pinnedMemQueue[pinnedMemQueueIndex].callbackFn = reqs->callbackFn;     
166     pinnedMemQueueIndex++;
167     if (pinnedMemQueueIndex == MAX_PINNED_REQ) {
168       printf("Error: pinned memory request buffer is overflowing\n"); 
169     }
170   }
171 }
172
173 void delayedFree(void *ptr){
174   if(currentDfr == MAX_DELAYED_FREE_REQS){
175     printf("Ran out of DFR queue space. Increase MAX_DELAYED_FREE_REQS\n");
176     exit(-1);
177   }
178   else{
179     delayedFreeReqs[currentDfr] = ptr;
180   }
181   currentDfr++;
182 }
183
184 void flushDelayedFrees(){
185   for(int i = 0; i < currentDfr; i++){
186     if(delayedFreeReqs[i] == NULL){
187       printf("recorded NULL ptr in delayedFree()");
188       exit(-1);
189     }
190     cudaFreeHost(delayedFreeReqs[i]);
191   }
192 }
193
194 /* flushPinnedMemQueue
195  *
196  * executes pending pinned memory allocation requests
197  *
198  */
199 void flushPinnedMemQueue() {
200
201   for (int i=0; i<pinnedMemQueueIndex; i++) {
202     pinnedMemReq *req = &pinnedMemQueue[i]; 
203     for (int j=0; j<req->nBuffers; j++) {
204       CUDA_SAFE_CALL_NO_SYNC(cudaMallocHost((void **) req->hostPtrs[j], 
205                                             req->sizes[j])); 
206     }
207     free(req->hostPtrs);
208     free(req->sizes);
209     CUDACallbackManager(pinnedMemQueue[i].callbackFn);    
210   }
211   pinnedMemQueueIndex = 0; 
212
213 }
214
215 /* allocateBuffers
216  *
217  * allocates a work request's data on the GPU
218  *
219  * used to allocate memory for work request data in advance in order
220  * to allow overlapping the work request's data transfer to the GPU
221  * with the execution of the previous kernel; the allocation needs to
222  * take place before the kernel starts executing in order to allow overlap
223  *
224  */
225
226 void allocateBuffers(workRequest *wr) {
227   dataInfo *bufferInfo = wr->bufferInfo; 
228
229   if (bufferInfo != NULL) {
230
231     for (int i=0; i<wr->nBuffers; i++) {
232       int index = bufferInfo[i].bufferID; 
233       int size = bufferInfo[i].size; 
234
235       // if index value is invalid, use an available ID  
236       if (index < 0 || index >= NUM_BUFFERS) {
237         int found = 0; 
238         for (int j=nextBuffer; j<NUM_BUFFERS*2; j++) {
239           if (devBuffers[j] == NULL) {
240             index = j;
241             found = 1; 
242             break;
243           }
244         }
245
246         /* if no index was found, try to search for a value at the
247          * beginning of the system addressed space 
248          */
249         
250         if (!found) {
251           for (int j=NUM_BUFFERS; j<nextBuffer; j++) {
252             if (devBuffers[j] == NULL) {        
253               index = j;
254               found = 1; 
255               break;
256             }
257           }
258         }
259
260         /* if no index was found, print an error */
261         if (!found) {
262           printf("Error: devBuffers is full \n");
263         }
264
265         nextBuffer = index+1; 
266         if (nextBuffer == NUM_BUFFERS * 2) {
267           nextBuffer = NUM_BUFFERS; 
268         }
269         
270         bufferInfo[i].bufferID = index; 
271
272       }      
273       
274       // allocate if the buffer for the corresponding index is NULL 
275       if (devBuffers[index] == NULL && size > 0) {
276 #ifdef GPU_PRINT_BUFFER_ALLOCATE
277         double mil = 1e3;
278         printf("*** ALLOCATE buffer 0x%x (%d) size %f kb\n", devBuffers[index], index, 1.0*size/mil);
279
280 #endif
281
282         CUDA_SAFE_CALL_NO_SYNC(cudaMalloc((void **) &devBuffers[index], size));
283 #ifdef GPU_DEBUG
284         printf("buffer %d allocated at time %.2f size: %d error string: %s\n", 
285                index, cutGetTimerValue(timerHandle), size, 
286                cudaGetErrorString( cudaGetLastError() ) );
287 #endif
288       }
289     }
290   }
291 }
292
293
294 /* setupData
295  *  copy data to the GPU before kernel execution 
296  */
297 void setupData(workRequest *wr) {
298   dataInfo *bufferInfo = wr->bufferInfo; 
299
300   if (bufferInfo != NULL) {
301     for (int i=0; i<wr->nBuffers; i++) {
302       int index = bufferInfo[i].bufferID; 
303       int size = bufferInfo[i].size; 
304       hostBuffers[index] = bufferInfo[i].hostBuffer; 
305       
306       /* allocate if the buffer for the corresponding index is NULL */
307       /*
308       if (devBuffers[index] == NULL) {
309         CUDA_SAFE_CALL_NO_SYNC(cudaMalloc((void **) &devBuffers[index], size));
310 #ifdef GPU_DEBUG
311         printf("buffer %d allocated %.2f\n", index,
312                cutGetTimerValue(timerHandle)); 
313 #endif
314       }
315       */
316       
317       if (bufferInfo[i].transferToDevice && size > 0) {
318         CUDA_SAFE_CALL_NO_SYNC(cudaMemcpyAsync(devBuffers[index], 
319           hostBuffers[index], size, cudaMemcpyHostToDevice, data_in_stream));
320 #ifdef GPU_DEBUG
321         printf("transferToDevice bufId: %d at time %.2f size: %d " 
322                "error string: %s\n", index, cutGetTimerValue(timerHandle), 
323                size, cudaGetErrorString( cudaGetLastError() )); 
324 #endif  
325         /*
326         CUDA_SAFE_CALL_NO_SYNC(cudaMemcpy(devBuffers[index], 
327           hostBuffers[index], size, cudaMemcpyHostToDevice));
328         */
329
330       }
331     }
332   }
333
334
335 /* copybackData
336  *  transfer data from the GPU to the CPU after a work request is done 
337  */ 
338 void copybackData(workRequest *wr) {
339   dataInfo *bufferInfo = wr->bufferInfo; 
340
341   if (bufferInfo != NULL) {
342     int nBuffers = wr->nBuffers; 
343     
344     for (int i=0; i<nBuffers; i++) {
345       int index = bufferInfo[i].bufferID; 
346       int size = bufferInfo[i].size; 
347       
348       if (bufferInfo[i].transferFromDevice && size > 0) {
349 #ifdef GPU_DEBUG
350         printf("transferFromDevice: %d at time %.2f size: %d "
351                "error string: %s\n", index, cutGetTimerValue(timerHandle), 
352                size, cudaGetErrorString( cudaGetLastError() )); 
353 #endif
354         
355         CUDA_SAFE_CALL_NO_SYNC(cudaMemcpyAsync(hostBuffers[index], 
356           devBuffers[index], size, cudaMemcpyDeviceToHost,
357           data_out_stream));
358         
359         /*
360         CUDA_SAFE_CALL_NO_SYNC(cudaMemcpy(hostBuffers[index], 
361           devBuffers[index], size, cudaMemcpyDeviceToHost));
362         */
363       }
364     }     
365   }
366 }
367
368 /* frees GPU memory for buffers specified by the user; also frees the
369  *  work request's bufferInfo array
370  */
371 void freeMemory(workRequest *wr) {
372   dataInfo *bufferInfo = wr->bufferInfo;   
373   int nBuffers = wr->nBuffers; 
374   if (bufferInfo != NULL) {
375     for (int i=0; i<nBuffers; i++) {    
376       int index = bufferInfo[i].bufferID; 
377       if (bufferInfo[i].freeBuffer) {
378 #ifdef GPU_PRINT_BUFFER_ALLOCATE
379         printf("*** FREE buffer 0x%x (%d)\n", devBuffers[index], index);
380 #endif
381
382 #ifdef GPU_DEBUG
383         printf("buffer %d freed at time %.2f error string: %s\n", 
384                index, cutGetTimerValue(timerHandle),  
385                cudaGetErrorString( cudaGetLastError() ));
386 #endif 
387         CUDA_SAFE_CALL_NO_SYNC(cudaFree(devBuffers[index])); 
388         devBuffers[index] = NULL; 
389       }
390     }
391     free(bufferInfo); 
392   }
393 }
394
395 /* 
396  * a switch statement defined by the user to allow the library to execute
397  * the correct kernel 
398  */ 
399 void kernelSelect(workRequest *wr);
400
401 /* initHybridAPI
402  *   initializes the work request queue, host/device buffer pointer
403  *   arrays, and CUDA streams
404  */
405 void initHybridAPI(int myPe) {
406
407   int deviceCount;
408   cudaGetDeviceCount(&deviceCount);
409
410   cudaSetDevice(myPe % deviceCount); 
411
412   initWRqueue(&wrQueue);
413
414   /* allocate host/device buffers array (both user and
415      system-addressed) */
416   hostBuffers = (void **) malloc(NUM_BUFFERS * 2 * sizeof(void *)); 
417   devBuffers = (void **) malloc(NUM_BUFFERS * 2 * sizeof(void *)); 
418
419   /* initialize device array to NULL */ 
420   for (int i=0; i<NUM_BUFFERS*2; i++) {
421     devBuffers[i] = NULL; 
422   }
423   
424   CUDA_SAFE_CALL_NO_SYNC(cudaStreamCreate(&kernel_stream)); 
425   CUDA_SAFE_CALL_NO_SYNC(cudaStreamCreate(&data_in_stream)); 
426   CUDA_SAFE_CALL_NO_SYNC(cudaStreamCreate(&data_out_stream)); 
427
428 #ifdef GPU_PROFILE
429   CUT_SAFE_CALL(cutCreateTimer(&timerHandle));
430   CUT_SAFE_CALL(cutStartTimer(timerHandle));
431 #endif
432
433   nextBuffer = NUM_BUFFERS;  
434
435 #ifdef GPU_TRACE
436   traceRegisterUserEvent("GPU Memory Setup", GPU_MEM_SETUP);
437   traceRegisterUserEvent("GPU Kernel Execution", GPU_KERNEL_EXEC);
438   traceRegisterUserEvent("GPU Memory Cleanup", GPU_MEM_CLEANUP);
439 #endif
440
441 #ifdef GPU_MEMPOOL
442   int nslots = GPU_MEMPOOL_NUM_SLOTS;
443   int *sizes;
444   sizes = (int *)malloc(sizeof(int)*nslots); 
445
446   memPoolBoundaries.reserve(GPU_MEMPOOL_NUM_SLOTS);
447   memPoolBoundaries.length() = GPU_MEMPOOL_NUM_SLOTS;
448
449   int bufSize = GPU_MEMPOOL_MIN_BUFFER_SIZE;
450   for(int i = 0; i < GPU_MEMPOOL_NUM_SLOTS; i++){
451     memPoolBoundaries[i] = bufSize;
452     bufSize = bufSize << 1;
453   }
454
455   //1K
456   sizes[0] = 512; 
457   //2K
458   sizes[1] = 512;
459   //4K
460   sizes[2] = 64;
461   //8K
462   sizes[3] = 64;
463   //16K
464   sizes[4] = 32;
465   //32K
466   sizes[5] = 32;
467   //64K
468   sizes[6] = 32;
469   //128K
470   sizes[7] = 32;
471   //256K
472   sizes[8] = 32;
473   //512K
474   sizes[9] = 32;
475   //1M
476   sizes[10] = 170;
477   //2M
478   sizes[11] = 16;
479   //4M
480   sizes[12] = 4;
481   //8M
482   sizes[13] = 2;
483   //16M
484   sizes[14] = 2; 
485
486   createPool(sizes, nslots, memPoolFreeBufs);
487   printf("[%d] done creating buffer pool\n", CmiMyPe());
488
489 #endif
490
491 #ifdef GPU_INSTRUMENT_WRS
492   initialized_instrument = false;
493 #endif
494 }
495
496 /* gpuProgressFn
497  *  called periodically to monitor work request progress, and perform
498  *  the prefetch of data for a subsequent work request
499  */
500 void gpuProgressFn() {
501   if (wrQueue == NULL) {
502     printf("Error: work request queue not initialized\n"); 
503     return; 
504   }
505   if (isEmpty(wrQueue)) {
506     flushPinnedMemQueue();    
507     flushDelayedFrees();
508     return;
509   } 
510   int returnVal; 
511   workRequest *head = firstElement(wrQueue); 
512   workRequest *second = secondElement(wrQueue);
513   workRequest *third = thirdElement(wrQueue); 
514
515   if (head->state == QUEUED) {
516 #ifdef GPU_PROFILE
517     gpuEvents[timeIndex].startTime = cutGetTimerValue(timerHandle); 
518     gpuEvents[timeIndex].eventType = DATA_SETUP; 
519     gpuEvents[timeIndex].ID = head->id; 
520     dataSetupIndex = timeIndex; 
521 #ifdef GPU_TRACE
522     gpuEvents[timeIndex].stage = GPU_MEM_SETUP; 
523     gpuEvents[timeIndex].cmistartTime = CmiWallTimer();
524 #endif
525     timeIndex++; 
526 #endif
527
528 #ifdef GPU_INSTRUMENT_WRS
529     head->startTime = CmiWallTimer(); 
530 #endif
531
532     allocateBuffers(head); 
533     setupData(head); 
534     head->state = TRANSFERRING_IN; 
535   }  
536   if (head->state == TRANSFERRING_IN) {
537     if ((returnVal = cudaStreamQuery(data_in_stream)) == cudaSuccess) {
538 #ifdef GPU_PROFILE
539       gpuEvents[dataSetupIndex].endTime = cutGetTimerValue(timerHandle);
540 #ifdef GPU_TRACE
541       gpuEvents[dataSetupIndex].cmiendTime = CmiWallTimer();
542       traceUserBracketEvent(gpuEvents[dataSetupIndex].stage, 
543                             gpuEvents[dataSetupIndex].cmistartTime, 
544                             gpuEvents[dataSetupIndex].cmiendTime); 
545 #endif 
546 #endif
547
548 #ifdef GPU_INSTRUMENT_WRS
549       {
550         if(initializedInstrument()){
551           double tt = CmiWallTimer()-(head->startTime);
552           int index = head->chareIndex;
553           char type = head->compType;
554           char phase = head->compPhase;
555
556           CkVec<RequestTimeInfo> &vec = avgTimes[index][type];
557           if(vec.length() <= phase){
558             vec.growAtLeast(phase);
559           }
560           vec[phase].transferTime += tt;
561         }
562       }
563 #endif
564
565       if (second != NULL /*&& (second->state == QUEUED)*/) {
566         allocateBuffers(second); 
567       }
568 #ifdef GPU_PROFILE
569       gpuEvents[timeIndex].startTime = cutGetTimerValue(timerHandle); 
570       gpuEvents[timeIndex].eventType = KERNEL_EXECUTION; 
571       gpuEvents[timeIndex].ID = head->id; 
572       runningKernelIndex = timeIndex; 
573 #ifdef GPU_TRACE
574       gpuEvents[timeIndex].stage = GPU_KERNEL_EXEC; 
575       gpuEvents[timeIndex].cmistartTime = CmiWallTimer();
576 #endif
577       timeIndex++; 
578 #endif
579 #ifdef GPU_INSTRUMENT_WRS
580       head->startTime = CmiWallTimer(); 
581 #endif
582
583       //flushPinnedMemQueue();
584       flushDelayedFrees();
585       kernelSelect(head); 
586
587       head->state = EXECUTING; 
588       if (second != NULL) {
589 #ifdef GPU_PROFILE
590         gpuEvents[timeIndex].startTime = cutGetTimerValue(timerHandle); 
591         gpuEvents[timeIndex].eventType = DATA_SETUP; 
592         gpuEvents[timeIndex].ID = second->id; 
593         dataSetupIndex = timeIndex; 
594 #ifdef GPU_TRACE
595         gpuEvents[timeIndex].stage = GPU_MEM_SETUP; 
596         gpuEvents[timeIndex].cmistartTime = CmiWallTimer();
597 #endif
598         timeIndex++; 
599 #endif
600
601 #ifdef GPU_INSTRUMENT_WRS
602         second->startTime = CmiWallTimer();
603 #endif
604         setupData(second); 
605         second->state = TRANSFERRING_IN;
606       }
607     }
608       /*
609 #ifdef GPU_DEBUG
610       printf("Querying memory stream returned: %d %.2f\n", returnVal, 
611              cutGetTimerValue(timerHandle));
612 #endif  
613       */
614   }
615   if (head->state == EXECUTING) {
616     if ((returnVal = cudaStreamQuery(kernel_stream)) == cudaSuccess) {
617 #ifdef GPU_PROFILE
618       gpuEvents[runningKernelIndex].endTime = cutGetTimerValue(timerHandle); 
619 #ifdef GPU_TRACE
620       gpuEvents[runningKernelIndex].cmiendTime = CmiWallTimer();
621       traceUserBracketEvent(gpuEvents[runningKernelIndex].stage, 
622                             gpuEvents[runningKernelIndex].cmistartTime, 
623                             gpuEvents[runningKernelIndex].cmiendTime); 
624 #endif
625 #endif
626 #ifdef GPU_INSTRUMENT_WRS
627       {
628         if(initializedInstrument()){
629           double tt = CmiWallTimer()-(head->startTime);
630           int index = head->chareIndex;
631           char type = head->compType;
632           char phase = head->compPhase;
633
634           CkVec<RequestTimeInfo> &vec = avgTimes[index][type];
635           if(vec.length() <= phase){
636             vec.growAtLeast(phase);
637           }
638           vec[phase].kernelTime += tt;
639         }
640       }
641 #endif
642
643       if (second != NULL && second->state == QUEUED) {
644 #ifdef GPU_PROFILE
645         gpuEvents[timeIndex].startTime = cutGetTimerValue(timerHandle); 
646         gpuEvents[timeIndex].eventType = DATA_SETUP; 
647         gpuEvents[timeIndex].ID = second->id; 
648         dataSetupIndex = timeIndex; 
649 #ifdef GPU_TRACE
650         gpuEvents[timeIndex].stage = GPU_MEM_SETUP; 
651         gpuEvents[timeIndex].cmistartTime = CmiWallTimer();
652 #endif
653         timeIndex++; 
654 #endif
655
656 #ifdef GPU_INSTRUMENT_MS
657         second->startTime = CmiWallTimer();
658 #endif
659         
660         allocateBuffers(second); 
661         setupData(second); 
662         second->state = TRANSFERRING_IN;        
663       } 
664       if (second != NULL && second->state == TRANSFERRING_IN) {
665         if (cudaStreamQuery(data_in_stream) == cudaSuccess) {
666 #ifdef GPU_PROFILE
667           gpuEvents[dataSetupIndex].endTime = cutGetTimerValue(timerHandle); 
668 #ifdef GPU_TRACE
669           gpuEvents[dataSetupIndex].cmiendTime = CmiWallTimer();
670           traceUserBracketEvent(gpuEvents[dataSetupIndex].stage, 
671                                 gpuEvents[dataSetupIndex].cmistartTime, 
672                                 gpuEvents[dataSetupIndex].cmiendTime); 
673 #endif
674 #endif
675 #ifdef GPU_INSTRUMENT_WRS
676           {
677             if(initializedInstrument()){
678               double tt = CmiWallTimer()-(head->startTime);
679               int index = second->chareIndex;
680               char type = second->compType;
681               char phase = second->compPhase;
682
683               CkVec<RequestTimeInfo> &vec = avgTimes[index][type];
684               if(vec.length() <= phase){
685                 vec.growAtLeast(phase);
686               }
687               vec[phase].transferTime += tt;
688             }
689           }
690 #endif
691
692           if (third != NULL /*&& (third->state == QUEUED)*/) {
693             allocateBuffers(third); 
694           }
695 #ifdef GPU_PROFILE
696           gpuEvents[timeIndex].startTime = cutGetTimerValue(timerHandle); 
697           gpuEvents[timeIndex].eventType = KERNEL_EXECUTION; 
698           gpuEvents[timeIndex].ID = second->id; 
699           runningKernelIndex = timeIndex; 
700 #ifdef GPU_TRACE
701           gpuEvents[timeIndex].stage = GPU_KERNEL_EXEC; 
702           gpuEvents[timeIndex].cmistartTime = CmiWallTimer();
703 #endif
704           timeIndex++; 
705 #endif
706 #ifdef GPU_INSTRUMENT_WRS
707           second->startTime = CmiWallTimer();
708 #endif
709           //        flushPinnedMemQueue();          
710           flushDelayedFrees();
711           kernelSelect(second); 
712           second->state = EXECUTING; 
713           if (third != NULL) {
714 #ifdef GPU_PROFILE
715             gpuEvents[timeIndex].startTime = cutGetTimerValue(timerHandle); 
716             gpuEvents[timeIndex].eventType = DATA_SETUP; 
717             gpuEvents[timeIndex].ID = third->id; 
718             dataSetupIndex = timeIndex; 
719 #ifdef GPU_TRACE
720             gpuEvents[timeIndex].stage = GPU_MEM_SETUP; 
721             gpuEvents[timeIndex].cmistartTime = CmiWallTimer();
722 #endif
723             timeIndex++; 
724 #endif
725
726 #ifdef GPU_INSTRUMENT_WRS
727             third->startTime = CmiWallTimer();
728 #endif
729             setupData(third); 
730             third->state = TRANSFERRING_IN;     
731           }
732         }
733       }
734 #ifdef GPU_PROFILE
735       gpuEvents[timeIndex].startTime = cutGetTimerValue(timerHandle); 
736       gpuEvents[timeIndex].eventType = DATA_CLEANUP; 
737       gpuEvents[timeIndex].ID = head->id; 
738       dataCleanupIndex = timeIndex;     
739 #ifdef GPU_TRACE
740       gpuEvents[timeIndex].stage = GPU_MEM_CLEANUP; 
741       gpuEvents[timeIndex].cmistartTime = CmiWallTimer();
742 #endif
743       timeIndex++; 
744 #endif
745 #ifdef GPU_INSTRUMENT_WRS
746       head->startTime = CmiWallTimer(); 
747 #endif
748       copybackData(head);
749       head->state = TRANSFERRING_OUT;
750     }
751       /*
752 #ifdef GPU_DEBUG
753       printf("Querying kernel completion returned: %d %.2f\n", returnVal,
754              cutGetTimerValue(timerHandle));
755 #endif  
756       */
757   }
758   if (head->state == TRANSFERRING_OUT) {
759     if (cudaStreamQuery(data_out_stream) == cudaSuccess && cudaStreamQuery(kernel_stream) == cudaSuccess){
760       freeMemory(head); 
761 #ifdef GPU_PROFILE
762       gpuEvents[dataCleanupIndex].endTime = cutGetTimerValue(timerHandle);
763 #ifdef GPU_TRACE
764       gpuEvents[dataCleanupIndex].cmiendTime = CmiWallTimer();
765       traceUserBracketEvent(gpuEvents[dataCleanupIndex].stage, 
766                             gpuEvents[dataCleanupIndex].cmistartTime, 
767                             gpuEvents[dataCleanupIndex].cmiendTime); 
768 #endif
769 #endif
770 #ifdef GPU_INSTRUMENT_WRS
771       {
772         if(initializedInstrument()){
773           double tt = CmiWallTimer()-(head->startTime);
774           int index = head->chareIndex;
775           char type = head->compType;
776           char phase = head->compPhase;
777
778           CkVec<RequestTimeInfo> &vec = avgTimes[index][type];
779           if(vec.length() <= phase){
780             vec.growAtLeast(phase);
781           }
782           vec[phase].cleanupTime += tt;
783           vec[phase].n++;
784         }
785       }
786 #endif
787
788       dequeue(wrQueue);
789       CUDACallbackManager(head->callbackFn);
790     }
791   }
792 }
793
794 #ifdef GPU_MEMPOOL
795 void releasePool(CkVec<BufferPool> &pools);
796 #endif
797 /* exitHybridAPI
798  *  cleans up and deletes memory allocated for the queue and the CUDA streams
799  */
800 void exitHybridAPI() {
801   deleteWRqueue(wrQueue); 
802   CUDA_SAFE_CALL_NO_SYNC(cudaStreamDestroy(kernel_stream)); 
803   CUDA_SAFE_CALL_NO_SYNC(cudaStreamDestroy(data_in_stream)); 
804   CUDA_SAFE_CALL_NO_SYNC(cudaStreamDestroy(data_out_stream)); 
805
806 #ifdef GPU_PROFILE
807   for (int i=0; i<timeIndex; i++) {
808     switch (gpuEvents[i].eventType) {
809     case DATA_SETUP:
810       printf("Kernel %d data setup", gpuEvents[i].ID); 
811       break;
812     case DATA_CLEANUP:
813       printf("Kernel %d data cleanup", gpuEvents[i].ID); 
814       break; 
815     case KERNEL_EXECUTION:
816       printf("Kernel %d execution", gpuEvents[i].ID); 
817       break;
818     default:
819       printf("Error, invalid timer identifier\n"); 
820     }
821     printf(" %.2f:%.2f\n", gpuEvents[i].startTime-gpuEvents[0].startTime, gpuEvents[i].endTime-gpuEvents[0].startTime); 
822   }
823
824   CUT_SAFE_CALL(cutStopTimer(timerHandle));
825   CUT_SAFE_CALL(cutDeleteTimer(timerHandle));  
826
827 #endif
828
829 #ifdef GPU_MEMPOOL
830   releasePool(memPoolFreeBufs);
831 #endif
832
833 }
834
835 #ifdef GPU_MEMPOOL
836 void releasePool(CkVec<BufferPool> &pools){
837   for(int i = 0; i < pools.length(); i++){
838     CUDA_SAFE_CALL_NO_SYNC(cudaFreeHost((void *)pools[i].head));
839   }
840   pools.free();
841 }
842
843 // Create a pool with nslots slots.
844 // There are nbuffers[i] buffers for each buffer size corresponding to entry i
845 // FIXME - list the alignment/fragmentation issues with either of two allocation schemes:
846 // if a single, large buffer is allocated for each subpool
847 // if multiple smaller buffers are allocated for each subpool
848 void createPool(int *nbuffers, int nslots, CkVec<BufferPool> &pools){
849   //pools  = (BufferPool *)malloc(nslots*sizeof(BufferPool));
850   pools.reserve(nslots);
851   pools.length() = nslots;
852
853   for(int i = 0; i < nslots; i++){
854     int bufSize = memPoolBoundaries[i];
855     int numBuffers = nbuffers[i];
856     pools[i].size = bufSize;
857     
858     CUDA_SAFE_CALL_NO_SYNC(cudaMallocHost((void **)(&pools[i].head), 
859                                           (sizeof(Header)+bufSize)*numBuffers));
860     if(pools[i].head == NULL){
861       abort();
862     }
863
864     Header *hd = pools[i].head;
865     Header *previous = NULL;
866     char *memory;
867
868     for(int j = 0; j < numBuffers; j++){
869       hd->slot = i;
870       hd->next = previous;
871       previous = hd;
872       hd++; // move ptr past header
873       memory = (char *)hd;
874       memory += bufSize;
875       hd = (Header *)memory;
876     }
877
878     pools[i].head = previous;
879 #ifdef GPU_MEMPOOL_DEBUG
880     pools[i].num = numBuffers;
881 #endif
882   }
883 }
884
885 int findPool(int size){
886   int boundaryArrayLen = memPoolBoundaries.length();
887   if(size <= memPoolBoundaries[0]){
888     return (0);
889   }
890   else if(size > memPoolBoundaries[boundaryArrayLen-1]){
891     // create new slot
892     memPoolBoundaries.push_back(size);
893
894     BufferPool newpool;
895     CUDA_SAFE_CALL_NO_SYNC(cudaMallocHost((void **)&newpool.head, size+sizeof(Header)));
896     newpool.size = size;
897 #ifdef GPU_MEMPOOL_DEBUG
898     newpool.num = 1;
899 #endif
900     memPoolFreeBufs.push_back(newpool);
901
902     Header *hd = newpool.head;
903     hd->next = NULL;
904     hd->slot = boundaryArrayLen;
905
906     return boundaryArrayLen;
907   }
908   for(int i = 0; i < GPU_MEMPOOL_NUM_SLOTS-1; i++){
909     if(memPoolBoundaries[i] < size && size <= memPoolBoundaries[i+1]){
910       return (i+1);
911     }
912   }
913   return -1;
914 }
915
916 void *getBufferFromPool(int pool, int size){
917   Header *ret;
918   if(pool < 0 || pool >= memPoolFreeBufs.length() || memPoolFreeBufs[pool].head == NULL){
919 #ifdef GPU_MEMPOOL_DEBUG
920     printf("(%d) pool %d size: %d, num: %d\n", CmiMyPe(), pool, size, memPoolFreeBufs[pool].num);
921 #endif
922     abort();
923   }
924   else{
925     ret = memPoolFreeBufs[pool].head;
926     memPoolFreeBufs[pool].head = ret->next;
927 #ifdef GPU_MEMPOOL_DEBUG
928     ret->size = size;
929     memPoolFreeBufs[pool].num--;
930 #endif
931     return (void *)(ret+1);
932   }
933   return NULL;
934 }
935
936 void returnBufferToPool(int pool, Header *hd){
937   hd->next = memPoolFreeBufs[pool].head;
938   memPoolFreeBufs[pool].head = hd;
939 #ifdef GPU_MEMPOOL_DEBUG
940   memPoolFreeBufs[pool].num++;
941 #endif
942 }
943
944 void *hapi_poolMalloc(int size){
945   int pool = findPool(size);
946   void *buf = getBufferFromPool(pool, size);
947 #ifdef GPU_MEMPOOL_DEBUG
948   printf("(%d) hapi_malloc size %d pool %d left %d\n", CmiMyPe(), size, pool, memPoolFreeBufs[pool].num);
949 #endif
950   return buf;
951 }
952
953 void hapi_poolFree(void *ptr){
954   Header *hd = ((Header *)ptr)-1;
955   int pool = hd->slot;
956   returnBufferToPool(pool, hd);
957 #ifdef GPU_MEMPOOL_DEBUG
958   int size = hd->size;
959   printf("(%d) hapi_free size %d pool %d left %d\n", CmiMyPe(), size, pool, memPoolFreeBufs[pool].num);
960 #endif
961 }
962
963
964 #endif
965
966 #ifdef GPU_INSTRUMENT_WRS
967 void hapi_initInstrument(int numChares, char types){
968   avgTimes.reserve(numChares);
969   avgTimes.length() = numChares;
970   for(int i = 0; i < numChares; i++){
971     avgTimes[i].reserve(types);
972     avgTimes[i].length() = types;
973   }
974   initialized_instrument = true;
975 }
976
977 bool initializedInstrument(){
978   return initialized_instrument;
979 }
980
981 RequestTimeInfo &hapi_queryInstrument(int chare, char type, char phase){
982   return avgTimes[chare][type][phase];
983 }
984
985 void hapi_clearInstrument(){
986   for(int chare = 0; chare < avgTimes.length(); chare++){
987     for(int type = 0; type < avgTimes[chare].length(); type++){
988       for(int phase = 0; phase < avgTimes[chare][type].length(); phase++){
989         avgTimes[chare][type][phase].transferTime = 0.0;
990         avgTimes[chare][type][phase].kernelTime = 0.0;
991         avgTimes[chare][type][phase].cleanupTime = 0.0;
992         avgTimes[chare][type][phase].n = 0;
993       }
994       avgTimes[chare][type].length() = 0;
995     }
996     avgTimes[chare].length() = 0;
997   }
998   avgTimes.length() = 0;
999   initialized_instrument = false;
1000 }
1001 #endif