其他分享
首页 > 其他分享> > c-Cuda:最小二乘求解,速度较差

c-Cuda:最小二乘求解,速度较差

作者:互联网

最近,我使用Cuda编写了一种称为“正交匹配追踪”的算法.在我丑陋的Cuda代码中,整个迭代需要60秒,而Eigen lib仅需要3秒…

在我的代码中,矩阵A为[640,1024],y为[640,1],在每一步中,我从A中选择一些向量以组成一个新的矩阵A_temp [640,itera],iter = 1:500.我在cpu中新建了一个数组MaxDex_Host []来告诉您选择哪一列.

我想使用最小二乘法从A_temp * x_temp = y获得x_temp [itera,1],我使用了cula API’culaDeviceSgels’和cublas矩阵向量乘法API.

因此,culaDeviceSgels将调用500次,而且我认为这将比Eigen lib的QR.Sovler更快.

我检查了Nisight的性能分析,发现cuestdestory需要很长时间.我在迭代之前先初始化cublas,然后在得到结果后对其进行销毁.所以我想知道custreamdestory和cublasdestory有什么不同?

主要问题是memcpy和函数’gemm_kernel1x1val’.我认为此功能来自“ culaDeviceSgels”

while(itera< 500):我使用cublasSgemv和cublasIsamax来获取MaxDex_Host [itera],然后

        MaxDex_Host[itera]=pos;
    itera++; 
    float* A_temp_cpu=new float[M*itera]; // matrix all in col-major
    for (int j=0;j<itera;j++) // to  get A_temp [M,itera] , the MaxDex_Host[] shows the positon of which column of A to chose , 
    {
        for (int i=0;i<M;i++) //M=640 , and A is 640*1024 ,itera is add 1 each step
        {
            A_temp_cpu[j*M+i]=A[MaxDex_Host[j]*M+i];
        }
    }
          // I must allocate one more array because culaDeviceSgels will decompose the one input Array ,  and I want to use A_temp after least-square solving.
    float* A_temp_gpu;
    float* A_temp2_gpu;  
    cudaMalloc((void**)&A_temp_gpu,Size_float*M*itera);
    cudaMalloc((void**)&A_temp2_gpu,Size_float*M*itera);
    cudaMemcpy(A_temp_gpu,A_temp_cpu,Size_float*M*itera,cudaMemcpyHostToDevice);
    cudaMemcpy(A_temp2_gpu,A_temp_gpu,Size_float*M*itera,cudaMemcpyDeviceToDevice);
    culaDeviceSgels('N',M,itera,1,A_temp_gpu,M,y_Gpu_temp,M);// the x_temp I want is in y_Gpu_temp's return value ,  stored in the y_Gpu_temp[0]——y_Gpu_temp[itera-1]
     float* x_temp;
    cudaMalloc((void**)&x_temp,Size_float*itera);
    cudaMemcpy(x_temp,y_Gpu_temp,Size_float*itera,cudaMemcpyDeviceToDevice);

Cuda的内存管理似乎太复杂了,是否还有其他方便的方法可以解决最小二乘问题?

解决方法:

我认为custreamdestory和gemm_kernel1x1val由您使用的API内部调用,因此与它们没有太多关系.

为了改进您的代码,我建议您执行以下操作.

>您可以通过保留矩阵A的设备副本来摆脱A_temp_cpu.然后,您可以通过内核分配将A的行复制到A_temp_gpu和A_temp2_gpu的行中.这样可以避免执行前两个cudaMemcpys.
>您可以使用itera而不是itera的最大可能值,在while循环外预分配A_temp_gpu和A_temp2_gpu.这将避免循环内的前两个cudaMalloc. x_temp也是如此.
>据我所知,culaDeviceSgels解决了线性方程组.我认为您也可以仅通过使用cuBLAS API来执行相同的操作.例如,您可以先通过cublasDgetrfBatched()执行LU分解,然后使用cublasStrsv()两次来求解两个出现的线性系统.您可能希望查看此解决方案是否导致更快的算法.

标签:cublas,c,cuda,matrix
来源: https://codeday.me/bug/20191009/1879799.html