开勇 OpenHero

OpenCL之我心得

http://bbs.hpctech.com

http://blog.csdn.net/openhero

http://www.comp.hkbu.edu.hk/~kyzhao/

OpenHero 开勇

2009/11/20

或许是为了最求片刻的宁静，选择了宁静的夜，或许在这个时候，只有实验室桌之上几台电脑嗡嗡的蜂鸣声能时刻的提醒我周围的存在。09年算是忙碌的一年，已经快到年底了，反观这一年自己做的事情，对于技术，或许更多的时间都和技术无关吧……追求片刻的宁静，寻求心灵的归属。最近几天黑白颠倒，每天只能睡4-5个小时，在安静的夜，或许能找到自己的真正的位置……

这两天在看OpenCL的一些特性，在反思这个过程，现在计算机技术的发展我想分为两个层面，一个是高层的，云；一个是底层的核心，纵核；在于云这个层面，我们可以采用集群，分布式的一系列的管理软件，任务控制软件来维护整个系统，而对于纵核系统，我们才刚刚起步，OpenCL作为一个纵核的标准，或许才算是一个起步，从OpenCL一出来就开始关注，到现在，OpenCL的特性里面学习了大型集群的思路，也体现了纵核的特点。我比较喜欢OpenCL里面的这个编程模型图：

OpenCL的一些概念我比较看好，例如在高层的多任务的配合，这里避开了在CUDA层面的时候，需要自己写一个多GPU管理的尴尬。

例如Nvidia提供的OpenCL的demo里面，oclMatrixMul的例子，我们可以看到在处理多卡这个环节，OpenCL可以直接利用多卡来分配各自的任务，然后再合起来，这样在一个节点上面，就可以同时利用所有的OpenCL的资源，可以把所有的GPU，CPU都可以利用起来，这样针对一个节点来说，就可以不用程序员人为的做太多的事情，就当成一个节点就可以了。

这个模型也体现了之前我在很多场合说的现在计算机集群的三个层面，

第一个是多机器的层面，我们可以用MPI这样的编程模型来理解这个层面，有多个计算机组成，可能包括同构和异构的集群，这个层面是在多个机器层面来说的。

第二个层面是针对每一个机器节点，我们可以理解为CPU，GPU多个处理器在一个运行系统中的，这个节点可以由系统统一管理。

第三个层面，针对于CPU和GPU本身，也是纵核的载体，这个层面。

我现在把OpenCL归并到第二个层面里面，OpenCL可以起到承上启下的作用。当然，如果OpenCL的底层支持库变换一下，也可以成为第一个层面，或者也可以成为第三个层面。但是现在我只把他归为第二个层面。

其实我们可以想象一个树的结构，一棵树，一棵草，一片树林，一片森林，这样的结构，或许也就是大自然的结构，所以我很喜欢上面的那个图。比较美妙：）ps：分形也十分的美妙……我很喜欢图形……或许是因为小学的时候，学过几年绘画吧……

在看OpenCL的时候，我们不能直接把它和CUDA来做比较，因为两者不是一个层面的东西。往往现在对CUDA的理解也是比较模糊的，用Nvidia自己的话来说，CUDA是一个平台，一个平台而已，它是一套ISA架构，就像X86一样的一套架构，我们现在用只不过是CUDA C，将来也可以是CUDA C++，CUDA java…不应该把CUDA和OpenCL用来做比较，两者是不一样的，OpenCL是一个更高的层面。

………………………………………………分割线

但是现在OpenCL还是在初步阶段，我这里有AMD的平台，有Nvidia的平台，各自都支持OpenCL，但是使用的底层库不一样，效果也不一样。IBM也有自己的OpenCL库，对于统一的支持，OpenCL还有一段路要走，一年，两年，还是三年？不得而知。

至少现在对于OpenCL编程，还是很痛苦的，大家的东西都不一样，而且是运行时编译执行的，这个速度肯定比原生态的要慢很多，在统一这个路上，OpenCL还有很长一段路要走。

……………………………………………………

对于OpenCL的使用，我更喜欢把它作为第二层来使用，作为一个多卡任务的管理者来使用，不过现在同一台机器上，还不能同时支持AMD和Nvidia同时使用OpenCL，所以，还需要等待。虽然我的机器上同时有AMD的卡，也有Nvidia的卡，但是两者都提供了各自的OpenCL，所以这块很麻烦~需要跑各自的程序的时候，还需要把DLL改成对应的才行。

