LINUX.ORG.RU

Любителям OpenCL и вычислений на GPU просьба помочь

 


1

2

Прошу прощения за простыню ниже, что поделать, раз тут не придумали код прятать под кат или в поля со скроллом. Загвоздка в том, что у меня есть небольшой массив четверок сингл флоатов (~ 10 млн.) Мне надо найти в нём сумму всех четверок, а также bounding box. Я не знаю, стоит ли для этого привлекать GPU. У меня зашибенная платформа: CPU AMD Ryzen 1600x и GPU AMD Radeon HD 6950. Вот я и не могу понять, то ли дело в старой видюхе, то ли юзкейс совсем не тот. Поэтому прошу собрать и запустить код ниже. Программа для GPU должна называеться «test.cl» (ну вы видите, там это захардкодено).

Собирается как-то так:

cc -I/usr/local/include -L/usr/local/lib -lOpenCL -o test test.c

Прошу привести вывод того, что напечатает. У меня вот скорость вычислений одного порядка на CPU и GPU, но время передачи данных велико и убивает весь профит. Ещё почему-то, если 2 раза подсчитать одно и то же, то во второй раз на GPU скорость растет на порядок. Только смысл в этом?

С НГ!

struct set_specs {
    float4 sum;
    float4 min;
    float4 max;
};

__kernel void pass1 (__global float4 *input,
                     __global struct set_specs *specs,
                     __local struct set_specs *tmp,
                     unsigned long n)
{
    size_t nitems = get_global_size (0);
    size_t loc_id = get_local_id (0);
    size_t glob_id = get_global_id (0);
    size_t group_id = get_group_id (0);

    float4 sum_item = (0, 0, 0, 0);
    float4 min_item = sum_item;
    float4 max_item = sum_item;

    size_t i, offset, mask;
    for (i = glob_id; i < n; i += nitems) {
        sum_item += input[i];
        min_item = min (min_item, input[i]);
        max_item = max (max_item, input[i]);
    }

    tmp[loc_id].sum = sum_item;
    tmp[loc_id].min = min_item;
    tmp[loc_id].max = max_item;
    barrier (CLK_LOCAL_MEM_FENCE);

    for (offset = 1; offset < get_local_size (0); offset <<= 1) {
        mask = (offset << 1) - 1;
        if ((loc_id & mask) == 0) {
            tmp[loc_id].sum = tmp[loc_id + offset].sum + tmp[loc_id].sum;
            tmp[loc_id].min = min (tmp[loc_id + offset].min, tmp[loc_id].min);
            tmp[loc_id].max = max (tmp[loc_id + offset].max, tmp[loc_id].max);
        }
        barrier(CLK_LOCAL_MEM_FENCE);
    }

    if (loc_id == 0) {
        specs[group_id].sum = tmp[0].sum;
        specs[group_id].min = tmp[0].min;
        specs[group_id].max = tmp[0].max;
    }
}

__kernel void pass2 (__global struct set_specs *input,
                     __global struct set_specs *specs,
                     __local float4 *sum4,
                     __local float4 *min4,
                     __local float4 *max4)
{
    size_t id = get_local_id (0);
    sum4[id] = input[id].sum;
    min4[id] = input[id].min;
    max4[id] = input[id].max;
    barrier(CLK_LOCAL_MEM_FENCE);

    size_t offset, mask;

    for (offset = 1; offset < get_local_size (0); offset <<= 1) {
        mask = (offset << 1) - 1;
        if ((id & mask) == 0) {
            sum4[id] = sum4[id] + sum4[id + offset];
            min4[id] = min (min4[id], min4[id + offset]);
            max4[id] = max (max4[id], max4[id + offset]);
        }
        barrier(CLK_LOCAL_MEM_FENCE);
    }

    if (id == 0) {
        specs->sum = sum4[0];
        specs->min = min4[0];
        specs->max = max4[0];
    }
}
#include <stdlib.h>
#include <stdio.h>
#include <fcntl.h>
#include <unistd.h>
#include <time.h>
#include <sys/types.h>
#include <sys/stat.h>
#include <sys/mman.h>
#include <xmmintrin.h>
#include "CL/cl.h"

struct set_specs {
    cl_float4 sum;
    cl_float4 min;
    cl_float4 max;
};

struct sse_specs {
    float sum[4];
    float min[4];
    float max[4];
};

