← Back to branch summary

~ubuntu-branches/ubuntu/precise/boinc/precise

« back to all changes in this revision

Viewing changes to samples/atiopencl/atiopencl.cpp

Tags: 6.12.8+dfsg-1
* New upstream release.
* Simplified debian/rules

Show diffs side-by-side

added added

removed removed

Lines of Context:
samples/atiopencl/atiopencl.cpp
 
1
// This file is part of BOINC.
 
2
// http://boinc.berkeley.edu
 
3
// Copyright (C) 2008 University of California
 
4
//
 
5
// BOINC is free software; you can redistribute it and/or modify it
 
6
// under the terms of the GNU Lesser General Public License
 
7
// as published by the Free Software Foundation,
 
8
// either version 3 of the License, or (at your option) any later version.
 
9
//
 
10
// BOINC is distributed in the hope that it will be useful,
 
11
// but WITHOUT ANY WARRANTY; without even the implied warranty of
 
12
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
 
13
// See the GNU Lesser General Public License for more details.
 
14
//
 
15
// You should have received a copy of the GNU Lesser General Public License
 
16
// along with BOINC.  If not, see <http://www.gnu.org/licenses/>.
 
17
//
 
18
// This program serves as both
 
19
// - An example BOINC-ATIOpenCL application, illustrating the use of the BOINC API
 
20
//   and ATIStream OpenCL API.
 
21
// - A program for testing various features of BOINC.
 
22
//
 
23
// The program reads the input nxn matrix from the "input" file, inverts the
 
24
// matrix NUM_ITERATIONS times and write to "output" file.
 
25
//
 
26
// command line options
 
27
// -run_slow: sleep 1 second after each character
 
28
// -cpu_time N: use about N CPU seconds after copying files
 
29
// -early_exit: exit(10) after 30 chars
 
30
// -early_crash: crash after 30 chars
 
31
//
 
32
// See http://boinc.berkeley.edu/trac/wiki/GPUApp for any compiling issues.
 
33
// Contributor: Tuan Le (tuanle86@berkeley.edu)
 
34
 
 
35
#include "atiopencl.hpp"
 
36
using std::string;
 