很希望将来OpenCL组委会自己提供统一的OpenCL，让后各自厂商提供各自的运行库比如NVOpenCL.dll AMDOpenCL.dll然后由OpenCL组委会提供的OpenCL.dll来管理这些厂商的DLL，这样对于开发人员来说，就会方便很多。我想，将来他们也会这样做。只是形式不一样而已。

OpenCL的特性今天说了一个，这个算是我对于OpenCL而言，最喜欢的一个特性吧。

以后有时间了，会继续分析OpenCL的一些特性，也喜欢随时欢迎大家一起讨论。

http://www.dearbook.com.cn/book/255585

http://www.china-pub.com/48582&ref=ps

http://product.dangdang.com/product.aspx?product_id=20705176

http://product.dangdang.com/product.aspx?product_id=20705187

http://www.amazon.cn/mn/detailApp/ref=sr_1_1?_encoding=UTF8&s=books&qid=1257088907&asin=B002U848C6&sr=1-1

历时应该有一年多时间，从最开始组织吆喝一群朋友帮忙写书，到中途一些事情的变故，再到后来张舒和艳利同学的不懈努力，终于经历一年多的时间把这本书完成了。

张舒和艳利在书的编写中作出的不懈努力值得我好好学习。张钰博同学已经去了UC 戴维斯 读博士，非常感谢晓炜同学那天能收留我在他们宿舍过夜，也是在那天晚上和张钰博同学聊到了CPU加速，聊到了GPU加速，聊到了CUDA，也才有了最终的这本书。很感谢张钰博同学能带我进入GPU的世界，也是沾了浙大的光……

对于这本书，对晓炜同学深表歉意，对boss来说也深表歉意。

这一年的光影仿佛历历在目，百感交集，有好多话想讲，却不知道从哪里开始。

生活，学习，工作……

继往开来，步步为营，踏踏实实，是对自己的督促，也是给朋友们的勉励。

CUDA让我结识了很多朋友，从全国的巡讲：

广州 （2009年3月28日，星期六） 华南理工大学 A4204/学术报告厅（圆满结束）
成都 （2009年4月11日，星期六） 电子科技大学 清水河基础教学楼B213(圆满结束)
武汉 （2009年4月23日，星期四） 华中科技大学 华中科技大学一号楼报告厅（圆满结束）
武汉 （2009年4月21日，星期二） 武汉大学 武汉大学莱卡厅（圆满结束）
合肥 （2009年4月25日，星期六） 合肥工业大学 南校区学术报告厅一楼大厅（圆满结束）
合肥 （2009年4月27日，星期一） 中国科技大学 西区活动中心学术报告厅（圆满结束）

到中科院的GPU研讨会：

再到Nvidia的GTC全球研讨会；

查看自己Email收发的CUDA相关的问题，超过几千封~回答过不止上千个的问题……

很感谢CSDN能给我提供这样的平台，也很感谢Nvidia的，邓老师，魏鸣，谢强，还有很多帮助过我的朋友，由衷的感谢……

中科院的葛蔚老师，张云泉老师，王龙博士，万宁博士，张先轶同学，陈实同学，刘伟峰...还有很多中科交流的时候认识的老师，同学……

清华的邓仰东老师……

还要感谢浪潮公司能给我和很多公司企业交流的机会,感谢王总，刘总对我的支持，还有team里面的兄弟的支持……；

在GTC遇到的联想的中国和俄罗斯大区的VP，魏总，也是我北理工的师兄，对我的鼓励，还有史博士，真的好感激……

吉星吉达的佟小龙……

全国各地好多的老师，同学，朋友……

还有国外的朋友……

还有找我写作业的国外的同学……让我知道了，不只是有国内的同学有找别人帮忙写作业的情况，其实全世界都有~~只是多少问题……

还有台湾国立中正大学电机所得 陈永维 同学，到香港来旅游的时候，还送我一盒月饼，真的好感谢……

还有参加GTC的时候送书给我的Steve Worley，虽然我知道他看不懂中文……

还有这一年来一直陪在我身边的朋友……

由于CUDA，我认识了好多好多的朋友，非常感谢，是你们的鼓励和支持，才让我感觉到了温暖……谢谢！