static struct {
    cl_context context;
    cl_kernel pass1, pass2;
    cl_program program;
    cl_command_queue queue;
    cl_mem output;

    size_t group_size;
    void *program_mapping;
    size_t program_size;
} opencl_context;

static void free_context ()
{
    if (opencl_context.output != NULL) clReleaseMemObject (opencl_context.output);
    if (opencl_context.program != NULL) clReleaseProgram (opencl_context.program);
    if (opencl_context.pass1 != NULL) clReleaseKernel (opencl_context.pass1);
    if (opencl_context.pass2 != NULL) clReleaseKernel (opencl_context.pass2);
    if (opencl_context.queue != NULL) clReleaseCommandQueue(opencl_context.queue);
    if (opencl_context.context != NULL) clReleaseContext(opencl_context.context);
    if (opencl_context.program_mapping != NULL)
        munmap (opencl_context.program_mapping, opencl_context.program_size);
}

static int init ()
{
    cl_context_properties properties[3];
    cl_uint num_of_platforms=0;
    cl_platform_id platform_id;
    cl_device_id device_id;
    cl_uint num_of_devices=0;

    // retreives a list of platforms available
    if (clGetPlatformIDs (1, &platform_id, &num_of_platforms)!= CL_SUCCESS) {
        fprintf(stderr, "Unable to get platform_id\n");
        goto bad;
    }

    // try to get a supported GPU device
    if (clGetDeviceIDs (platform_id, CL_DEVICE_TYPE_GPU, 1, &device_id,
                        &num_of_devices) != CL_SUCCESS) {
        fprintf(stderr, "Unable to get device_id\n");
        goto bad;
    }

    // context properties list - must be terminated with 0
    properties[0]= CL_CONTEXT_PLATFORM;
    properties[1]= (cl_context_properties) platform_id;
    properties[2]= 0;

    int fd = open ("test.cl", O_RDONLY);
    if (fd == -1) {
        fprintf (stderr, "Cannot open OpenCL program\n");
        goto bad;
    }

    struct stat sb;
    fstat (fd, &sb);
    opencl_context.program_size = sb.st_size;
    printf ("Program size %lu bytes\n", opencl_context.program_size);

    opencl_context.program_mapping = mmap (NULL, opencl_context.program_size, PROT_READ, MAP_PRIVATE, fd, 0);
    close (fd);
    
    opencl_context.context = clCreateContext (properties, 1, &device_id, NULL, NULL, NULL);
    if (opencl_context.context == NULL) {
        fprintf (stderr, "Cannot create context\n");
        goto bad;
    }

    opencl_context.queue = clCreateCommandQueue (opencl_context.context, device_id, 0, NULL);
    if (opencl_context.queue == NULL) {
        fprintf (stderr, "Cannot create command queue\n");
        goto bad;
    }

    opencl_context.program = clCreateProgramWithSource (opencl_context.context, 1, (const char **)
                                                        &opencl_context.program_mapping, NULL, NULL);
    if (opencl_context.program == NULL) {
        fprintf (stderr, "Cannot create program\n");
        goto bad;
    }

    if (clBuildProgram (opencl_context.program, 0, NULL, "-cl-fast-relaxed-math -cl-unsafe-math-optimizations",
                        NULL, NULL) != CL_SUCCESS) {
        fprintf(stderr, "Error building program\n");
        char buffer[4096];
        size_t length;
        clGetProgramBuildInfo(opencl_context.program, device_id, CL_PROGRAM_BUILD_LOG,
                              sizeof(buffer), buffer, &length);
        fprintf(stderr, "%s\n",buffer);
        free_context ();
        exit (1);
    }

    // specify which kernel from the program to execute
    opencl_context.pass1 = clCreateKernel (opencl_context.program, "pass1", NULL);
    if (opencl_context.pass1 == NULL) {
        fprintf (stderr, "Cannot create kernel\n");
        goto bad;
    }
    if (clGetKernelWorkGroupInfo (opencl_context.pass1, device_id,
                                  CL_KERNEL_PREFERRED_WORK_GROUP_SIZE_MULTIPLE, sizeof (size_t),
                                  &opencl_context.group_size, NULL) != CL_SUCCESS) {
        fprintf(stderr, "Error getting group size program\n");
        goto bad;
    }
    printf ("Group size = %lu\n", opencl_context.group_size);

    opencl_context.pass2 = clCreateKernel (opencl_context.program, "pass2", NULL);
    if (opencl_context.pass2 == NULL) {
        fprintf (stderr, "Cannot create kernel\n");
        goto bad;
    }

    // create buffers for the input and ouput
    opencl_context.output = clCreateBuffer(opencl_context.context, CL_MEM_READ_WRITE,
                                           sizeof (struct set_specs) * opencl_context.group_size, NULL, NULL);

    return 1;
bad:
    free_context ();
    return 0;
}