37
 
 
38
int main(int argc, char * argv[]) {
 
39
    int i, retval, lastInversion=0, checkpointExists=0, matrixSize=0;
 
40
    double fd;
 
41
    char input_path[512], output_path[512], chkpt_path[512], buf[256];
 
42
    MFILE out;
 
43
    FILE* state, *infile;
 
44
 
 
45
    generate_random_input_file(MATRIX_SIZE); //call this if you don't want to
 
46
                                             //construct the input file manually
 
47
    
 
48
        for (i=0; i<argc; i++) {
 
49
        if (!strcmp(argv[i], "-early_exit")) early_exit = true;
 
50
        if (!strcmp(argv[i], "-early_crash")) early_crash = true;
 
51
        if (!strcmp(argv[i], "-early_sleep")) early_sleep = true;
 
52
        if (!strcmp(argv[i], "-run_slow")) run_slow = true;
 
53
        if (!strcmp(argv[i], "-cpu_time")) {
 
54
            cpu_time = atof(argv[++i]);
 
55
        }
 
56
    }
 
57
        
 
58
    retval = boinc_init();
 
59
    if (retval) {
 
60
        fprintf(stderr, "%s boinc_init returned %d\n",
 
61
                boinc_msg_prefix(buf, sizeof(buf)), retval );
 
62
        exit(retval);
 
63
    }
 
64
    
 
65
    // open the input file (resolve logical name first)
 
66
    //
 
67
    boinc_resolve_filename(INPUT_FILENAME, input_path, sizeof(input_path));
 
68
    infile = boinc_fopen(input_path, "r");
 
69
    if (!infile) {
 
70
        fprintf(stderr,
 
71
            "%s Couldn't find input file in boinc\\win_build, resolved name %s.\n",
 
72
            boinc_msg_prefix(buf, sizeof(buf)), input_path
 
73
        );
 
74
        getchar();
 
75
        exit(-1);
 
76
    }
 
77
    
 
78
    boinc_resolve_filename(OUTPUT_FILENAME, output_path, sizeof(output_path));
 
79
    
 
80
    // See if there's a valid checkpoint file.
 
81
    // If so retrieve the current matrix and inversion number
 
82
    //
 
83
    boinc_resolve_filename(CHECKPOINT_FILE, chkpt_path, sizeof(chkpt_path));
 
84
    state = boinc_fopen(chkpt_path, "r");
 
85
    if (state) {
 
86
        printf("Checkpoint file is detected. Read from checkpoint file ... \n");
 
87
        checkpointExists=fscanf(state, "%d", &lastInversion); 
 
88
        if (checkpointExists == 1) {
 
89
            isStateFileInUse=true;
 
90
            printf("Last inversion # is : %d\n",lastInversion); 
 
91
            fscanf(state,"%d",&matrixSize);
 
92
            width=height=matrixSize;
 
93
            printf("Initialize host ....\n");
 
94
            initialize_host(state);
 
95
        }
 
96
        fclose(state);
 
97
    } else {
 
98
        printf("There's no valid checkpoint file!\n");
 
99
    }
 
100
 
 
101
    retval = out.open(output_path, "wb");
 
102
 
 
103
    if (retval) {
 
104
        fprintf(stderr, "%s APP: matrix_inversion output open failed:\n",
 
105
            boinc_msg_prefix(buf, sizeof(buf))
 
106
        );
 
107
        fprintf(stderr, "%s resolved name %s, retval %d\n",
 
108
            boinc_msg_prefix(buf, sizeof(buf)), output_path, retval
 
109
        );
 
110
        perror("open");
 
111
        exit(1);
 
112
    }
 
113
 
 
114
#ifdef APP_GRAPHICS
 
115
    // create shared mem segment for graphics, and arrange to update it
 
116
    //
 
117
    shmem = (UC_SHMEM*)boinc_graphics_make_shmem("matrix_inversion", sizeof(UC_SHMEM));
 
118
    if (!shmem) {
 
119
        fprintf(stderr, "%s failed to create shared mem segment\n",
 
120
            boinc_msg_prefix(buf, sizeof(buf))
 
121
        );
 
122
    }
 
123
    update_shmem();
 
124
    boinc_register_timer_callback(update_shmem);
 
125
#endif
 
126
 
 
127
    if (checkpointExists != 1) { //checkpoint file is not found.
 
128
        matrixSize=get_matrix_size(infile);
 
129
        printf("Matrix Size: width = height = %d\n",matrixSize);
 
130
        width=height=matrixSize;
 
131
        // Initialize Host application
 
132
        printf("Initialize host ....\n");
 
133
        if (initialize_host(infile)==1) {
 
134
            return 1;   
 
135
        }
 
136
        out.printf("\n----------------- Before being inversed ----------------\n\n");
 
137
        printf("Computation is running ... Inverse the matrix %d times. Start at inversion #1\n",
 
138
               NUM_ITERATIONS);
 
139
    } else {
 
140
        out.printf("\n----------------- Last checkpointed inversion #%d ----------------\n\n",
 
141
                   lastInversion);
 
142
        printf("Computation is resumed ... Inverse the matrix %d more times. Start at inversion #%d\n",
 
143
               NUM_ITERATIONS-lastInversion,lastInversion+1);
 
144
    }
 
145
 
 
146
    // Initialize OpenCL resources
 
147
    if (initialize_cl()==1) {
 
148
        return 1;
 
149
    }
 
150
 
 
151
    print_to_file(&out,input,matrixSize);
 
152
 
 
153
    for (int i=lastInversion+1;i<=NUM_ITERATIONS;++i) {
 
154
        //the invert function will trigger kernel calls.
 
155
        invert(input,output,matrixSize);
 
156
        printf("Finish inversion #%d\n",i);
 
157
        for (int j=0;j<matrixSize*matrixSize;++j) {
 
158
            input[j]=output[j]; //change the input for the next iteration
 
159
        }
 
160
        if (run_slow) {
 
161
            boinc_sleep(1.);
 
162
        }
 
163
 
 
164
        if (early_exit && i>30) {
 
165
            exit(-10);
 
166
        }
 
167
 
 
168
        if (early_crash && i>30) {
 
169
            boinc_crash();
 
170
        }
 
171
 
 
172
        if (early_sleep && i>30) {
 
173
            g_sleep = true;
 
174
            while (1) boinc_sleep(1);
 
175
        }
 
176
                
 
177
        if (boinc_time_to_checkpoint()) {
 
178
            printf("Perform checkpointing at inversion # %d\n",i);
 
179
            //we'll need to write the current matrix to the state file.
 
180
            retval = do_checkpoint(out, i, input, matrixSize); 
 
181
            if (retval) {
 
182
                fprintf(stderr, "%s APP: matrix_inversion checkpoint failed %d\n",
 
183
                    boinc_msg_prefix(buf, sizeof(buf)), retval
 
184
                );
 
185
                exit(retval);
 
186
            }
 
187
            boinc_checkpoint_completed();
 
188
        }
 
189
        fd = i/NUM_ITERATIONS;
 
190
        if (cpu_time) fd /= 2;
 
191
        boinc_fraction_done(fd);
 
192
    }
 
193
 
 
194
    out.printf("\n\n----------------- Final inversion #%d ----------------\n\n",
 
195
               NUM_ITERATIONS);
 
196
    print_to_file(&out,output,matrixSize);
 
197
 
 
198
    retval = out.flush(); //force the output file to be closed.
 
199
    if (retval) {
 
200
        fprintf(stderr, "%s APP: matrix_inversion flush failed %d\n",
 
201
            boinc_msg_prefix(buf, sizeof(buf)), retval
 
202
        );
 
203
        exit(1);
 
204
    }
 
205
 
 
206
    // Releases OpenCL resources 
 
207
    if (cleanup_cl()==1) {
 
208
        printf("Error!");
 
209
        return 1;
 
210
    }
 
211
 
 
212
    // Release host resources
 
213
    cleanup_host();
 
214
 
 
215
    // burn up some CPU time if needed
 
216
    //
 
217
    if (cpu_time) {
 
218
        printf("\nBurning up some CPU time ... \n");
 
219
        double start = dtime();
 
220
        for (int i=0; ; i++) {
 
221
            double e = dtime()-start;
 
222
            if (e > cpu_time) break;
 
223
            fd = .5 + .5*(e/cpu_time);
 
224
            boinc_fraction_done(fd);
 
225
 
 
226
            if (boinc_time_to_checkpoint()) {
 
227
                retval = do_checkpoint(out, NUM_ITERATIONS, input, matrixSize);
 
228
                if (retval) {
 
229
                    fprintf(stderr, "%s APP: maxtrix_inversion checkpoint failed %d\n",
 
230
                        boinc_msg_prefix(buf, sizeof(buf)), retval
 
231
                    );
 
232
                    exit(1);
 
233
                }
 
234
                boinc_checkpoint_completed();
 
235
            }
 
236
            comp_result = do_a_giga_flop(i);
 
237
        }
 
238
    }
 
239
    boinc_fraction_done(1);
 
240
 
 
241
#ifdef APP_GRAPHICS
 
242
    update_shmem();
 
243
#endif
 
244
 
 
245
    printf("\nDone! Please press ENTER to exit. ");
 
246
    getchar();
 
247
    boinc_finish(0);
 
248
}
 
249
 
 
250
#ifdef _WIN32
 
