nw.c 15 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420
  1. #ifdef RD_WG_SIZE_0_0
  2. #define BLOCK_SIZE RD_WG_SIZE_0_0
  3. #elif defined(RD_WG_SIZE_0)
  4. #define BLOCK_SIZE RD_WG_SIZE_0
  5. #elif defined(RD_WG_SIZE)
  6. #define BLOCK_SIZE RD_WG_SIZE
  7. #else
  8. #define BLOCK_SIZE 16
  9. #endif
  10. #define LIMIT -999
  11. #include <stdio.h>
  12. #include <string.h>
  13. #include <stdlib.h>
  14. #include <iostream>
  15. #include <string>
  16. #include <sys/time.h>
  17. #ifdef NV //NVIDIA
  18. #include <oclUtils.h>
  19. #else
  20. #include <CL/cl.h>
  21. #endif
  22. //global variables
  23. int blosum62[24][24] = {
  24. { 4, -1, -2, -2, 0, -1, -1, 0, -2, -1, -1, -1, -1, -2, -1, 1, 0, -3, -2, 0, -2, -1, 0, -4},
  25. {-1, 5, 0, -2, -3, 1, 0, -2, 0, -3, -2, 2, -1, -3, -2, -1, -1, -3, -2, -3, -1, 0, -1, -4},
  26. {-2, 0, 6, 1, -3, 0, 0, 0, 1, -3, -3, 0, -2, -3, -2, 1, 0, -4, -2, -3, 3, 0, -1, -4},
  27. {-2, -2, 1, 6, -3, 0, 2, -1, -1, -3, -4, -1, -3, -3, -1, 0, -1, -4, -3, -3, 4, 1, -1, -4},
  28. { 0, -3, -3, -3, 9, -3, -4, -3, -3, -1, -1, -3, -1, -2, -3, -1, -1, -2, -2, -1, -3, -3, -2, -4},
  29. {-1, 1, 0, 0, -3, 5, 2, -2, 0, -3, -2, 1, 0, -3, -1, 0, -1, -2, -1, -2, 0, 3, -1, -4},
  30. {-1, 0, 0, 2, -4, 2, 5, -2, 0, -3, -3, 1, -2, -3, -1, 0, -1, -3, -2, -2, 1, 4, -1, -4},
  31. { 0, -2, 0, -1, -3, -2, -2, 6, -2, -4, -4, -2, -3, -3, -2, 0, -2, -2, -3, -3, -1, -2, -1, -4},
  32. {-2, 0, 1, -1, -3, 0, 0, -2, 8, -3, -3, -1, -2, -1, -2, -1, -2, -2, 2, -3, 0, 0, -1, -4},
  33. {-1, -3, -3, -3, -1, -3, -3, -4, -3, 4, 2, -3, 1, 0, -3, -2, -1, -3, -1, 3, -3, -3, -1, -4},
  34. {-1, -2, -3, -4, -1, -2, -3, -4, -3, 2, 4, -2, 2, 0, -3, -2, -1, -2, -1, 1, -4, -3, -1, -4},
  35. {-1, 2, 0, -1, -3, 1, 1, -2, -1, -3, -2, 5, -1, -3, -1, 0, -1, -3, -2, -2, 0, 1, -1, -4},
  36. {-1, -1, -2, -3, -1, 0, -2, -3, -2, 1, 2, -1, 5, 0, -2, -1, -1, -1, -1, 1, -3, -1, -1, -4},
  37. {-2, -3, -3, -3, -2, -3, -3, -3, -1, 0, 0, -3, 0, 6, -4, -2, -2, 1, 3, -1, -3, -3, -1, -4},
  38. {-1, -2, -2, -1, -3, -1, -1, -2, -2, -3, -3, -1, -2, -4, 7, -1, -1, -4, -3, -2, -2, -1, -2, -4},
  39. { 1, -1, 1, 0, -1, 0, 0, 0, -1, -2, -2, 0, -1, -2, -1, 4, 1, -3, -2, -2, 0, 0, 0, -4},
  40. { 0, -1, 0, -1, -1, -1, -1, -2, -2, -1, -1, -1, -1, -2, -1, 1, 5, -2, -2, 0, -1, -1, 0, -4},
  41. {-3, -3, -4, -4, -2, -2, -3, -2, -2, -3, -2, -3, -1, 1, -4, -3, -2, 11, 2, -3, -4, -3, -2, -4},
  42. {-2, -2, -2, -3, -2, -1, -2, -3, 2, -1, -1, -2, -1, 3, -3, -2, -2, 2, 7, -1, -3, -2, -1, -4},
  43. { 0, -3, -3, -3, -1, -2, -2, -3, -3, 3, 1, -2, 1, -1, -2, -2, 0, -3, -1, 4, -3, -2, -1, -4},
  44. {-2, -1, 3, 4, -3, 0, 1, -1, 0, -3, -4, 0, -3, -3, -2, 0, -1, -4, -3, -3, 4, 1, -1, -4},
  45. {-1, 0, 0, 1, -3, 3, 4, -2, 0, -3, -3, 1, -1, -3, -1, 0, -1, -3, -2, -2, 1, 4, -1, -4},
  46. { 0, -1, -1, -1, -2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -2, 0, 0, -2, -1, -1, -1, -1, -1, -4},
  47. {-4, -4, -4, -4, -4, -4, -4, -4, -4, -4, -4, -4, -4, -4, -4, -4, -4, -4, -4, -4, -4, -4, -4, 1}
  48. };
  49. // local variables
  50. static cl_context context;
  51. static cl_command_queue cmd_queue;
  52. static cl_device_type device_type;
  53. static cl_device_id * device_list;
  54. static cl_int num_devices;
  55. static int initialize(int use_gpu)
  56. {
  57. cl_int result;
  58. size_t size;
  59. // modification to handle the case in which we have more than one OpenCL platform available on the system.
  60. cl_uint platformCount;
  61. // create OpenCL context
  62. clGetPlatformIDs(0, NULL, &platformCount);
  63. cl_platform_id *platforms_ids;
  64. platforms_ids = (cl_platform_id*) malloc(sizeof(cl_platform_id) * platformCount);
  65. if (clGetPlatformIDs(platformCount, platforms_ids, NULL) != CL_SUCCESS) { printf("ERROR: clGetPlatformIDs(1,*,0) failed\n"); return -1; }
  66. cl_context_properties ctxprop[] = { CL_CONTEXT_PLATFORM, (cl_context_properties)platforms_ids[1], 0};
  67. device_type = use_gpu ? CL_DEVICE_TYPE_GPU : CL_DEVICE_TYPE_CPU;
  68. context = clCreateContextFromType( ctxprop, device_type, NULL, NULL, NULL );
  69. if( !context ) { printf("ERROR: clCreateContextFromType(%s) failed\n", use_gpu ? "GPU" : "CPU"); return -1; }
  70. // get the list of GPUs
  71. result = clGetContextInfo( context, CL_CONTEXT_DEVICES, 0, NULL, &size );
  72. num_devices = (int) (size / sizeof(cl_device_id));
  73. printf("num_devices = %d\n", num_devices);
  74. if( result != CL_SUCCESS || num_devices < 1 ) { printf("ERROR: clGetContextInfo() failed\n"); return -1; }
  75. device_list = new cl_device_id[num_devices];
  76. if( !device_list ) { printf("ERROR: new cl_device_id[] failed\n"); return -1; }
  77. result = clGetContextInfo( context, CL_CONTEXT_DEVICES, size, device_list, NULL );
  78. if( result != CL_SUCCESS ) { printf("ERROR: clGetContextInfo() failed\n"); return -1; }
  79. // create command queue for the first device
  80. cmd_queue = clCreateCommandQueue( context, device_list[0], 0, NULL );
  81. if( !cmd_queue ) { printf("ERROR: clCreateCommandQueue() failed\n"); return -1; }
  82. return 0;
  83. }
  84. static int shutdown()
  85. {
  86. // release resources
  87. if( cmd_queue ) clReleaseCommandQueue( cmd_queue );
  88. if( context ) clReleaseContext( context );
  89. if( device_list ) delete device_list;
  90. // reset all variables
  91. cmd_queue = 0;
  92. context = 0;
  93. device_list = 0;
  94. num_devices = 0;
  95. device_type = 0;
  96. return 0;
  97. }
  98. int maximum( int a,
  99. int b,
  100. int c){
  101. int k;
  102. if( a <= b )
  103. k = b;
  104. else
  105. k = a;
  106. if( k <=c )
  107. return(c);
  108. else
  109. return(k);
  110. }
  111. void usage(int argc, char **argv)
  112. {
  113. fprintf(stderr, "Usage: %s <max_rows/max_cols> <penalty> \n", argv[0]);
  114. fprintf(stderr, "\t<dimension> - x and y dimensions\n");
  115. fprintf(stderr, "\t<penalty> - penalty(positive integer)\n");
  116. fprintf(stderr, "\t<file> - filename\n");
  117. exit(1);
  118. }
  119. double gettime() {
  120. struct timeval t;
  121. gettimeofday(&t,NULL);
  122. return t.tv_sec+t.tv_usec*1e-6;
  123. }
  124. int main(int argc, char **argv){
  125. printf("WG size of kernel = %d \n", BLOCK_SIZE);
  126. int max_rows, max_cols, penalty;
  127. char * tempchar;
  128. // the lengths of the two sequences should be able to divided by 16.
  129. // And at current stage max_rows needs to equal max_cols
  130. if (argc == 4)
  131. {
  132. max_rows = atoi(argv[1]);
  133. max_cols = atoi(argv[1]);
  134. penalty = atoi(argv[2]);
  135. tempchar = argv[3];
  136. }
  137. else{
  138. usage(argc, argv);
  139. }
  140. if(atoi(argv[1])%16!=0){
  141. fprintf(stderr,"The dimension values must be a multiple of 16\n");
  142. exit(1);
  143. }
  144. max_rows = max_rows + 1;
  145. max_cols = max_cols + 1;
  146. int *reference;
  147. int *input_itemsets;
  148. int *output_itemsets;
  149. reference = (int *)malloc( max_rows * max_cols * sizeof(int) );
  150. input_itemsets = (int *)malloc( max_rows * max_cols * sizeof(int) );
  151. output_itemsets = (int *)malloc( max_rows * max_cols * sizeof(int) );
  152. srand(7);
  153. //initialization
  154. for (int i = 0 ; i < max_cols; i++){
  155. for (int j = 0 ; j < max_rows; j++){
  156. input_itemsets[i*max_cols+j] = 0;
  157. }
  158. }
  159. for( int i=1; i< max_rows ; i++){ //initialize the cols
  160. input_itemsets[i*max_cols] = rand() % 10 + 1;
  161. }
  162. for( int j=1; j< max_cols ; j++){ //initialize the rows
  163. input_itemsets[j] = rand() % 10 + 1;
  164. }
  165. for (int i = 1 ; i < max_cols; i++){
  166. for (int j = 1 ; j < max_rows; j++){
  167. reference[i*max_cols+j] = blosum62[input_itemsets[i*max_cols]][input_itemsets[j]];
  168. }
  169. }
  170. for( int i = 1; i< max_rows ; i++)
  171. input_itemsets[i*max_cols] = -i * penalty;
  172. for( int j = 1; j< max_cols ; j++)
  173. input_itemsets[j] = -j * penalty;
  174. int sourcesize = 1024*1024;
  175. char * source = (char *)calloc(sourcesize, sizeof(char));
  176. if(!source) { printf("ERROR: calloc(%d) failed\n", sourcesize); return -1; }
  177. // read the kernel core source
  178. char * kernel_nw1 = "nw_kernel1";
  179. char * kernel_nw2 = "nw_kernel2";
  180. FILE * fp = fopen(tempchar, "rb");
  181. if(!fp) { printf("ERROR: unable to open '%s'\n", tempchar); return -1; }
  182. fread(source + strlen(source), sourcesize, 1, fp);
  183. fclose(fp);
  184. int nworkitems, workgroupsize = 0;
  185. nworkitems = BLOCK_SIZE;
  186. if(nworkitems < 1 || workgroupsize < 0){
  187. printf("ERROR: invalid or missing <num_work_items>[/<work_group_size>]\n");
  188. return -1;
  189. }
  190. // set global and local workitems
  191. size_t local_work[3] = { (workgroupsize>0)?workgroupsize:1, 1, 1 };
  192. size_t global_work[3] = { nworkitems, 1, 1 }; //nworkitems = no. of GPU threads
  193. int use_gpu = 0;
  194. // OpenCL initialization
  195. if(initialize(use_gpu)) return -1;
  196. // compile kernel
  197. cl_int err = 0;
  198. const char * slist[2] = { source, 0 };
  199. cl_program prog = clCreateProgramWithSource(context, 1, slist, NULL, &err);
  200. if(err != CL_SUCCESS) { printf("ERROR: clCreateProgramWithSource() => %d\n", err); return -1; }
  201. char clOptions[110];
  202. // sprintf(clOptions,"-I../../src");
  203. sprintf(clOptions," ");
  204. #ifdef BLOCK_SIZE
  205. sprintf(clOptions + strlen(clOptions), " -DBLOCK_SIZE=%d", BLOCK_SIZE);
  206. #endif
  207. err = clBuildProgram(prog, 0, NULL, clOptions, NULL, NULL);
  208. /*{ // show warnings/errors
  209. static char log[65536]; memset(log, 0, sizeof(log));
  210. cl_device_id device_id = 0;
  211. err = clGetContextInfo(context, CL_CONTEXT_DEVICES, sizeof(device_id), &device_id, NULL);
  212. clGetProgramBuildInfo(prog, device_id, CL_PROGRAM_BUILD_LOG, sizeof(log)-1, log, NULL);
  213. if(err || strstr(log,"warning:") || strstr(log, "error:")) printf("<<<<\n%s\n>>>>\n", log);
  214. }*/
  215. if(err != CL_SUCCESS) { printf("ERROR: clBuildProgram() => %d\n", err); return -1; }
  216. cl_kernel kernel1;
  217. cl_kernel kernel2;
  218. kernel1 = clCreateKernel(prog, kernel_nw1, &err);
  219. kernel2 = clCreateKernel(prog, kernel_nw2, &err);
  220. if(err != CL_SUCCESS) { printf("ERROR: clCreateKernel() 0 => %d\n", err); return -1; }
  221. clReleaseProgram(prog);
  222. // creat buffers
  223. cl_mem input_itemsets_d;
  224. cl_mem output_itemsets_d;
  225. cl_mem reference_d;
  226. input_itemsets_d = clCreateBuffer(context, CL_MEM_READ_WRITE, max_cols * max_rows * sizeof(int), NULL, &err );
  227. if(err != CL_SUCCESS) { printf("ERROR: clCreateBuffer input_item_set (size:%d) => %d\n", max_cols * max_rows, err); return -1;}
  228. reference_d = clCreateBuffer(context, CL_MEM_READ_WRITE, max_cols * max_rows * sizeof(int), NULL, &err );
  229. if(err != CL_SUCCESS) { printf("ERROR: clCreateBuffer reference (size:%d) => %d\n", max_cols * max_rows, err); return -1;}
  230. output_itemsets_d = clCreateBuffer(context, CL_MEM_READ_WRITE, max_cols * max_rows * sizeof(int), NULL, &err );
  231. if(err != CL_SUCCESS) { printf("ERROR: clCreateBuffer output_item_set (size:%d) => %d\n", max_cols * max_rows, err); return -1;}
  232. //write buffers
  233. err = clEnqueueWriteBuffer(cmd_queue, input_itemsets_d, 1, 0, max_cols * max_rows * sizeof(int), input_itemsets, 0, 0, 0);
  234. if(err != CL_SUCCESS) { printf("ERROR: clEnqueueWriteBuffer bufIn1 (size:%d) => %d\n", max_cols * max_rows, err); return -1; }
  235. err = clEnqueueWriteBuffer(cmd_queue, reference_d, 1, 0, max_cols * max_rows * sizeof(int), reference, 0, 0, 0);
  236. if(err != CL_SUCCESS) { printf("ERROR: clEnqueueWriteBuffer bufIn2 (size:%d) => %d\n", max_cols * max_rows, err); return -1; }
  237. int worksize = max_cols - 1;
  238. printf("worksize = %d\n", worksize);
  239. //these two parameters are for extension use, don't worry about it.
  240. int offset_r = 0, offset_c = 0;
  241. int block_width = worksize/BLOCK_SIZE ;
  242. clSetKernelArg(kernel1, 0, sizeof(void *), (void*) &reference_d);
  243. clSetKernelArg(kernel1, 1, sizeof(void *), (void*) &input_itemsets_d);
  244. clSetKernelArg(kernel1, 2, sizeof(void *), (void*) &output_itemsets_d);
  245. clSetKernelArg(kernel1, 3, sizeof(cl_int) * (BLOCK_SIZE + 1) *(BLOCK_SIZE+1), (void*)NULL );
  246. clSetKernelArg(kernel1, 4, sizeof(cl_int) * BLOCK_SIZE * BLOCK_SIZE, (void*)NULL );
  247. clSetKernelArg(kernel1, 5, sizeof(cl_int), (void*) &max_cols);
  248. clSetKernelArg(kernel1, 6, sizeof(cl_int), (void*) &penalty);
  249. clSetKernelArg(kernel1, 8, sizeof(cl_int), (void*) &block_width);
  250. clSetKernelArg(kernel1, 9, sizeof(cl_int), (void*) &worksize);
  251. clSetKernelArg(kernel1, 10, sizeof(cl_int), (void*) &offset_r);
  252. clSetKernelArg(kernel1, 11, sizeof(cl_int), (void*) &offset_c);
  253. clSetKernelArg(kernel2, 0, sizeof(void *), (void*) &reference_d);
  254. clSetKernelArg(kernel2, 1, sizeof(void *), (void*) &input_itemsets_d);
  255. clSetKernelArg(kernel2, 2, sizeof(void *), (void*) &output_itemsets_d);
  256. clSetKernelArg(kernel2, 3, sizeof(cl_int) * (BLOCK_SIZE + 1) *(BLOCK_SIZE+1), (void*)NULL );
  257. clSetKernelArg(kernel2, 4, sizeof(cl_int) * BLOCK_SIZE *BLOCK_SIZE, (void*)NULL );
  258. clSetKernelArg(kernel2, 5, sizeof(cl_int), (void*) &max_cols);
  259. clSetKernelArg(kernel2, 6, sizeof(cl_int), (void*) &penalty);
  260. clSetKernelArg(kernel2, 8, sizeof(cl_int), (void*) &block_width);
  261. clSetKernelArg(kernel2, 9, sizeof(cl_int), (void*) &worksize);
  262. clSetKernelArg(kernel2, 10, sizeof(cl_int), (void*) &offset_r);
  263. clSetKernelArg(kernel2, 11, sizeof(cl_int), (void*) &offset_c);
  264. printf("Processing upper-left matrix\n");
  265. for( int blk = 1 ; blk <= worksize/BLOCK_SIZE ; blk++){
  266. global_work[0] = BLOCK_SIZE * blk;
  267. local_work[0] = BLOCK_SIZE;
  268. clSetKernelArg(kernel1, 7, sizeof(cl_int), (void*) &blk);
  269. err = clEnqueueNDRangeKernel(cmd_queue, kernel1, 2, NULL, global_work, local_work, 0, 0, 0);
  270. if(err != CL_SUCCESS) { printf("ERROR: 1 clEnqueueNDRangeKernel()=>%d failed\n", err); return -1; }
  271. }
  272. clFinish(cmd_queue);
  273. printf("Processing lower-right matrix\n");
  274. for( int blk = worksize/BLOCK_SIZE - 1 ; blk >= 1 ; blk--){
  275. global_work[0] = BLOCK_SIZE * blk;
  276. local_work[0] = BLOCK_SIZE;
  277. clSetKernelArg(kernel2, 7, sizeof(cl_int), (void*) &blk);
  278. err = clEnqueueNDRangeKernel(cmd_queue, kernel2, 2, NULL, global_work, local_work, 0, 0, 0);
  279. if(err != CL_SUCCESS) { printf("ERROR: 2 clEnqueueNDRangeKernel()=>%d failed\n", err); return -1; }
  280. }
  281. // Lingjie Zhang modified at Nov 1, 2015
  282. // clFinish(cmd_queue);
  283. // fflush(stdout);
  284. //end Lingjie Zhang modification
  285. err = clEnqueueReadBuffer(cmd_queue, input_itemsets_d, 1, 0, max_cols * max_rows * sizeof(int), output_itemsets, 0, 0, 0);
  286. clFinish(cmd_queue);
  287. //#define TRACEBACK
  288. #ifdef TRACEBACK
  289. FILE *fpo = fopen("result.txt","w");
  290. fprintf(fpo, "print traceback value GPU:\n");
  291. for (int i = max_rows - 2, j = max_rows - 2; i>=0, j>=0;){
  292. int nw, n, w, traceback;
  293. if ( i == max_rows - 2 && j == max_rows - 2 )
  294. fprintf(fpo, "%d ", output_itemsets[ i * max_cols + j]); //print the first element
  295. if ( i == 0 && j == 0 )
  296. break;
  297. if ( i > 0 && j > 0 ){
  298. nw = output_itemsets[(i - 1) * max_cols + j - 1];
  299. w = output_itemsets[ i * max_cols + j - 1 ];
  300. n = output_itemsets[(i - 1) * max_cols + j];
  301. }
  302. else if ( i == 0 ){
  303. nw = n = LIMIT;
  304. w = output_itemsets[ i * max_cols + j - 1 ];
  305. }
  306. else if ( j == 0 ){
  307. nw = w = LIMIT;
  308. n = output_itemsets[(i - 1) * max_cols + j];
  309. }
  310. else{
  311. }
  312. //traceback = maximum(nw, w, n);
  313. int new_nw, new_w, new_n;
  314. new_nw = nw + reference[i * max_cols + j];
  315. new_w = w - penalty;
  316. new_n = n - penalty;
  317. traceback = maximum(new_nw, new_w, new_n);
  318. if(traceback == new_nw)
  319. traceback = nw;
  320. if(traceback == new_w)
  321. traceback = w;
  322. if(traceback == new_n)
  323. traceback = n;
  324. fprintf(fpo, "%d ", traceback);
  325. if(traceback == nw )
  326. {i--; j--; continue;}
  327. else if(traceback == w )
  328. {j--; continue;}
  329. else if(traceback == n )
  330. {i--; continue;}
  331. else
  332. ;
  333. }
  334. fclose(fpo);
  335. #endif
  336. printf("Computation Done\n");
  337. // OpenCL shutdown
  338. if(shutdown()) return -1;
  339. clReleaseMemObject(input_itemsets_d);
  340. clReleaseMemObject(output_itemsets_d);
  341. clReleaseMemObject(reference_d);
  342. free(reference);
  343. free(input_itemsets);
  344. free(output_itemsets);
  345. }