static void find_specs (cl_mem input, size_t n, struct set_specs *out)
{
    size_t local_size = opencl_context.group_size;
    size_t global_size = local_size * local_size;
    cl_ulong len_arg;

    len_arg = n;
    clSetKernelArg(opencl_context.pass1, 0, sizeof(cl_mem), &input);
    clSetKernelArg(opencl_context.pass1, 1, sizeof(cl_mem), &opencl_context.output);
    clSetKernelArg(opencl_context.pass1, 2, local_size * sizeof(struct set_specs), NULL);
    clSetKernelArg(opencl_context.pass1, 3, sizeof(cl_ulong), &len_arg);
    clEnqueueNDRangeKernel(opencl_context.queue, opencl_context.pass1, 1, NULL, &global_size,
                           &local_size, 0, NULL, NULL);

    clSetKernelArg(opencl_context.pass2, 0, sizeof(cl_mem), &opencl_context.output);
    clSetKernelArg(opencl_context.pass2, 1, sizeof(cl_mem), &opencl_context.output);
    clSetKernelArg(opencl_context.pass2, 2, sizeof(cl_float4) * local_size, NULL);
    clSetKernelArg(opencl_context.pass2, 3, sizeof(cl_float4) * local_size, NULL);
    clSetKernelArg(opencl_context.pass2, 4, sizeof(cl_float4) * local_size, NULL);
    clEnqueueNDRangeKernel(opencl_context.queue, opencl_context.pass2, 1, NULL, &local_size,
                           &local_size, 0, NULL, NULL);

    clEnqueueReadBuffer(opencl_context.queue, opencl_context.output, CL_TRUE, 0,
                        sizeof(struct set_specs), out, 0, NULL, NULL);

}

#define N 11000000

typedef float vect4[4] __attribute__((aligned (16)));

static void traditional_specs (const vect4 *dots, size_t n, struct sse_specs *specs)
{
    __v4sf sum = _mm_set_ps1 (0);
    __v4sf min = sum;
    __v4sf max = sum;
    __v4sf dot;
    size_t i;

    for (i=0; i<n; i++) {
        dot = _mm_load_ps (dots[i]);
        sum += dot;
        min = _mm_min_ps (min, dot);
        max = _mm_max_ps (max, dot);
    }
    _mm_store_ps (specs->sum, sum);
    _mm_store_ps (specs->min, min);
    _mm_store_ps (specs->max, max);
}

static long gettime()
{
    struct timespec tv;

    clock_gettime (CLOCK_REALTIME, &tv);
    return tv.tv_sec * 1000000000L + tv.tv_nsec;
}