251
int WINAPI WinMain(HINSTANCE hInst, HINSTANCE hPrevInst, LPSTR Args, int WinMode) {
 
252
    LPSTR command_line;
 
253
    char* argv[100];
 
254
    int argc;
 
255
 
 
256
    command_line = GetCommandLine();
 
257
    argc = parse_command_line( command_line, argv );
 
258
    return main(argc, argv);
 
259
}
 
260
#endif
 
261
 
 
262
/*** BOINC FUNCTION DEFINITIONS ***/
 
263
 
 
264
/* Do a billion floating-point ops */
 
265
static double do_a_giga_flop(int foo) {
 
266
    double x = 3.14159*foo;
 
267
    int i;
 
268
    for (i=0; i<500000000; i++) {
 
269
        x += 5.12313123;
 
270
        x *= 0.5398394834;
 
271
    }
 
272
    return x;
 
273
}
 
274
 
 
275
/* Save the computation state into checkpoint file */
 
276
int do_checkpoint(MFILE& mf, int n, cl_float *input, int matrixSize) {
 
277
    int retval;
 
278
    string resolved_name;
 
279
 
 
280
    FILE* f = fopen("temp", "w");
 
281
    if (!f) return 1;
 
282
    fprintf(f, "%d", n); //write inversion number
 
283
    fprintf(f, " ");
 
284
    fprintf(f, "%d", matrixSize); //write matrixSize
 
285
    fprintf(f, " ");
 
286
    for (int i=0;i<matrixSize*matrixSize;++i) {
 
287
        fprintf(f, " ");
 
288
        fprintf(f, "%f", input[i]);
 
289
    }
 
290
    fclose(f);
 
291
    retval = mf.flush();
 
292
    if (retval) return retval;
 
293
    boinc_resolve_filename_s(CHECKPOINT_FILE, resolved_name);
 
294
    retval = boinc_rename("temp", resolved_name.c_str());
 
295
    if (retval) return retval;
 
296
    return 0; //return 0 to indicate success.
 
297
}
 
298
 
 
299
/*** FUNCTION DEFINITIONS ***/
 
300
 
 
301
/* Create an input file filled with random data of type cl_float. */
 
302
void generate_random_input_file(int n) {
 
303
    FILE *infile;
 
304
 
 
305
    infile=fopen(INPUT_FILENAME,"w");
 
306
    cl_float *input = (cl_float *)malloc(sizeof(cl_float)*(n*n));
 
307
    srand(n);
 
308
    for( int i = 0; i < n; i++ ) {
 
309
        for (int j = 0; j < n; j++) {
 
310
            input[i*n+j] = 2.0*(rand()%32768)/32768.0 - 1.0;
 
311
        }
 
312
        input[i*n+i] += sqrt((float)n);
 
313
    }
 
314
    int j=0;
 
315
    for (int i=0;i<n*n;++i) {
 
316
        fprintf(infile,"%15f",input[i]);
 
317
        if (j+1==n) {
 
318
            fprintf(infile,"\n");
 
319
            j=0;
 
320
        } else {
 
321
            ++j;
 
322
        }
 
323
    }
 
324
    fclose(infile);
 
325
    free(input);
 
326
}
 
327
 
 
328
/*
 
329
 * Parse the input file and determine the size of the matrix.
 
330
 * This is an nxn matrix. Note: if width<> height, the matrix is
 
331
 * non-invertible.
 
332
 */
 
333
int get_matrix_size(FILE *infile) {
 
334
    int w=0;
 
335
    char c;
 
336
        
 
337
    fseek(infile,0,SEEK_SET);
 
338
    while (true) {
 
339
        do {
 
340
            c=fgetc(infile);
 
341
            if (c == EOF || c == '\n') {
 
342
                goto exitLoop;
 
343
            }
 
344
        } while (isspace(c));
 
345
 
 
346
        if (isdigit(c) || c=='.' || c=='-') {
 
347
            ++w;
 
348
        }
 
349
 
 
350
        do {
 
351
            c=fgetc(infile);
 
352
            if (c == EOF || c == '\n') {
 
353
                goto exitLoop;
 
354
            }
 
355
        } while (isdigit(c) || c=='.' || c=='-');
 
356
 
 
357
        if (c==EOF || c == '\n') {
 
358
            break;
 
359
        }
 
360
    }
 
361
    exitLoop:
 
362
    return w;
 
363
}
 
364
 
 
365
/*
 
366
 * \brief Host Initialization 
 
367
 *        Allocate and initialize memory 
 
368
 *        on the host. Print input array. 
 
369
 */
 
370
int initialize_host(FILE *infile) {
 
371
    input  = NULL;
 
372
    output = NULL;
 
373
 
 
374
    if (width!=height) {
 
375
        printf("Error: non nxn matrix cannot be invertiable.\n");
 
376
        return 1;
 
377
    }
 
378
 
 
379
    /////////////////////////////////////////////////////////////////
 
380
    // Allocate and initialize memory used by host 
 
381
    /////////////////////////////////////////////////////////////////
 
382
    cl_uint sizeInBytes = width * height * sizeof(cl_float);
 
383
    input = (cl_float *) malloc(sizeInBytes);
 
384
    if (input == NULL) {
 
385
        printf("Error: Failed to allocate input memory on host\n");
 
386
        return 1; 
 
387
    }
 
388
 
 
389
    output = (cl_float *) malloc(sizeInBytes);
 
390
    if(output == NULL) {
 
391
        printf("Error: Failed to allocate output memory on host\n");
 
392
        return 1; 
 
393
    }
 
394
 
 
395
    //fillRandom(input,width,height);
 
396
    fetch_elements_into_host_memory(infile,input);
 
397
    return 0;
 
398
}
 
399
 
 
400
/*
 
401
 * Read the float values from input file into "input" array.
 
402
 */
 
