自己编写了一个CUDA求解线性方程组的程序,能够正确实现,但是速度很不理想。我用的是联想T430笔记本自带的显卡NVS5400,只比CPU快4倍多一点,考虑到CPU是双核四线程,如果开四个线程粗粒度并行的话,GPU的速度只和CPU大致相当,比较失望。现在把最耗时的函数贴出来,求各位大佬给出一些优化策略来,不胜感激。
typedef double floatype;
//对 Ax=Y线性方程式进行变换,将系数矩阵变成三角矩阵
//pfDeviceA:系数矩阵;pfDeviceY:等式右边的Y向量;nVectorSize:方程组维数;nMainElem:对角线主元;
//nBlockThrdNum:每个块的线程数量;nGridThrdNum:每个Grid的线程数量
global void cudaEliminate( floatype * pfDeviceA, floatype * pfDeviceY, int nVectorSize,
int nMainElem, int nBlockThrdNum, int nGridThrdNum )
{
int nTid = blockIdx.x * nBlockThrdNum + threadIdx.x;
while( nMainElem + nTid < nVectorSize - 1 )
{
int nOriginRowStart = nMainElem * nVectorSize;
int nEliminatedRowStart = (nMainElem+nTid+1)*nVectorSize;
floatype fFactor = pfDeviceA[ nEliminatedRowStart + nMainElem ]/pfDeviceA[ nOriginRowStart + nMainElem ];
for( int i=nMainElem; i<nVectorSize; i++ )
{
pfDeviceA[ nEliminatedRowStart + i ] -= fFactor * pfDeviceA[ nOriginRowStart + i ];
}
pfDeviceY[ nMainElem+nTid+1 ] -= fFactor * pfDeviceY[ nMainElem ];
nTid += nGridThrdNum;
}
}
(1)求上面函数的优化方法。
(2)上面的for循环如果用“动态并行”是不是很适合?速度大概能提高多少?
(3)求一个寻找矩阵中一列的最大值及其所在位置的实现方法,考虑到耗时占比比较少,可以不用归约,能并行实现就行。
(4)为什么我将floatype定义为float速度反而慢了四分之一?
LZ您好:
因为我并不熟悉线性代数相关的数值方法,因此无法告诉您您的方法是否合适而得当,亦无法给出通盘建议,下面仅就局部实现细节简要说下:
1:(本条略过,原因如上)。
2:目测您的for循环仅仅是读取,相减,写入的工作,这个直接用多个线程展开实现即可。
3:只考虑一列的话,你基本上只能用规约了,但是如果要处理很多列的话,你可以一个线程循环解决一列,相邻的一组线程解决相邻的多列,这样如果是行优先存储的矩阵,同时也是合并访问的,效果较好,实现也较为方便。
4:可能和您实现有关,尚不清楚。
大致如上,未及指出,请其他熟悉线性代数数值方法的网友/斑竹/总版主/NV原厂支持予以补充和讨论。
祝您编码顺利~
[
果真如此!我将它展开后再调整了一下block size,比CPU的速度由4倍多提升到了14倍!现将展开后的代码以及调用它的代码粘贴如下:
global void cudaEliminate( floatype * pfDeviceA, floatype * pfDeviceY, int nVectorSize,
int nMainElem, int nBlockThrdNum, int nGridThrdNum )
{
int nTid = blockIdx_x * nBlockThrdNum + threadIdx_x;
while( nTid < (nVectorSize-nMainElem)*(nVectorSize-nMainElem-1) )
{
int nSubBlockWidth = nVectorSize - nMainElem;
int nRowOfThisThrd = nMainElem + nTid/nSubBlockWidth + 1;
int nColOfThisThrd = nMainElem + nTid%nSubBlockWidth;//nTid-(nTid/nSubBlockWidth)*nSubBlockWidth + 1;//
floatype fFactor = pfDeviceA[ nRowOfThisThrd * nVectorSize + nMainElem ]
/ pfDeviceA[ nMainElem * nVectorSize + nMainElem ];
if( 0 == nTid % nSubBlockWidth )//if( nTid == (nTid/nSubBlockWidth)*nSubBlockWidth )//
{
pfDeviceY[ nRowOfThisThrd ] -= fFactor * pfDeviceY[ nMainElem ];
}
else
{
pfDeviceA[ nRowOfThisThrd * nVectorSize + nColOfThisThrd ]
-= fFactor * pfDeviceA[ nMainElem * nVectorSize + nColOfThisThrd ];
}
nTid += nGridThrdNum;
}
}///
调用部分:
bool CudaSolveEquations( floatype * pfDeviceA, floatype * pfDeviceY, int nVectorSize )
{
int i, nBlockThrdNum, nBlockNum;
for( i=0; i<nVectorSize-1; i++ )
{
//nBlockNum = GetBlockNumAndBlockThrdNum( nVectorSize-1, & nBlockThrdNum );
nBlockNum = GetBlockNumAndBlockThrdNum( (nVectorSize-i)(nVectorSize-i-1), & nBlockThrdNum );
cudaEliminate<<<nBlockNum, nBlockThrdNum>>>( pfDeviceA, pfDeviceY, nVectorSize, i, nBlockThrdNum, nBlockNum * nBlockThrdNum );
cudaDeviceSynchronize();
}
…
(1)继续求优化措施。
(2)内存合并访问是否和显卡的位宽有关?比如显卡位宽为64位,如果顺序读入8个双精度数,是否就可以一次读入?
(3)怎样去掉那恼人的C4819警告?
(4)回seumx:你说的方法对我来说就是天书,我目前不准备用混合编程,不过还是要谢谢你了。
合并访问与显存位宽无关,需要相邻线程访问一段连续的(最好是满足对齐要求的)显存。
system
6
接着问问题:
(1)解线性方程组的时候当维数低于1500时工作正常,超过后出现unknown error错误,问一下有经验的版主,一般情况下是哪里出了问题?
(2)维数为1000的双精度,与CPU的答案大概只有前4、5位数字是相同的,单精度的就更低了,正常么?
system
7
LZ您好:
1:这个一般是访存越界了,您可以上一下nsight或者cuda memchecker以便为您快速指出越界位置。
2:这个不好判定,因为浮点数是不满足实数的交换律和结合律的,不同的计算过程(比如CPU上的循环VS GPU上的规约)可能会产生不同的结果,并且这个也不好说哪些是更为准确的,因为CPU结果也只是“一种浮点运算顺序”下的结果,同时,还可能因为不同的中间过程精度,比如CPU上如果使用X87指令,那么会有80bit的精度,到最后才截断,但是如果SIMD指令则没有这么好的中间精度;对于乘加运算,GPU有FMA计算,保留了较好的中间精度,而CPU上一般只能分开计算(最新的CPU才支持FMA的指令集)。
当然如果是您代码实现考虑不周也是有可能的。
大致如此,祝您编码顺利~
system
8
[
汇报一下:收到别的帖子启发,我用每个block处理一行数据,那个乘数因子fFactor申请为共享内存,速度进一步提升了70%-80%,比CPU快了20倍多。
这个核函数就是一个高斯消元过程:一个n*n的方阵,将第二行减去乘了乘数因子的第一行,使得第二行的第一个元素为0;将第三行减去乘了乘数因子的第一行,使得第三行的第一个元素为0;以此类推,最终使得第一列除了第一行其他元素都为0.
由于这个核函数是个主耗时函数,且高斯消元的算法复杂度为N的立方,还是请各位大侠继续提出优化措施,我将感激不尽。将第一行映射为纹理内存是否可行?怎么做?