int main ()
{
    if (!init ()) return 1;

    cl_mem input;
    int i;
    long time;

    input = clCreateBuffer(opencl_context.context, CL_MEM_READ_WRITE,
                           sizeof (cl_float4) * N, NULL, NULL);
    cl_float4 *buffer = clEnqueueMapBuffer (opencl_context.queue, input, CL_TRUE,
                                            CL_MAP_WRITE, 0, sizeof (cl_float4) * N,
                                            0, NULL, NULL, NULL);

    time = gettime();
    for (i=0; i<N; i++) {
        buffer[i].x = i;
        buffer[i].y = i;
        buffer[i].z = i;
        buffer[i].w = i;
    }
    time = gettime() - time;
    clEnqueueUnmapMemObject (opencl_context.queue, input, buffer, 0, NULL, NULL);
    printf ("Data loading time = %li\n", time);

    struct set_specs output;

    time = gettime();
    find_specs (input, N, &output);
    time = gettime() - time;
    clReleaseMemObject (input);

    printf ("sum = <%f, %f, %f, %f>, time=%li\n", output.sum.x, output.sum.y,
            output.sum.z, output.sum.w, time);
    printf ("min = <%f, %f, %f, %f>, time=%li\n", output.min.x, output.min.y,
            output.min.z, output.min.w, time);
    printf ("max = <%f, %f, %f, %f>, time=%li\n", output.max.x, output.max.y,
            output.max.z, output.max.w, time);


    vect4 *buffer2 = aligned_alloc (16, sizeof (vect4) * N);
    struct sse_specs output2 __attribute__((aligned (16)));
    for (i=0; i<N; i++) {
        buffer2[i][0] = i;
        buffer2[i][1] = i;
        buffer2[i][2] = i;
        buffer2[i][3] = i;
    }

    time = gettime();
    traditional_specs (buffer2, N, &output2);
    time = gettime() - time;
    printf ("sum = <%f, %f, %f, %f>, time=%li\n", output2.sum[0], output2.sum[1],
            output2.sum[2], output2.sum[3], time);
    printf ("min = <%f, %f, %f, %f>, time=%li\n", output2.min[0], output2.min[1],
            output2.min[2], output2.min[3], time);
    printf ("max = <%f, %f, %f, %f>, time=%li\n", output2.max[0], output2.max[1],
            output2.max[2], output2.max[3], time);
    
    free_context ();
    free (buffer2);
    return 0;

Код толком не читал, но не уверен, что эффективно реализован алгоритм на opencl. В качестве примера, как можно делать reduce, посмотри сорцы pyopencl, как этот алгоритм там реализован. Вообще, вычисления на gpu — это про эффективные операции с памятью, остальное второстепенно.

lu4nik ★★★ ()

Intel(R) Core(TM) i7-3612QM CPU @ 2.10GHz
Intel® HD Graphics 4000

Program size 2334 bytes
Group size = 16
Data loading time = 73320111
sum = <60499974684672.000000, 60499974684672.000000, 60499974684672.000000, 60499974684672.000000>, time=37664331
min = <0.000000, 0.000000, 0.000000, 0.000000>, time=37664331
max = <10999999.000000, 10999999.000000, 10999999.000000, 10999999.000000>, time=37664331
sum = <59260360720384.000000, 59260360720384.000000, 59260360720384.000000, 59260360720384.000000>, time=67240326
min = <0.000000, 0.000000, 0.000000, 0.000000>, time=67240326
max = <10999999.000000, 10999999.000000, 10999999.000000, 10999999.000000>, time=67240326
./test  0,60s user 0,08s system 94% cpu 0,721 total

ox55ff ★★ ()
Ответ на: комментарий от deadplace

Ну хз. У меня видюшка встроенная в проц. И её память находится в RAM. Т.е. в моём случае пересылки - это тупо копирование в пределах RAM. И это если нет оптимизации, чтобы не копировать, а использовать данные по месту. Я вопрос не изучал.

ox55ff ★★ ()
float4 sum_item = (0, 0, 0, 0);

Тут потенциальная бага. Правильно - либо с преобразованием типа и размножением

float4 sum_item = 0;
либо векторный литерал
float4 sum_item = (float4)(0, 0, 0, 0);
Так, как сейчас - это оператор «запятая» из C, в данном случае (со скобками) последний нолик приведётся к float и разиножится до float4. Первые три нолика никуда не идут. Ко всем указателям на входные данные по хорошему надо приписать const и restrict. AMD-шному компилятору можно помочь, задав перед каждым ядром __attribute__((reqd_work_group_size(64, 1, 1))).

У меня AMD R9 380

gcc -O3 test.c -o test -lOpenCL
$ ./test 
Program size 2335 bytes
Group size = 64
Data loading time = 31778963
sum = <60499991461888.000000, 60499991461888.000000, 60499991461888.000000, 60499991461888.000000>, time=190229101
min = <0.000000, 0.000000, 0.000000, 0.000000>, time=190229101
max = <10999999.000000, 10999999.000000, 10999999.000000, 10999999.000000>, time=190229101
sum = <59260360720384.000000, 59260360720384.000000, 59260360720384.000000, 59260360720384.000000>, time=14654692
min = <0.000000, 0.000000, 0.000000, 0.000000>, time=14654692
max = <10999999.000000, 10999999.000000, 10999999.000000, 10999999.000000>, time=14654692
Если данные изначально не на GPU и нет сложных вычислений - нет никакого смысла качать их на GPU, считать что-то там, а потом выкачивать на хост. Быстрее посчитать на процессоре. Так для всех memory-bound задач. (11000000⋅4⋅4/1024/1024/1024)/(14654692e−9)=11.2 ГБ/с на процессоре. Это уже быстрее в разы, чем передать по PCI в одну сторону. Другое дело - если бы был APU с доступом к той же памяти, но тут другой случай.

tim239 ()
Ответ на: комментарий от ox55ff

Ну хз. У меня видюшка встроенная в проц. И её память находится в RAM. Т.е. в моём случае пересылки - это тупо копирование в пределах RAM. И это если нет оптимизации, чтобы не копировать, а использовать данные по месту.

Это не так. Доставание данных из оперативки в кеш ГПУ – очень дорогая операция.

Вообще, судя по «~ 10 млн», 10'000'000 х 16 ÷ 1024 ÷ 1024 == 153 MiB, оно должно сносно работать на современных ГПУ (как у ОП), но код не читал :)