403
void fetch_elements_into_host_memory(FILE *infile, cl_float *input) {
 
404
    float num=0;
 
405
    int i=0;
 
406
    if (!isStateFileInUse) {
 
407
        fseek(infile,0,SEEK_SET);
 
408
    }
 
409
    while (fscanf(infile,"%f",&num)==1) {
 
410
        input[i]=num;
 
411
        ++i;
 
412
    }
 
413
}
 
414
 
 
415
/*
 
416
 * Converts the contents of a file into a string
 
417
 */
 
418
char * convert_to_string(const char *fileName) {
 
419
    int count=0;
 
420
    char *s;
 
421
    char c;
 
422
    int i=0;
 
423
 
 
424
    // look for "atiopencl_kernels.cl" in "boinc/samples/atiopencl/debug" or
 
425
    // in "boinc/samples/atiopencl/release". Note that "atiopencl_kernels.cl"
 
426
    // is automatically copied to these directories along the building process.
 
427
    FILE *infile=fopen(fileName,"r");
 
428
    if (!infile) { //not found. This typically happens on Linux or Mac.
 
429
        //look for "atiopencl_kernels.cl" in "boinc/sample/atiopencl/" instead.
 
430
        infile = fopen(KERNELS_FILEPATH,"r"); 
 
431
        if (!infile) {
 
432
            printf("File open Error!");
 
433
            exit(0);
 
434
        }
 
435
    }
 
436
    fseek(infile,0,SEEK_SET);
 
437
    while (fgetc(infile)!=EOF) count++;
 
438
    s=(char *) malloc(sizeof(char)*(count+1)); //add 1 for string terminator.
 
439
    fseek(infile,0,SEEK_SET);   
 
440
    while ((c=fgetc(infile))!=EOF) {
 
441
        s[i++]=c;
 
442
    }
 
443
    s[i]='\0';
 
444
    return s;
 
445
}
 
446
 
 
447
/*
 
448
 * \brief OpenCL related initialization 
 
449
 *        Create Context, Device list, Command Queue
 
450
 *        Load CL file, compile, link CL source 
 
451
 *                Build program and kernel objects
 
452
 */
 
453
 
 
454
 // Note: OpenCL memory buffer objects will be created in invert
 
455
 //       function before kernel calls are made.
 
456
int initialize_cl(void) {
 
457
    cl_int status = 0;
 
458
    size_t deviceListSize;
 
459
 
 
460
    localThreads[0]  = LOCAL_WORK_SIZE;
 
461
    globalThreads[0] = GLOBAL_WORK_SIZE;
 
462
 
 
463
    /*
 
464
     * Have a look at the available platforms and pick either
 
465
     * the AMD one if available or a reasonable default.
 
466
     */
 
467
 
 
468
    cl_uint numPlatforms;
 
469
    cl_platform_id platform = NULL;
 
470
    status = clGetPlatformIDs(0, NULL, &numPlatforms);
 
471
    if(status != CL_SUCCESS) {
 
472
        printf("Error: Getting Platforms. (clGetPlatformsIDs)\n");
 
473
        return 1;
 
474
    }
 
475
    
 
476
    if (numPlatforms > 0) {
 
477
        cl_platform_id* platforms = (cl_platform_id *)
 
478
                                                    malloc(sizeof(cl_platform_id)*numPlatforms);
 
479
        status = clGetPlatformIDs(numPlatforms, platforms, NULL);
 
480
        if (status != CL_SUCCESS) {
 
481
            printf("Error: Getting Platform Ids. (clGetPlatformsIDs)\n");
 
482
            return 1;
 
483
        }
 
484
        for (unsigned int i=0; i < numPlatforms; ++i) {
 
485
            char pbuff[100];
 
486
            status = clGetPlatformInfo(platforms[i],
 
487
                                       CL_PLATFORM_VENDOR,
 
488
                                       sizeof(pbuff),
 
489
                                       pbuff,
 
490
                                       NULL);
 
491
            if (status != CL_SUCCESS) {
 
492
                printf("Error: Getting Platform Info.(clGetPlatformInfo)\n");
 
493
                return 1;
 
494
            }
 
495
            platform = platforms[i];
 
496
            if (!strcmp(pbuff, "Advanced Micro Devices, Inc.")) {
 
497
                break;
 
498
            }
 
499
        }
 
500
        delete platforms;
 
501
    }
 
502
 
 
503
    if(NULL == platform) {
 
504
        printf("NULL platform found so Exiting Application.");
 
505
        return 1;
 
506
    }
 
507
 
 
508
    /*
 
509
     * If we could find our platform, use it. Otherwise use just available platform.
 
510
     */
 
511
    cl_context_properties cps[3] = { CL_CONTEXT_PLATFORM,
 
512
                                             (cl_context_properties)platform,
 
513
                                                                         0 
 
514
                                       };
 
515
 
 
516
    /////////////////////////////////////////////////////////////////
 
517
    // Create an OpenCL context
 
518
    /////////////////////////////////////////////////////////////////
 
519
    context = clCreateContextFromType(cps, CL_DEVICE_TYPE_ALL, NULL, NULL, &status);
 
520
    if (status != CL_SUCCESS) {  
 
521
        printf("Error: Creating Context. (clCreateContextFromType)\n");
 
522
        return 1; 
 
523
    }
 
524
 
 
525
    /* First, get the size of device list data */
 
526
    status = clGetContextInfo(context, CL_CONTEXT_DEVICES, 0, NULL, &deviceListSize);
 
527
    if (status != CL_SUCCESS) {  
 
528
        printf("Error: Getting Context Info (device list size, clGetContextInfo)\n");
 
529
        return 1;
 
530
    }
 
531
 
 
532
    /////////////////////////////////////////////////////////////////
 
533
    // Detect OpenCL devices
 
534
    /////////////////////////////////////////////////////////////////
 
535
    devices = (cl_device_id *)malloc(deviceListSize);
 
536
    if (devices == 0) {
 
537
        printf("Error: No devices found.\n");
 
538
        return 1;
 
539
    }
 
540
 
 
541
    /* Now, get the device list data */
 
542
    status = clGetContextInfo(context, CL_CONTEXT_DEVICES, deviceListSize, devices, NULL);
 
543
    if (status != CL_SUCCESS) { 
 
544
        printf("Error: Getting Context Info (device list, clGetContextInfo)\n");
 
545
        return 1;
 
546
    }
 
547
 
 
548
    /////////////////////////////////////////////////////////////////
 
549
    // Create an OpenCL command queue
 
550
    /////////////////////////////////////////////////////////////////
 
551
    commandQueue = clCreateCommandQueue(context, devices[0], 0, &status);
 
552
    if(status != CL_SUCCESS) { 
 
553
        printf("Creating Command Queue. (clCreateCommandQueue)\n");
 
554
        return 1;
 
555
    }
 
556
    
 
557
    /////////////////////////////////////////////////////////////////
 
558
    // Load CL file, build CL program object, create CL kernel object
 
559
    /////////////////////////////////////////////////////////////////
 
560
    source = convert_to_string(KERNELS_FILENAME);
 
561
    size_t sourceSize[]    = { strlen(source) };
 
562
    program = clCreateProgramWithSource(context, 1, &source, sourceSize, &status);
 
563
    if (status != CL_SUCCESS) { 
 
564
        printf("Error: Loading Binary into cl_program (clCreateProgramWithBinary)\n");
 
565
        return 1;
 
566
    }
 
567
 
 
568
    /* create a cl program executable for all the devices specified */
 
569
    status = clBuildProgram(program, 1, devices, NULL, NULL, NULL);
 
570
    if (status != CL_SUCCESS)  {
 
571
        printf("Error: Building Program (clBuildProgram)\n");
 
572
        return 1; 
 
573
    }
 
574
 
 
575
    /* get a kernel object handle for a kernel with the given name */
 
576
    GEStep1A_kernel = clCreateKernel(program, "GEStep1A", &status);
 
577
    if (status != CL_SUCCESS) {  
 
578
        printf("Error: clCreateKernel (GEStep1A)\n");
 
579
        return 1;
 
580
    }
 
581
 
 
582
    GEStep2_kernel = clCreateKernel(program, "GEStep2", &status);
 
583
    if (status != CL_SUCCESS) {
 
584
        printf("Error: clCreateKernel (GEStep2)\n");
 
585
        return 1;
 
586
    }
 
587
 
 
588
        GEStep3_kernel = clCreateKernel(program, "GEStep3", &status);
 
589
    if (status != CL_SUCCESS) {
 
590
        printf("Error: clCreateKernel (GEStep3)\n");
 
591
        return 1;
 
592
    }
 
593
 
 
594
    return 0;
 
595
}
 
