FPGA開発日記

FPGAというより、コンピュータアーキテクチャかもね! カテゴリ別記事インデックス https://sites.google.com/site/fpgadevelopindex/

GPGPU with OpenCL 高速化技法の勉強

GPGPUを手に入れたので、OpenCLの本を買って勉強している。とりあえずこれを買った(3500円もして痛い出費だったが...)

改訂新版 OpenCL入門 1.2対応 マルチコアCPU・GPUのための並列プログラミング

改訂新版 OpenCL入門 1.2対応 マルチコアCPU・GPUのための並列プログラミング

GPGPUの基礎がまったく分かっていないのであれだが、いろんな情報から鑑みるに、どうやらメモリ転送のボトルネックを潰すところが重要になるように思える。 思えば、これ系のアクセラレータというのは世の中に腐るほど存在して、それぞれの使い分けは、どのプログラミングモデルがどのアプリケーションに向いているか、を判断する必要があるように思える。 それはアルゴリズムとしてと適正だけではなく、たとえば、メモリアクセスを頻繁におこすのか、どのようなパタンでアクセスするのか、に寄るような気がしている。

たとえば、AVX/SSEのようなCPUのより近い部分で処理するようなアーキテクチャは、CPUとのやりとりが頻繁に発生するようなアプリケーションで有利だし、FPGAでARMとハードウェアでアクセラレーションするときは、 FPGAの特性を生かすことができるような、より演算単位の細かい細粒度なアプリケーションが向いているように思える(あとはデータパスさえ用意すれば流せるような単純なアプリケーションなど)

そういう意味でGPGPUはメモリアクセス(というかGPU向けにデータを転送する部分)がアプリケーションの設計の上で肝になると思っており、 大量のデータをどのように効率的に流すかを考える必要がある。

という訳で、とりあえず上記の本のチュートリアルも大分すすんだので、試しに配列中の浮動小数点の積和を取っていくプログラムを書いてみた。 以下超適当なアプリケーションで、各エレメントで演算した結果を最後に足し算しているだけなので、正確に全体の配列の積和を取っている訳ではない。

  • 本当にやりたかったこと
answer = 0;
for (i = 0; i < MAX_ARRAY; i++) {
    answer += A[i] * B[i];
}
  • 並列化にあたり考えるのが面倒なため
answer_array[MAX_PE] = {0};
answer = 0;
for (pe = 0; pe < MAX_PE; pe++) {
    base = pe * (MAX_ARRAY/MAX_PE);
    for (i = 0; i < MAX_ARRAY; i++) {
        answer_array[pe] += A[base + i] * B[base + i];
    }
    answer += answer_array[pe];
}
#pragma warning(disable:4996)

#include <stdio.h>
#include <stdlib.h>
#include <time.h>

#ifdef __APPLE__
#include <OpenCL/opencl.h>
#else  // __APPLE__
#include <CL/cl.h>
#endif // __APPLE__

#define MAX_ARRAY (10000)
#define MAX_LOOP  (10000)
#define MAX_SOURCE_SIZE (0x10000)
#define NUM_PE 100