/********************************************************************
    created:    2009/09/19
    created:    19:9:2009   12:00
    filename:     vector_matrix_multiplication.cu
    file base:    vector_matrix_multiplication
    file ext:    CUDA
    author:        zhao.kaiyong(at)gmail.com
    purpose:    vector matrix multiplication
    copyright:    everyone can use this code, please specify source.
                任意使用，请注明出处;
http://www.hpctech.com
http://openhero.net
*********************************************************************/

#include <stdio.h>
#include <stdlib.h>
#include <cuda_runtime.h>
#include <cutil.h>

/************************************************************************/
/* Init CUDA                                                            */
/************************************************************************/
#if __DEVICE_EMULATION__

bool InitCUDA(void){return true;}

#else
bool InitCUDA(void)
{
    int count = 0;
    int i = 0;

    cudaGetDeviceCount(&count);
    if(count == 0) {
        fprintf(stderr, "There is no device.\n");
        return false;
    }

    for(i = 0; i < count; i++) {
        cudaDeviceProp prop;
        if(cudaGetDeviceProperties(&prop, i) == cudaSuccess) {
            if(prop.major >= 1) {
                break;
            }
        }
    }
    if(i == count) {
        fprintf(stderr, "There is no device supporting CUDA.\n");
        return false;
    }
    cudaSetDevice(i);

    printf("CUDA initialized.\n");
    return true;
}

#endif
/************************************************************************/
/* Example                                                              */
/************************************************************************/
__global__ static void vector_matrix_mult_kernel(float* A, long wA, float* B, long wB, float* C)
{
    __shared__ float subA[64];

    A = A + threadIdx.x;
    B = B + blockIdx.x * 64 + threadIdx.x;
    C = C + blockIdx.x * 64 + threadIdx.x;

    float subC = 0.0;

    for (int i = 0; i < wA; i+=64)
    {
        subA[threadIdx.x] = A[i];
        __syncthreads();

#pragma unroll
        for (int j = 0; j < 64; j++, B += wB)
        {
            subC += subA[j] * B[0];
        }
        __syncthreads();
    }

    C[0] = subC;
}

//__global__ static void vector_matrix_mult_kernel_32T(float* A, long wA, float* B, long wB, float* C)
//{

//}

#define  RUN_TEST

#define  WA 32
#define  WB 64
/************************************************************************/
/* HelloCUDA                                                            */
/************************************************************************/
int main(int argc, char* argv[])
{

    if(!InitCUDA()) {
        return 0;
    }

    srand(2009);

    long wA = 64* WA;
    long wB = 64* WB;

    long size_A = wA;
    long size_B = wA*wB;
    long size_C = wB;

    float    *hA = (float*)malloc(sizeof(float) * size_A);
    float    *hB = (float*)malloc(sizeof(float) * size_B);
    float    *hC = (float*)malloc(sizeof(float) * size_C);
    float    *testhC = (float*)malloc(sizeof(float) * size_C);

    for (int i = 0; i < size_A; i++)
    {
        hA[i] = (float)rand()/(float)RAND_MAX;
    }

    for (int i = 0; i < size_B; i++)
    {
        hB[i] = (float)rand()/(float)RAND_MAX;
    }

    float    *dA    = 0;
    float    *dB = 0;
    float    *dC = 0;

    CUDA_SAFE_CALL( cudaMalloc((void**) &dA, sizeof(float) * size_A));
    CUDA_SAFE_CALL( cudaMalloc((void**) &dB, sizeof(float) * size_B));
    CUDA_SAFE_CALL( cudaMalloc((void**) &dC, sizeof(float) * size_C));

    CUDA_SAFE_CALL( cudaMemcpy(dA, hA, sizeof(float)*size_A, cudaMemcpyHostToDevice));
    CUDA_SAFE_CALL( cudaMemcpy(dB, hB, sizeof(float)*size_B, cudaMemcpyHostToDevice));

    unsigned int timer = 0;
    CUT_SAFE_CALL( cutCreateTimer( &timer));
    CUT_SAFE_CALL( cutStartTimer( timer));

    dim3 threads = 64;
    dim3 blocks = wB/64;
    vector_matrix_mult_kernel<<<blocks, threads>>>(dA, wA, dB, wB, dC);

    CUT_CHECK_ERROR("Kernel execution failed\n");

    CUDA_SAFE_CALL( cudaMemcpy(hC, dC, sizeof(float) * size_C, cudaMemcpyDeviceToHost));
    CUT_SAFE_CALL( cutStopTimer( timer));
    printf("Processing time: %f (ms)\n", cutGetTimerValue( timer));
    CUT_SAFE_CALL( cutResetTimer( timer));

    for (int i = 0; i < wB; i++)
    {
        float subC = 0.0;
        for (int j = 0; j <wA; j++)
        {
            subC += hA[j] * hB[j*wB + i];
        }
        testhC[i] = subC;
    }

    CUT_SAFE_CALL( cutStopTimer( timer));
    printf("Processing time: %f (ms)\n", cutGetTimerValue( timer));

    CUT_SAFE_CALL( cutDeleteTimer( timer));

#ifdef RUN_TEST

    CUTBoolean res = cutCompareL2fe(testhC, hC, size_C, 1e-6f);
    printf("Test %s \n", (1 == res) ? "PASSED":"FAILED");
#endif

    CUDA_SAFE_CALL( cudaFree(dA));
    CUDA_SAFE_CALL( cudaFree(dB));
    CUDA_SAFE_CALL( cudaFree(dC));
    free(hA);
    free(hB);
    free(hC);
    CUT_EXIT(argc, argv);

    return 0;
}