596
 
 
597
/*
 
598
 * \brief Release OpenCL resources (Context, Memory etc.) 
 
599
 */
 
600
int cleanup_cl(void) {
 
601
    cl_int status;
 
602
 
 
603
    status = clReleaseKernel(GEStep1A_kernel);
 
604
    if (status != CL_SUCCESS) {
 
605
        printf("Error: In clReleaseKernel (GEStep1A_kernel)\n");
 
606
        return 1; 
 
607
    }
 
608
 
 
609
    status = clReleaseKernel(GEStep2_kernel);
 
610
    if (status != CL_SUCCESS) {
 
611
        printf("Error: In clReleaseKernel (GEStep2_kernel)\n");
 
612
        return 1; 
 
613
        }
 
614
 
 
615
    status = clReleaseKernel(GEStep3_kernel);
 
616
    if (status != CL_SUCCESS) {
 
617
        printf("Error: In clReleaseKernel (GEStep3_kernel)\n");
 
618
        return 1; 
 
619
    }
 
620
 
 
621
    status = clReleaseProgram(program);
 
622
    if (status != CL_SUCCESS) {
 
623
        printf("Error: In clReleaseProgram\n");
 
624
        return 1; 
 
625
    }
 
626
 
 
627
    status = clReleaseMemObject(inputBuffer);
 
628
    if (status != CL_SUCCESS) {
 
629
        printf("Error: In clReleaseMemObject (inputBuffer)\n");
 
630
        return 1; 
 
631
    }
 
632
 
 
633
    status = clReleaseCommandQueue(commandQueue);
 
634
    if (status != CL_SUCCESS) {
 
635
        printf("Error: In clReleaseCommandQueue\n");
 
636
        return 1;
 
637
    }
 
638
 
 
639
    status = clReleaseContext(context);
 
640
    if (status != CL_SUCCESS) {
 
641
        printf("Error: In clReleaseContext\n");
 
642
        return 1;
 
643
    }
 
644
 
 
645
    return 0;
 
646
}
 
647
 
 
648
/* 
 
649
 * \brief Releases program's resources 
 
650
 */
 
651
void cleanup_host(void) {
 
652
    if (input != NULL) {
 
653
        free(input);
 
654
        input = NULL;
 
655
    }
 
656
 
 
657
    if (output != NULL) {
 
658
        free(output);
 
659
        output = NULL;
 
660
    }
 
661
 
 
662
    if (devices != NULL) {
 
663
        free(devices);
 
664
        devices = NULL;
 
665
    }
 
666
 
 
667
    if (source != NULL) {
 
668
        free((char *)source);
 
669
        source = NULL;
 
670
    }
 
671
}
 
672
 
 
673
/*
 
674
 * Write the result to output file
 
675
 */
 
