OpenCL最小线程并行粒度

介绍

由于OpenCL是为各类处理器设备而打造的开发标准的计算语言。因此跟CUDA不太一样的是，其对设备特征查询的项更上层，而没有提供一些更为底层的特征查询。比如，你用OpenCL的设备查询API只能获取最大work group size，但无法获取到最小线程并行粒度。

但是，由于最小线程并行粒度对于OpenCL应用领域最广的GPU而言确实是一个比较重要的参数。如果你的work group的work item的个数是最小线程并行粒度的倍数，那么你的OpenCL kernel程序往往会达到很高的计算效率，同时也能基于这个模型来做一些Memory Bank Confliction的避免措施。因此，我这里提供了一个比较简单的OpenCL kernel来获取当前GPU或其它处理器的最小线程并行粒度。

基本理论

我们知道，一个计算设备由若干个Compute Unit构i成，而一个Compute Unit中包含了多个Processing Element，一个Compute Unit中的所有Processing Element对于一条算术逻辑指令而言是同时进行操作的。而不同的Compute Unit之间也可以是同时进行操作。因此，GPU的并行可以划分为两个层次——一层是Compute Unit内的所有Processing Element的并行操作；另一层是各个Compute Unit的并行操作。

上面是物理层面，如果对于OpenCL逻辑层面，我们可以认为，一个work group的最大work item个数是指一个compute unit最多能调度、分配的线程数。这个数值一般就是一个CU内所包含的PE的个数的倍数。比如，如果一个GPU有2个CU，每个CU含有8个PE，而Max work group size是512，那么说明一个CU至少可以分配供512个线程并发操作所需要的各种资源。由于一个GPU根据一条算术逻辑指令能对所有PE发射若干次作为一个“原子的”发射操作，因此，这一个对程序员而言作为“原子的”发射操作启动了多少个线程，那么我们就可以认为是该GPU的最小并行线程数。如果一款GPU的最小线程并行数是32，那么该GPU将以32个线程作为一组原子的线程组。这意味着，如果遇到分支，那么一组32个线程组中的所有线程都将介入这个分支，对于不满足条件的线程，则会等到这32个线程中其它线程都完成分支处理之后再一起执行下面的指令。

如果我将work group size指定为64，并且在kernel程序里加一个判断，如果pid小于32做操作A，否则做操作B，那么pid为0~31的线程组会执行操作A，而pid为32到63的线程组不会受到阻塞，而会立马执行操作B。此时，两组线程将并发操作（注意，这里是并发，而不是并行。因为上面讲过，GPU一次发射32个线程的话，那么对于多个32线程组将会调度发射指令）。

根据这个特性，我们就可以写一个OpenCL kernel程序来判别当前GPU的最小并行线程粒度。

我们首先会将work group size定为最大能接受的尺寸。然后，我们将这个work group平均划分为两组，对它们进行测试。我们在中间定义了一个local memory的变量，每个线程都能访问它，不过我们只让pid为0以及pid为[max_work_group_size / 2]的线程去访问它，以不受太多干扰。如果这个标志在线程组0执行时被线程组1改变，那么我们就知道这个粒度并非是最小的，然后对前一组再平均划分为2，递归操作。如果在执行线程组0之后标志没有被更改，那么说明这整个线程组是一个原子的线程组，也就是我们所要的最小并行的线程粒度。

在内核程序中，我们还传了一个用于延迟的循环次数，使得非原子的线程组能够被并发执行。

OpenCL代码片断

执行环境为：Windows 7 64bit AMD-APU A6-3420M Visual Studio 2008 ATI Strom APP SDK.

主机端的部分代码片断：

        cl_context context = nullptr;       // OpenCL contextcl_command_queue commandQueue = nullptr;cl_program program = nullptr;       // OpenCL kernel program object that'll be running on the compute devicecl_mem outputMemObj = nullptr;      // output memory object for outputcl_kernel kernel = nullptr;         // kernel objectconst int deviceIndex = 0;context = clCreateContext(NULL,1, &devices[deviceIndex],NULL,NULL,NULL);commandQueue = clCreateCommandQueue(context, devices[deviceIndex], 0, NULL);// Read the kernel code to the bufferFILE *fp = fopen("cl_kernel.cl", "rb");if(fp == nullptr){puts("The kernel file not found!");goto RELEASE_RESOURCES;}fseek(fp, 0, SEEK_END);size_t kernelLength = ftell(fp);fseek(fp, 0, SEEK_SET);char *kernelCodeBuffer = (char*)malloc(kernelLength + 1);fread(kernelCodeBuffer, 1, kernelLength, fp);kernelCodeBuffer[kernelLength] = '\0';fclose(fp);const char *aSource = kernelCodeBuffer;program = clCreateProgramWithSource(context, 1, &aSource, &kernelLength, NULL);status = clBuildProgram(program, 1, &devices[deviceIndex], NULL, NULL, NULL);cl_int outputArg = 0;outputMemObj = clCreateBuffer(context, CL_MEM_WRITE_ONLY, sizeof(outputArg), NULL, NULL);kernel = clCreateKernel(program,"QueryMinimumGranularity", NULL);cl_int inputArg = 1000;status = clSetKernelArg(kernel, 0, sizeof(inputArg), &inputArg);status = clSetKernelArg(kernel, 1, sizeof(outputMemObj), &outputMemObj);size_t groupSize;clGetDeviceInfo(devices[deviceIndex], CL_DEVICE_MAX_WORK_GROUP_SIZE, sizeof(groupSize), &groupSize, NULL);size_t global_work_size[1] = { groupSize };size_t local_work_size[1] = { groupSize };status = clEnqueueNDRangeKernel(commandQueue, kernel, 1, NULL, global_work_size, local_work_size, 0, NULL, NULL);clFinish(commandQueue);     // Force wait until the OpenCL kernel is completedstatus = clEnqueueReadBuffer(commandQueue, outputMemObj, CL_TRUE, 0, sizeof(outputArg), &outputArg, 0, NULL, NULL);char chBuffer[256];wchar_t wsBuffer[256];sprintf(chBuffer, "The minimum granularity is: %d", outputArg);MBString2WCString(wsBuffer, chBuffer, false);MessageBox(hWnd, wsBuffer, L"Notice", MB_OK);

kernel代码：

__kernel void QueryMinimumGranularity(int nLoop, __global int *pOut)
{__local volatile int flag;int index = get_global_id(0);int totalItems = get_global_size(0);do{int halfIndex = totalItems / 2;if(index == 0)flag = 1;barrier(CLK_LOCAL_MEM_FENCE);if(index < halfIndex){for(int i = 0; i < nLoop; i++){if(flag == -1)break;}if(flag != -1){if(index == 0){*pOut = totalItems;flag = 2;}}}else{if(index == halfIndex){if(flag != 2){//while(flag != 1);flag = -1;}}}barrier(CLK_LOCAL_MEM_FENCE);if(flag == 2)break;totalItems /= 2;}while(totalItems > 0);
}

关于Image Engineering & Computer Vision的更多讨论与交流，敬请关注本博客和新浪微博songzi_tea.