http://cuda.csdn.net/Contest/pro/nvidia_results.html

历时5个月的NVIDIA CUDA 在8月落下帷幕，大赛举办期间。有300多个学校700多名学生报名。

最终结果

特别贡献奖

奖 别 姓 名 学 校 参赛作品

特别贡献奖 彭江锋 施少怀 漆舒汉 杨植群 华南理工大学/软件学院 基于GPU的多模式网页精确匹配系统

自选赛奖项

奖 别 姓 名 学 校 参赛作品

一等奖 陈实富 中国科学院/深圳先进技术研究院 Glirt:基于CUDA的多模态三维医学图像配准

二等奖 陈伯君 北京大学/工学院 基于CUDA的DSMC高性能并行计算

二等奖 喻勤 电子科技大学 /电子工程学院 基于CUDA的高性能SAR成像模拟

三等奖 胡晓菡 上海交通大学 /生命科学技术学院 CUDA平台下的复杂疾病全基因组基因相互作用计算

三等奖 陆扬 中国科学技术大学 /工程科学学院 使用CUDA模拟顶点、几何着色器-原型演示

命题赛奖项

奖 别 姓 名 学 校 参赛作品

一等奖 翟艳堂 中国科学院 /计算机网络信息中心 稀疏大矩阵与矢量的乘积运算

二等奖 包南森 上海大学/计算机学院 粒子邻居搜索

二等奖 尚书杰 谢迪 浙江大学 /计算机学院 稀疏大矩阵与矢量的乘积运算

三等奖 褚艳利 邵帅 年华 西安电子科技大学 /计算机学院 基于CSR表示的稀疏矩阵-矢量乘法的CUDA实现

三等奖 张广勇 内蒙古大学 /计算机学院 基于GPU 平台的稀疏矩阵与矢量乘积

优秀作品奖

姓 名 学 校 参赛作品

朱小松 大连理工大学 /土木水利学院 基于CUDA的有限差分法求解简单热传递问题

冯涛 沈乐 清华大学 /工程物理系 CUDA 平台下三维反应堆中子通量分布和有效增

尹康学 长安大学 /信息学院 稀疏大矩阵与矢量相乘基于CUDA的一种可行算法

朱德东 复旦大学 /软件学院 基于CUDA的多分辨率图像融合算法

杨植群 华南理工大学 /计算机科学与工程学院 粒子邻居搜索

石丹 范正娟 赵安元 四川大学 /软件学院 CUDA平台下的超声图像斑点噪声抑制

陈晓熹 厦门大学 /信息科学与技术学院 CuParcone

李斌 段元泽 云南大学 /资源环境与地球科学学院 即时计算的动态分形屏保

优秀作品奖得主将会获得由NVIDIA 中国公司颁发的荣誉证书

恭喜以上获奖选手，请保持联系方式畅通，近期会与您联系相关事宜。

10. CUDA cosnstant使用（一）------GPU的革命

序言：最近的事情无比的多，差点就找不到回家的路了，都快忘记出发的起点的时候，冷静下来，侧夜未眠，事情再多，都要一件一件的做好，做不好的，就不接，我想别人也可以接受的，我的个人的能力也是有限的，尽人事，听天命。有的时候，更多的时候是在听天命，而不是在尽人事。不过有时候累得不能再累的时候，真的好想什么都不做，就想回家好好的躺上几天，什么都不用思考，安安静静的就好。回过神来，一件事情，一件事情的做，做好一件，在做另一件，给自己信心，只有自己能把一件一件的事情都步步的做好的时候，就是给自己最大的信心。信心不是别人给的，是自己给自己的。