676
void print_to_file(MFILE *out, float *h_odata, int n) {
 
677
    int count=0;
 
678
    int move=0;
 
679
    int num_elements=n*n;
 
680
    while (num_elements>0) {
 
681
        out->printf("%15f    ",h_odata[move]);
 
682
        ++count;
 
683
        ++move;
 
684
        if (count==n) {
 
685
            out->printf("\n");
 
686
            count=0;
 
687
        }
 
688
        --num_elements;
 
689
    }
 
690
}
 
691
 
 
692
/*
 
693
 * \brief Run OpenCL program 
 
694
 *                
 
695
 *        Bind host variables to kernel arguments 
 
696
 *                Run the CL kernel
 
697
 */
 
698
int run_GEStep1A_kernel(cl_float * AI, int i, int n2, int lda2) {
 
699
    cl_int status;
 
700
    cl_event events[2];
 
701
 
 
702
    /* 
 
703
     * the input array to the kernel. This array will eventually be modified
 
704
     * to the inverted array.
 
705
     */
 
706
    status = clSetKernelArg(GEStep1A_kernel, 0, sizeof(cl_mem), (void *)&inputBuffer);
 
707
    if (status != CL_SUCCESS) { 
 
708
        printf("Error: Setting kernel argument. (input)\n");
 
709
        return 1;
 
710
    }
 
711
 
 
712
    /*i*/
 
713
    status = clSetKernelArg(GEStep1A_kernel, 1, sizeof(int), (void *)&i);
 
714
    if (status != CL_SUCCESS) { 
 
715
        printf("Error: Setting kernel argument. (i)\n");
 
716
        return 1;
 
717
    }
 
718
 
 
719
    /*n2*/
 
720
    status = clSetKernelArg(GEStep1A_kernel, 2, sizeof(int), (void *)&n2);
 
721
    if (status != CL_SUCCESS) { 
 
722
        printf("Error: Setting kernel argument. (n2)\n");
 
723
        return 1;
 
724
    }
 
725
 
 
726
    /*lda2*/
 
727
    status = clSetKernelArg(GEStep1A_kernel, 3, sizeof(int), (void *)&lda2);
 
728
    if (status != CL_SUCCESS) { 
 
729
        printf("Error: Setting kernel argument. (lda2)\n");
 
730
        return 1;
 
731
    }
 
732
 
 
733
    /* 
 
734
     * Enqueue a kernel run call.
 
735
     */
 
736
    status = clEnqueueNDRangeKernel(commandQueue,
 
737
                                    GEStep1A_kernel,
 
738
                                    1,
 
739
                                    NULL,
 
740
                                    globalThreads,
 
741
                                    localThreads,
 
742
                                    0,
 
743
                                    NULL,
 
744
                                    &events[0]);
 
745
    if (status != CL_SUCCESS) { 
 
746
        printf("Error: Enqueueing kernel onto command queue. (clEnqueueNDRangeKernel)\n");
 
747
        return 1;
 
748
    }
 
749
 
 
750
    /* wait for the kernel call to finish execution */
 
751
    status = clWaitForEvents(1, &events[0]);
 
752
    if (status != CL_SUCCESS) { 
 
753
        printf("Error: Waiting for kernel run to finish. (clWaitForEvents)\n");
 
754
        return 1;
 
755
    }
 
756
 
 
757
    status = clReleaseEvent(events[0]);
 
758
    if (status != CL_SUCCESS) { 
 
759
        printf("Error: Release event object. (clReleaseEvent)\n");
 
760
        return 1;
 
761
    }
 
762
 
 
763
    /* Enqueue readBuffer*/  //Note: we are reading back from inputBuffer since AI is modified directly in kernel
 
764
    status = clEnqueueReadBuffer(commandQueue,
 
765
                                 inputBuffer,
 
766
                                 CL_TRUE,
 
767
                                 0,
 
768
                                 globalThreads[0] * sizeof(cl_float),
 
769
                                 AI,
 
770
                                 0,
 
771
                                 NULL,
 
772
                                 &events[1]);
 
773
 
 
774
    if(status != CL_SUCCESS) { 
 
775
        printf("Error: clEnqueueReadBuffer failed. (clEnqueueReadBuffer)\n");
 
776
        return 1;
 
777
    }
 
778
 
 
779
    /* Wait for the read buffer to finish execution */
 
780
    status = clWaitForEvents(1, &events[1]);
 
781
    if (status != CL_SUCCESS) {
 
782
        printf("Error: Waiting for read buffer call to finish. (clWaitForEvents)\n");
 
783
        return 1;
 
784
    }
 
785
 
 
786
    status = clReleaseEvent(events[1]);
 
787
    if (status != CL_SUCCESS) { 
 
788
        printf("Error: Release event object. (clReleaseEvent)\n");
 
789
        return 1;
 
790
    }
 
791
    return 0;
 
792
}
 