Stil ★★★★★ ()

Ещё почему-то, если 2 раза подсчитать одно и то же, то во второй раз на GPU скорость растет на порядок. Только смысл в этом?

А если написать(версию под cpu) нормально, то её скорость так же вырастит на порядок. Такие дела.

vcerloman ()
Ответ на: комментарий от tim239

Запустил из-под CodeXL. Из 190 мс работы на GPU ядра занимают 2,5 мс - причём второе совсем мелкое. Быстрее CPU. Копирование же 167,8 МБ по PCI занимает 15 мс, т.е уже больше чем посчитать на процессоре. Остальное видимо набегает из-за ленивых инициализаций, если всё уже «прогрето» и такие операции (копирование данных + запуск двух ядер) идут одна за другой то время выполнения одной будет стремиться к 17,5 мс. А вообще второе ядро подозрительное - если размер группы 64 то у нас 11000000/64 = 171875 групп, и при такой двухпроходной редукции должны просматриваться в цикле данные от всех групп, а учитываются только input[0]..input[63].

tim239 ()
Ответ на: комментарий от tim239

Спасибо за столь подробный анализ. Касательно алгоритма, там всё верно. Первое ядро для каждого из 64*64 воркайтема делает редукцию, а потом каждый воркгроуп редуцирует 64 эл-та в 1. В итоге всегда имеем 64 элемента на выходе. Второе ядро делает из 64-х эл-тов один.

anonymous ()
Ответ на: комментарий от anonymous

Ну это вы, анонимусы, такое сосете, я тут почетный 5-звездочный клоун и каку в рот не беру :)

А чего это у тебя бомбануло за чужую системку?

RedEyedMan4 ★★★★★ ()
Последнее исправление: RedEyedMan4 (всего исправлений: 2)

А как дела сейчас обстоят с OpenMP 4x вроде как #pragma omp target map() позволяет делегировать вычисления на gpu и потом забрать от туда данные.

Deleted ()
Ответ на: комментарий от anonymous

Даже так! Хотя мне openmp не нравится из-за своих неродных языку конструкций

Есть ещё openacc для прямого переноса вычислений на gpu или иные ускорители и тоже работает через #pragma.

неродных языку конструкций

Вполне родные иструкции и что самое важное переносимое, если система не поддерживает что-то из то это просто игнорируется и выполняется в одном потоке если поддерживается то получаем все плюшки ускорения пусть даже с некоторым падением производительности, ну, это издержки портируемости и абстракции. Да и вообще как по мне сейчас надо брат за правило хорошего тона касательно С использовать openmp+openacc везде, там где компилятор не умеет или система не может в ускорители директивы просто будут проигнорированы, а спец функции заменены заглушками. Ускорение банальных for за счёт многопоточности и/или простое использование simd без привязки к нему же вполне замечательно будет смотреться. Но всё же самое главное, это то что любой существующий код можно забустить там где это возможно не меняя по сути самого кода, а лишь в некоторых местах помазать #pragma omp /#pragma acc )))

Deleted ()
Ответ на: комментарий от Deleted

Мне GCD больше нравится, которая реализует синхронные и асинхронные очереди. Пофиг, что этой штуки нет под линукс. OpenMP не особо копал, а в GCD очень хорошие механизмы синхронизации.

anonymous ()
Ответ на: комментарий от anonymous

Ну касательно open/acc/mp у меня самого ещё много вопросов, я конечно за них топлю тут, но есть моменты мне не понятные, надо тред поднять на это счёт

Deleted ()