正文： 书接上回《9.CUDA shared mem使用》讲了shared memory的使用，最近有几个朋友都在问我cosntant的使用的问题，这次首先先讲一下cosntant的使用，下一章节才讲一下cosntant的使用中性能的体现；

下面是一个简单的代码:

/********************************************************************

* cosntant_test.cu

* This is a example of the CUDA program.

* author: zhao.kaiyong(at)gmail.com

*********************************************************************/

#include <stdio.h>

#include <stdlib.h>

#include <cuda_runtime.h>

#include <cutil.h>

/************************************************************************/

/* Init CUDA */

/************************************************************************/

#if __DEVICE_EMULATION__

bool InitCUDA(void){return true;}

#else

bool InitCUDA(void)

{

int count = 0;

int i = 0;

cudaGetDeviceCount(&count);

if(count == 0) {

fprintf(stderr, "There is no device.\n");

return false;

}

for(i = 0; i < count; i++) {

cudaDeviceProp prop;

if(cudaGetDeviceProperties(&prop, i) == cudaSuccess) {

if(prop.major >= 1) {

break;

}

}

}

if(i == count) {

fprintf(stderr, "There is no device supporting CUDA.\n");

return false;

}

cudaSetDevice(i);

printf("CUDA initialized.\n");

return true;

}

#endif

/************************************************************************/

/* Example */

/************************************************************************/

__constant__ char p_HelloCUDA[11];// = "Hello CUDA!";

__constant__ int t_HelloCUDA[11]={0,1,2,3,4,5,6,7,8,9,10};

__constant__ int num = 11;

__global__ static void HelloCUDA(char* result)

{

int i = 0;

for(i = 0; i < num; i++) {

result[i] = p_HelloCUDA[i]+t_HelloCUDA[i];

}

}

/************************************************************************/

/* HelloCUDA */

/************************************************************************/

int main(int argc, char* argv[])

{

if(!InitCUDA()) {

return 0;

}

char helloCUDA[] = "Hdjik CUDA!";

char *device_result = 0;

char host_result[12] ={0};

CUDA_SAFE_CALL( cudaMalloc((void**) &device_result, sizeof(char) * 11));

CUDA_SAFE_CALL( cudaMemcpyToSymbol(p_HelloCUDA, helloCUDA, sizeof(char)*11) );

unsigned int timer = 0;

CUT_SAFE_CALL( cutCreateTimer( &timer));

CUT_SAFE_CALL( cutStartTimer( timer));

HelloCUDA<<<1, 1, 0>>>(device_result);

CUT_CHECK_ERROR("Kernel execution failed\n");

CUDA_SAFE_CALL( cudaThreadSynchronize() );

CUT_SAFE_CALL( cutStopTimer( timer));

printf("Processing time: %f (ms)\n", cutGetTimerValue( timer));

CUT_SAFE_CALL( cutDeleteTimer( timer));

CUDA_SAFE_CALL( cudaMemcpy(&host_result, device_result, sizeof(char) * 11, cudaMemcpyDeviceToHost));

printf("%s\n", host_result);

CUDA_SAFE_CALL( cudaFree(device_result));

CUT_EXIT(argc, argv);

return 0;

}

这里写了两种使用cosntant的方法：

1. 一种方法是直接在定义的时候初始化constant：

__constant__ int t_HelloCUDA[11]={0,1,2,3,4,5,6,7,8,9,10};

__constant__ int num = 11;

在kernel里面直接使用就可以了；

2. 第二种方法是定义一个cosntant数组，然后使用函数初始化它；

__constant__ char p_HelloCUDA[11];// = "Hello CUDA!";

CUDA_SAFE_CALL( cudaMemcpyToSymbol(p_HelloCUDA, helloCUDA, sizeof(char)*11) );

前面一句是定义constant，后面一句是初始化这个常量，在kernel里面使用的时候，就按照定义的方式使用就可以了；

当然也有朋友想定义自己的结构体，当然是可以的，只是在初始化的时候，copy相应的结构体就可以了。这个只是一个初步的使用方法，希望对家有用。

下一次会对cosntant的使用中体现出来的性能优势做一个简单的分析。