793
 
 
794
int run_GEStep2_kernel(cl_float * AI, cl_float diag, int i, int n2, int lda2) {
 
795
    cl_int status;
 
796
    cl_event events[2];
 
797
 
 
798
    /* 
 
799
     * the input array to the kernel. This array will eventually be modified 
 
800
     * to the inverted array.  
 
801
     */
 
802
    status = clSetKernelArg(GEStep2_kernel, 0, sizeof(cl_mem), (void *)&inputBuffer);
 
803
    if (status != CL_SUCCESS) { 
 
804
        printf("Error: Setting kernel argument. (AI)\n");
 
805
        return 1;
 
806
    }
 
807
 
 
808
    /*diag*/
 
809
    status = clSetKernelArg(GEStep2_kernel, 1, sizeof(cl_float), (void *)&diag);
 
810
    if (status != CL_SUCCESS) { 
 
811
        printf("Error: Setting kernel argument. (diag)\n");
 
812
        return 1;
 
813
    }
 
814
 
 
815
    /*i*/
 
816
    status = clSetKernelArg(GEStep2_kernel, 2, sizeof(int), (void *)&i);
 
817
    if (status != CL_SUCCESS) { 
 
818
        printf("Error: Setting kernel argument. (i)\n");
 
819
        return 1;
 
820
    }
 
821
 
 
822
    /*n2*/
 
823
    status = clSetKernelArg(GEStep2_kernel, 3, sizeof(int), (void *)&n2);
 
824
    if (status != CL_SUCCESS) { 
 
825
        printf("Error: Setting kernel argument. (n2)\n");
 
826
        return 1;
 
827
    }
 
828
 
 
829
    /*lda2*/
 
830
    status = clSetKernelArg(GEStep2_kernel, 4, sizeof(int), (void *)&lda2);
 
831
    if (status != CL_SUCCESS) {
 
832
        printf("Error: Setting kernel argument. (lda2)\n");
 
833
        return 1;
 
834
    }
 
835
 
 
836
    /* 
 
837
     * Enqueue a kernel run call.
 
838
     */
 
839
    status = clEnqueueNDRangeKernel(commandQueue,
 
840
                                    GEStep2_kernel,
 
841
                                    1,
 
842
                                                                        NULL,
 
843
                                    globalThreads,
 
844
                                    localThreads,
 
845
                                    0,
 
846
                                    NULL,
 
847
                                    &events[0]);
 
848
    if (status != CL_SUCCESS) { 
 
849
        printf("Error: Enqueueing kernel onto command queue. (clEnqueueNDRangeKernel)\n");
 
850
        return 1;
 
851
    }
 
852
 
 
853
    /* wait for the kernel call to finish execution */
 
854
    status = clWaitForEvents(1, &events[0]);
 
855
    if (status != CL_SUCCESS) { 
 
856
        printf("Error: Waiting for kernel run to finish. (clWaitForEvents)\n");
 
857
        return 1;
 
858
    }
 
859
 
 
860
    status = clReleaseEvent(events[0]);
 
861
    if (status != CL_SUCCESS) { 
 
862
        printf("Error: Release event object. (clReleaseEvent)\n");
 
863
        return 1;
 
864
    }
 
865
 
 
866
    /* Enqueue readBuffer*/
 
867
        //Note: we are reading back from inputBuffer since AI is modified directly in kernel
 
868
    status = clEnqueueReadBuffer(commandQueue,
 
869
                                 inputBuffer,
 
870
                                 CL_TRUE,
 
871
                                 0,
 
872
                                 globalThreads[0] * sizeof(cl_float),
 
873
                                 AI,
 
874
                                 0,
 
875
                                 NULL,
 
876
                                 &events[1]);
 
877
    if (status != CL_SUCCESS) { 
 
878
        printf("Error: clEnqueueReadBuffer failed. (clEnqueueReadBuffer)\n");
 
879
        return 1;
 
880
    }
 
881
    
 
882
    /* Wait for the read buffer to finish execution */
 
883
    status = clWaitForEvents(1, &events[1]);
 
884
    if (status != CL_SUCCESS) { 
 
885
        printf("Error: Waiting for read buffer call to finish. (clWaitForEvents)\n");
 
886
        return 1;
 
887
    }
 
888
 
 
889
    status = clReleaseEvent(events[1]);
 
890
    if (status != CL_SUCCESS) { 
 
891
        printf("Error: Release event object. (clReleaseEvent)\n");
 
892
        return 1;
 
893
    }
 
894
    return 0;
 
895
}
 
896
 
 
897
int run_GEStep3_kernel(cl_float * AI, int i, int n2, int lda2) {
 
898
    cl_int status;
 
899
    cl_event events[2];
 
900
 
 
901
    /* 
 
902
     * The input array to the kernel. This array will eventually be modified
 
903
     * to the inverted array.
 
904
     */
 
905
    status = clSetKernelArg(GEStep3_kernel, 0, sizeof(cl_mem), (void *)&inputBuffer);
 
906
    if (status != CL_SUCCESS) { 
 
907
        printf("Error: Setting kernel argument. (input)\n");
 
908
        return 1;
 
909
    }
 
910
 
 
911
    /*i*/
 
912
    status = clSetKernelArg(GEStep3_kernel, 1, sizeof(int), (void *)&i);
 
913
    if (status != CL_SUCCESS) { 
 
914
        printf("Error: Setting kernel argument. (i)\n");
 
915
        return 1;
 
916
    }
 
917
 
 
918
    /*n2*/
 
919
    status = clSetKernelArg(GEStep3_kernel, 2, sizeof(int), (void *)&n2);
 
920
    if (status != CL_SUCCESS) { 
 
921
        printf("Error: Setting kernel argument. (n2)\n");
 
922
        return 1;
 
923
    }
 
924
 
 
925
        /*lda2*/
 
926
    status = clSetKernelArg(GEStep3_kernel, 3, sizeof(int), (void *)&lda2);
 
927
    if (status != CL_SUCCESS) { 
 
928
        printf("Error: Setting kernel argument. (lda2)\n");
 
929
        return 1;
 
930
    }
 
931
 
 
932
    /* 
 
933
     * Enqueue a kernel run call.
 
934
     */
 
935
    status = clEnqueueNDRangeKernel(commandQueue,
 
936
                                    GEStep3_kernel,
 
937
                                    1,
 
938
                                    NULL,
 
939
                                    globalThreads,
 
940
                                    localThreads,
 
941
                                    0,
 
942
                                    NULL,
 
943
                                    &events[0]);
 
944
    if (status != CL_SUCCESS) { 
 
945
        printf("Error: Enqueueing kernel onto command queue. (clEnqueueNDRangeKernel)\n");
 
946
        return 1;
 
947
    }
 
948
 
 
949
    /* wait for the kernel call to finish execution */
 
950
    status = clWaitForEvents(1, &events[0]);
 
951
    if (status != CL_SUCCESS) { 
 
952
        printf("Error: Waiting for kernel run to finish. (clWaitForEvents)\n");
 
953
        return 1;
 
954
    }
 
955
 
 
956
    status = clReleaseEvent(events[0]);
 
957
    if (status != CL_SUCCESS) { 
 
958
        printf("Error: Release event object. (clReleaseEvent)\n");
 
959
        return 1;
 
960
    }
 
961
 
 
962
    /* Enqueue readBuffer*/
 
963
        //Note: we are reading back from inputBuffer since AI is modified directly in kernel
 
964
    status = clEnqueueReadBuffer(commandQueue,
 
965
                                 inputBuffer,
 
966
                                 CL_TRUE,
 
967
                                 0,
 
968
                                 globalThreads[0] * sizeof(cl_float),
 
969
                                 AI,
 
970
                                 0,
 
971
                                 NULL,
 
972
                                 &events[1]);
 
973
    if (status != CL_SUCCESS) { 
 
974
        printf("Error: clEnqueueReadBuffer failed. (clEnqueueReadBuffer)\n");
 
975
        return 1;
 
976
    }
 
977
 
 
978
    /* Wait for the read buffer to finish execution */
 
979
    status = clWaitForEvents(1, &events[1]);
 
980
    if (status != CL_SUCCESS) { 
 
981
        printf("Error: Waiting for read buffer call to finish. (clWaitForEvents)\n");
 
982
        return 1;
 
983
    }
 
984
 
 
985
    status = clReleaseEvent(events[1]);
 
986
    if(status != CL_SUCCESS) { 
 
987
        printf("Error: Release event object. (clReleaseEvent)\n");
 
988
        return 1;
 
989
    }
 
990
 
 
991
        return 0;
 
992
}
 