int main()
{
    int i, k, loop;
    float A[MAX_ARRAY], B[MAX_ARRAY], answer_cpu[NUM_PE], answer = 0.0f;
    /* Intializes random number generator */
    time_t t;
    srand((unsigned)time(&t));


    //% making initial random pattern
    for (i = 0; i < MAX_ARRAY; i++) {
        A[i] = (float)rand() / (float)RAND_MAX - 0.5f;
        B[i] = (float)rand() / (float)RAND_MAX - 0.5f;
        // printf("%7.5f %7.5f\n", A[i], B[i]);
    }

    printf("CPU start... ");

    for (loop = 0; loop < MAX_LOOP; loop++) {
        for (i = 0; i < NUM_PE; i++) {
            int base = i * (MAX_ARRAY / NUM_PE);
            answer_cpu[i] = 0.0f;
            for (k = 0; k < MAX_ARRAY / NUM_PE; k++) {
                answer_cpu[i] += A[k + base] * B[k + base];
            }
            // printf("[%d]=%7.5f\n", i, answer_cpu[i]);
            answer += answer_cpu[i];
        }
        // printf("%7.5f\n", answer);
    }

    printf("Answer=%7.5f\n", answer);

    //* GPU program start.

    cl_platform_id platform_id = NULL;
    cl_device_id device_id = NULL;
    cl_context context = NULL;
    cl_command_queue command_queue = NULL;
    cl_mem Amobj = NULL;
    cl_mem Bmobj = NULL;
    cl_mem AnsObj = NULL;
    cl_mem Cmobj = NULL;
    cl_program program = NULL;
    cl_kernel kernel = NULL;
    cl_uint ret_num_devices;
    cl_uint ret_num_platforms;
    cl_int ret;

    FILE *fp;
    const char fileName[] = "./fmaf_performance.cl";
    size_t source_size;
    char *source_str;

    /* loading source code including kernel */
    fp = fopen(fileName, "r");
    if (!fp){
        fprintf(stderr, "Failed to load kernel.\n");
        exit(1);
    }
    source_str = (char *)malloc(MAX_SOURCE_SIZE);
    source_size = fread(source_str, 1, MAX_SOURCE_SIZE, fp);
    fclose(fp);

    /* getting platform device information */
    ret = clGetPlatformIDs(1, &platform_id, &ret_num_platforms);
    ret = clGetDeviceIDs(platform_id, CL_DEVICE_TYPE_DEFAULT, 1, &device_id, &ret_num_devices);

    /* creating OpenCL context */
    context = clCreateContext(NULL, 1, &device_id, NULL, NULL, &ret);

    /* creating command queue */
    command_queue = clCreateCommandQueue(context, device_id, 0, &ret);

    /* Creating buffer object */
    Amobj  = clCreateBuffer(context, CL_MEM_READ_ONLY,  MAX_ARRAY * sizeof (float), NULL, NULL);
    Bmobj  = clCreateBuffer(context, CL_MEM_READ_ONLY,  MAX_ARRAY * sizeof (float), NULL, NULL);
    AnsObj = clCreateBuffer(context, CL_MEM_READ_ONLY,  sizeof(float), NULL, NULL);
    Cmobj  = clCreateBuffer(context, CL_MEM_READ_WRITE, NUM_PE * sizeof(float), NULL, NULL);

    size_t global_item_size = NUM_PE;
    size_t local_item_size = 1;
    int ArraySize = MAX_ARRAY / global_item_size;

    float gpu_answer = 0.0f;

    /* creating kernel program from read source */
    program = clCreateProgramWithSource(context, 1, (const char **)&source_str, (const size_t *)&source_size, &ret);
    ret = clBuildProgram(program, 1, &device_id, NULL, NULL, NULL);

    /* creating kernel with data prallel OpenCL */
    kernel = clCreateKernel(program, "fmaf_performance", &ret);

    printf("GPU start... ");
    for (loop = 0; loop < MAX_LOOP; loop++) {


        /* transfer data into memory buffer */
        ret = clEnqueueWriteBuffer(command_queue, Amobj, CL_TRUE, 0, MAX_ARRAY * sizeof(float), A, 0, NULL, NULL);
        ret = clEnqueueWriteBuffer(command_queue, Bmobj, CL_TRUE, 0, MAX_ARRAY * sizeof(float), B, 0, NULL, NULL);
        ret = clEnqueueWriteBuffer(command_queue, AnsObj, CL_TRUE, 0, sizeof(float), &ArraySize, 0, NULL, NULL);

        /* Setting OpenCL kernel argument */
        ret = clSetKernelArg(kernel, 0, sizeof(cl_mem), (void *)&Amobj);
        ret = clSetKernelArg(kernel, 1, sizeof(cl_mem), (void *)&Bmobj);
        ret = clSetKernelArg(kernel, 2, sizeof(cl_mem), (void *)&AnsObj);
        ret = clSetKernelArg(kernel, 3, sizeof(cl_mem), (void *)&Cmobj);

        /* Execute OpenCL kernel with Data parallel */
        ret = clEnqueueNDRangeKernel(command_queue, kernel, 1, NULL,
            &global_item_size, &local_item_size, 0, NULL, NULL);

        float answer_array[NUM_PE];
        ret = clEnqueueReadBuffer(command_queue, Cmobj, CL_TRUE, 0, sizeof(float) * NUM_PE, answer_array, 0, NULL, NULL);

        for (i = 0; i < NUM_PE; i++){
            gpu_answer += answer_array[i];
        }
    }

    printf("Answer=%7.5f\n", gpu_answer);

    ret = clFlush(command_queue);
    ret = clFinish(command_queue);
    ret = clReleaseKernel(kernel);
    ret = clReleaseProgram(program);
    ret = clReleaseMemObject(Amobj);
    ret = clReleaseMemObject(Bmobj);
    ret = clReleaseMemObject(Cmobj);
    ret = clReleaseCommandQueue(command_queue);
    ret = clReleaseContext(context);
    free(source_str);

    return 0;
}

コードとしてはだいぶヒドイ。しかもGPU版とCPU版で答えあってないし... 最初にCPUで計算し、次にGPUで計算する。面倒なのでMAX_LOOPで同じ値を何度も計算して、足しこんでいる。 こういう楽をしたときに、VisualStudioが最適化してループを消していないかどうか不安だったりする ...(確かめていない)

ぶっちゃけ、この程度だと全然速度は変わらない。むしろオーバヘッドがあってちょっと遅い気がする。 もうすこし複雑な、というか演算インテンシブなアプリケーションを考えたほうがいいのかな (FFTとか)