GPGPU.ORG.CN

根据第二节的介绍，matlab中可以通过mex文件的方式编译C/C++代码，但是对于.cu文件则无能为力。为了解决这个问题，NV的工程师们开发了用于编译cu文件的脚本，下面进行详细介绍。
测试环境：vs2005,matlab 7.6(r2008a),cuda 2.3

首先在NV的网站http://developer.nvidia.com/object/matlab_cuda.html下载matlab的插件包，并解压，共有4个重要文件：

nvmex.m
nvmex_helper.m
nvmexopts.bat
.\bin\nvmex.pl

首先将nvmex.pl拷贝到matlab的安装目录中的bin目录下，例如笔者的系统中就拷贝到C:\MATLAB\bin中。
下面编写cu文件，我们仍以两个矩阵相加为例，具体流程见代码注释。从代码中可看到，同时使用了matlab函数和cuda函数，并在GPU中完成了两个矩阵相加。完成addMatrix.cu的编写后，需要对其进行编译。

首先将nvmex.m,nvmex_helper.m,nvmexopts.bat拷贝到addMatrix.cu文件所在目录，并设置为matlab的当前目录。

1、在matlab命令窗口中运行: nvmex -setup，选择编译器，此处选择vs2005。
2、编辑nvmexopts.bat，修改其中的VSINSTALLDIR选项，指定为vs2005的安装目录，如笔者的系统中：C:\Program Files\Microsoft Visual Studio 8
3、编译命令：
nvmex -f nvmexopts.bat addMatrix.cu -IC:\cuda\include -LC:\cuda\lib –lcudart

编译成功后，在matlab中即可像常规函数一样使用addMatrix函数。

可能会遇到的问题：
1、目录设置：需要仔细按照上述步骤设置好相关目录
2、编译时出现未定义的变量及结构，如_wchar等，则可更新nvmex.pl和nvmexopts.bat来解决。下载链接

/********************************************************************
filename:     addMatrix.cu
file ext:       cu
author:        wy@gpgpu.org.cn
purpose:  test matlab with cuda by nvidia’s nvmex
*********************************************************************/

#include "mex.h"
#include "cuda.h"

/* Kernel*/
__global__ void addMatrix(float* pfC, float* pfA, float* pfB, int M, int N)
{
    int x = blockIdx.x*blockDim.x + threadIdx.x;
    int y = blockIdx.y*blockDim.y + threadIdx.y;
    if (x >= M ||
        y >= N)
    {
        return;
    }
    int index = y*M + x;
    pfC[index] = pfA[index] + pfB[index];
}

void mexFunction( int nlhs, mxArray *plhs[],
                  int nrhs, const mxArray *prhs[])
{
    int i, mA, nA, mB, nB;
    double *A, *B, *C;
    /* 输入变量检查 */
    if (nrhs != 2)
        mexErrMsgTxt("输入参数必须为2个");
    if (nlhs != 1)
        mexErrMsgTxt("输出参数必须为1个");
    mA = mxGetM(prhs[0]);
    nA = mxGetN(prhs[0]);
    mB = mxGetM(prhs[1]);
    nB = mxGetN(prhs[1]);
    if (mA != mB ||
        nA != nB)
        mexErrMsgTxt("输入矩阵尺寸必须相同");

    /*为输出结果分配内存*/
    plhs[0]=mxCreateDoubleMatrix(mA,nA,mxREAL);
    /*获取数据指针*/
    A = mxGetPr(prhs[0]);
    B = mxGetPr(prhs[1]);
    C = mxGetPr(plhs[0]);

    float* d_pfA = NULL;
    float* d_pfB = NULL;
    float* d_pfC = NULL;
    /* 在GPU中分配显存 */
    cudaMalloc( (void **) &d_pfA,sizeof(float)*mA*nA);
    cudaMalloc( (void **) &d_pfB,sizeof(float)*mA*mA);
    cudaMalloc( (void **) &d_pfC,sizeof(float)*mA*mA);

    /*创建一个单精度float缓存区*/
    float* pfSwap = (float *) mxMalloc(sizeof(float)*mA*nA);

    /* 检查输入矩阵A的精度，并假设A, B, C精度相同*/
    mxClassID category = mxGetClassID(prhs[0]);
    if (category == mxSINGLE_CLASS) 
    {
        cudaMemcpy(d_pfA, A, sizeof(float)*mA*nA, cudaMemcpyHostToDevice);
        cudaMemcpy(d_pfB, B, sizeof(float)*mA*nA, cudaMemcpyHostToDevice);
    }
    else if (category == mxDOUBLE_CLASS)
    {
        for (i=0; i < mA*nA; i++)
        {
            pfSwap[i] = (float) A[i];
        }
        cudaMemcpy(d_pfA, pfSwap, sizeof(float)*mA*nA, cudaMemcpyHostToDevice);

        for (i=0; i < mA*nA; i++)
        {
            pfSwap[i] = (float) B[i];
        }
        cudaMemcpy(d_pfB, pfSwap, sizeof(float)*mA*nA, cudaMemcpyHostToDevice);
    }
    /*此处也可以将块的配置作为参数传入*/
    dim3 db,dg;
    db.x = 128;
    db.y = 4;
    dg.x = mA/db.x;
    dg.y = mB/db.y;
    if (mA % db.x > 0)
        dg.x += 1;
    if (nA % db.y > 0)
        dg.y += 1;

     /*调用Kernel函数*/
    addMatrix<<<dg, db>>>(d_pfC, d_pfA, d_pfB, mA, nA);
    /*将计算结果拷贝到临时分配的内存中*/
    cudaMemcpy(pfSwap, d_pfC, sizeof(float)*mA*nA, cudaMemcpyDeviceToHost);
    /*将计算结果赋值给matlab缓存中*/
    for (i=0; i < mA*nA; i++)
    {
        C[i] = (double)pfSwap[i];
    }
    mxFree(pfSwap);
    cudaFree(d_pfA);
    cudaFree(d_pfB);
    cudaFree(d_pfC);
  return;
}