993
 
 
994
void invertge(cl_float * AI_d, int lda, int n) {
 
995
    int lda2 = lda * 2;
 
996
    // perform elementary row operations till A in AI becomes identity matrix
 
997
    for (int i = 0; i < n; i++) {
 
998
        // execute kernel
 
999
        run_GEStep1A_kernel(AI_d,i,n*2, lda2);
 
1000
    }
 
1001
 
 
1002
    for (int i = n-1; i >= 0; i--) {
 
1003
        cl_float diag = 1.0;
 
1004
        diag=AI_d[i*lda2+i];
 
1005
        // execute kernels
 
1006
        run_GEStep2_kernel(AI_d,diag,i,n*2, lda2);
 
1007
        run_GEStep3_kernel(AI_d,i,n*2, lda2);
 
1008
    }
 
1009
}
 
1010
 
 
1011
/* inverts nxn matrix input and stores the result in output */
 
1012
void invert(cl_float * input, cl_float *output, int n) {
 
1013
    fprintf(stderr,"starting inversion n = %d ", n);
 
1014
    volatile clock_t gputime;
 
1015
    gputime=clock();
 
1016
 
 
1017
    int lda = ((n+15)&~15|16);
 
1018
    cl_float * AI_d = (cl_float *)malloc(sizeof(cl_float)*n*lda*2);
 
1019
    memset(AI_d,0,sizeof(cl_float)*n*lda*2);
 
1020
    for (int i = 0; i < n; i++) {
 
1021
        memcpy(&AI_d[lda*i*2], &input[n*i], sizeof(cl_float)*n);
 
1022
        AI_d[lda*i*2+n+i] = 1;
 
1023
    }
 
1024
 
 
1025
    cl_int status;
 
1026
 
 
1027
    /////////////////////////////////////////////////////////////////
 
1028
    // Create OpenCL memory buffer
 
1029
    /////////////////////////////////////////////////////////////////
 
1030
    inputBuffer = clCreateBuffer(context,
 
1031
                                 CL_MEM_READ_WRITE | CL_MEM_USE_HOST_PTR,
 
1032
                                 sizeof(cl_float) * globalThreads[0],
 
1033
                                 AI_d,
 
1034
                                 &status);
 
1035
    if (status != CL_SUCCESS) { 
 
1036
        printf("Error: clCreateBuffer (inputBuffer)\n");
 
1037
        exit(0);
 
1038
    }
 
1039
    // Note: there's no output buffer. In kernel, AI_d is modified directly.
 
1040
        // Thus, we should read the result back to host from inputBuffer as well.
 
1041
 
 
1042
    invertge(AI_d, lda, n);
 
1043
    gputime=clock()-gputime;fprintf(stderr, " %7.1f ms ",gputime/1.e3f);
 
1044
    fprintf(stderr, " %7.2f Gflops", 1e-3*(3.0)*n*n*n/3.0/gputime);
 
1045
 
 
1046
#ifdef VERIFY   
 
1047
    // let's verify that
 
1048
    cl_float error=0.0;
 
1049
 
 
1050
    // multiply inverse*xcopy, should be Identity matrix
 
1051
    for (int k = 0; k < n; k++) {
 
1052
        for (int j = 0; j < n; j++) {
 
1053
            cl_float sum = 0;
 
1054
            for (int i = 0; i < n; i++) {
 
1055
                sum += AI[j*lda*2+n+i]*A[i*n+k];
 
1056
                }
 
1057
            if (j!=k) {
 
1058
                error += sum * sum;
 
1059
            } else {
 
1060
                error += (1.0-sum) * (1.0-sum);
 
1061
            }
 
1062
        }
 
1063
    }
 
1064
    fprintf(stderr, " %6.2f SSE", error);
 
1065
#endif  
 
1066
 
 
1067
    //copy the result to output
 
1068
    for (int i = 0; i < n; i++) {
 
1069
        memcpy(&output[n*i], &AI_d[lda*i*2+n], sizeof(cl_float)*n);
 
1070
    }
 
1071
    free(AI_d);
 
1072
    fprintf(stderr," done!\n");
 
1073
}