我正在使用Anaconda分发的Python和Numba,我编写了以下Python函数,它将稀疏矩阵A
(以CSR格式存储)乘以密集向量x
:
@jit
def csrMult( x, Adata, Aindices, Aindptr, Ashape ):
numRowsA = Ashape[0]
Ax = numpy.zeros( numRowsA )
for i in range( numRowsA ):
Ax_i = 0.0
for dataIdx in range( Aindptr[i], Aindptr[i+1] ):
j = Aindices[dataIdx]
Ax_i += Adata[dataIdx] * x[j]
Ax[i] = Ax_i
return Ax
这A
是一个大的scipy
稀疏矩阵,
>>> A.shape
( 56469, 39279 )
# having ~ 142,258,302 nonzero entries (so about 6.4% )
>>> type( A[0,0] )
dtype( 'float32' )
并且x
是一个numpy
数组。以下是调用上述函数的代码片段:
x = numpy.random.randn( A.shape[1] )
Ax = A.dot( x )
AxCheck = csrMult( x, A.data, A.indices, A.indptr, A.shape )
注意@jit
-decorator告诉Numba为该csrMult()
函数进行即时编译。
在我的实验中,我的功能csrMult()
大约是该方法的两倍scipy
.dot()
。对于Numba来说,这是一个令人印象深刻的结果。
然而,仍然MATLAB执行该矩阵-向量乘法约快6倍比csrMult()
。我相信这是因为MATLAB在执行稀疏矩阵向量乘法时使用多线程。
for
使用Numba时如何并行化外环?
Numba曾经拥有一个prange()
函数,这使得并行化并行化并行化for
循环变得简单。不幸的是,Numba不再prange()
[ 实际上,这是错误的,请看下面的编辑 ]。 那么,现在将这个for
-loop 并行化的正确方法是什么,Numba的prange()
功能已经消失了?
当prange()
从Numba中删除时,Numba的开发人员有什么选择?
编辑1: 我更新到Numba的最新版本,即.35,然后
prange()
又回来了!它没有包含在版本.33中,我使用过的版本。 这是个好消息,但不幸的是,当我尝试并行使用for循环时,我收到一条错误消息prange()
。这是来自Numba文档的并行for循环示例(请参见第1.9.2节“显式并行循环”),下面是我的新代码:
from numba import njit, prange
@njit( parallel=True )
def csrMult_numba( x, Adata, Aindices, Aindptr, Ashape):
numRowsA = Ashape[0]
Ax = np.zeros( numRowsA )
for i in prange( numRowsA ):
Ax_i = 0.0
for dataIdx in range( Aindptr[i],Aindptr[i+1] ):
j = Aindices[dataIdx]
Ax_i += Adata[dataIdx] * x[j]
Ax[i] = Ax_i
return Ax
当我使用上面给出的代码片段调用此函数时,收到以下错误:
AttributeError:nopython失败(转换为parfors)'SetItem'对象没有属性'get_targets'
prange
崩溃,我的问题是:并行化这个Python -loop 的正确方法(使用prange
或替代方法)是for
什么?
如下所述,在C ++中并行化一个类似的for循环并获得8x加速,在20 -omp-threads 上运行是微不足道的。必须有一种方法可以使用Numba,因为for循环是令人尴尬的并行(并且因为稀疏矩阵 - 向量乘法是科学计算中的基本操作)。
编辑2: 这是我的C ++版本
csrMult()
。for()
在C ++版本中并行化循环使得代码在我的测试中快了大约8倍。这告诉我,在使用Numba时,Python版本应该可以实现类似的加速。
void csrMult(VectorXd& Ax, VectorXd& x, vector<double>& Adata, vector<int>& Aindices, vector<int>& Aindptr)
{
// This code assumes that the size of Ax is numRowsA.
#pragma omp parallel num_threads(20)
{
#pragma omp for schedule(dynamic,590)
for (int i = 0; i < Ax.size(); i++)
{
double Ax_i = 0.0;
for (int dataIdx = Aindptr[i]; dataIdx < Aindptr[i + 1]; dataIdx++)
{
Ax_i += Adata[dataIdx] * x[Aindices[dataIdx]];
}
Ax[i] = Ax_i;
}
}
}
发布于 2018-09-06 15:05:03
Numba已经更新并且prange()
现在正在使用! (我正在回答我自己的问题。)
在2017年12月12日的博客文章中讨论了对Numba并行计算功能的改进。以下是博客的相关片段:
很久以前(超过20个版本!),Numba过去常常支持一个习惯用来为所谓的循环写并行
prange()
。在2014年对代码库进行重大重构之后,必须删除此功能,但这是自那时以来最常请求的Numba功能之一。在英特尔开发人员对数组表达式进行并行化之后,他们意识到恢复prange
相当容易
使用Numba版本0.36.1,我可以for
使用以下简单代码并行化我令人尴尬的并行循环:
@numba.jit(nopython=True, parallel=True)
def csrMult_parallel(x,Adata,Aindices,Aindptr,Ashape):
numRowsA = Ashape[0]
Ax = np.zeros(numRowsA)
for i in numba.prange(numRowsA):
Ax_i = 0.0
for dataIdx in range(Aindptr[i],Aindptr[i+1]):
j = Aindices[dataIdx]
Ax_i += Adata[dataIdx]*x[j]
Ax[i] = Ax_i
return Ax
在我的实验中,并行化for
-loop使得函数的执行速度比我在问题开始时发布的版本快八倍,该版本已经使用了Numba,但是没有并行化。此外,在我的实验中,并行化版本比Ax = A.dot(x)
使用scipy的稀疏矩阵向量乘法函数的命令快约5倍。 Numba压碎了scipy,我终于有了一个python稀疏矩阵向量乘法程序,它和MATLAB一样快。
发布于 2018-09-06 15:56:26
感谢您的量子更新,丹尼尔。
以下几行可能难以接受,但请相信我,还有更多的事情需要考虑。我曾研究过具有大规模矩阵的HPC /并行计算问题N [TB]; N > 10
及其稀疏的伴随,因此一些经验可能对您的进一步观点有用。
将一段代码并行化的愿望听起来像是一种越来越多的当代重新表达的法术力。问题不在于代码,而在于此类移动的成本。
经济是头号问题。Amdahl定律最初由Gene Amdahl制定,没有考虑过程[PAR]
- 设置+过程[PAR]
- 终结和终止的成本,这些成本确实必须在每个真实世界的实施中支付。
该开销严格阿姆达尔定律描述的这些未避免不良反应的规模,并有助于理解有一个选择只引进并行化之前评估一些新的方面(在这样的一个可接受的成本,因为它是非常,确实很很容易支付超过一个人可能从中获得的收益 - 从降级的处理性能的天真失望是故事的更容易的部分)。
如果愿意更好地理解这个主题并预先计算实际的“ 最小 ”-subProblem-“ 大小 ”,请随时阅读更多关于开销严格的Amdahl定律重新制定的帖子,其中的总和[PAR]
将为至少从现实世界的工具中得到证明,以便将subProblem的并行拆分引入 N_trully_[PAR]_processes
(不是任何“只是” - [CONCURRENT]
但是真实[PARALLEL]
- 这些方式不相等)。
Python是一个伟大的原型生态系统,而且numba
,numpy
和其他编译扩展有助于提高性能,远远超过原生GIL步进python(co - )处理通常提供的。
在这里,你尝试强制numba.jit()
安排工作几乎仅仅通过它的自动-用于自由,jit()
-时间词法分析器(你要把你的代码),它应该既“懂”你的全球目标(什么做的),和还提出了一些矢量化技巧(如何最好地组装一堆CPU指令,以实现这种代码执行的最大效率)。
这听起来很容易,但事实并非如此。
Travis Oliphant的团队在工具方面取得了巨大进步numba
,但是.jit()
当试图转换代码并组装更高效的流程时,让我们切合实际,不要期望在-lexer +代码分析中实现任何形式的自动化操作。机器指令实现高级任务的目标。
@guvectorize
?这里?真的吗?由于[PSPACE]
尺寸调整,你可能会立即忘记要求numba
以某种方式有效地“填充”GPU引擎的数据,其内存占用量远远落后于GPU-GDDR尺寸(根本不是说太“浅”的GPU-这种数学上的内核大小 - “微小”处理只是乘法,可能在[PAR]
,但后来加入[SEQ]
)。
(重新) - 使用数据加载GPU需要大量时间。如果支付了这一点,那么GPU内存延迟对于“微小的”-GPU内核经济性来说也不是很友好 - 你的GPU-SMX代码执行必须支付〜350-700 [ns]
只是为了获取一个数字(很可能不是自动的)重新调整以便在后续步骤中使用最佳的合并SM缓存重用,您可能会注意到,您永远不会让我重复它,永远不要重复使用单个矩阵单元,因此缓存本身不会在350~700 [ns]
每个矩阵单元下提供任何东西),而智能纯粹的numpy
矢量化代码1 [ns]
甚至可以在最大的[PSPACE]
足迹上处理少于每个单元格的矩阵矢量产品。
这是比较的标准。
(分析将更好地显示这里的事实,但这个原则事先是众所周知的,没有测试如何将一些TB
数据移动到GPU结构上只是为了实现这一点。)
给定矩阵的内存规模,A
预期的效果越差,矩阵表示的存储的稀疏组织最有可能破坏大多数(如果不是全部)可能的性能增益,这些性能可通过numba
密集矩阵表示上的矢量化技巧实现。因为高效内存获取高速缓存行重用的可能性几乎为零,并且稀疏性也将破坏任何简单的方法来实现矢量化操作的紧凑映射,并且这些很难保持能够轻松转换为高级CPU硬件矢量处理资源。
Ax = np.zeros_like( A[:,0] )
并将其作为另一个参数传递到numba.jit()
代码的编译部分,以避免重复支付额外的[PTIME,PSPACE]
成本(再次)创建(再次)新的内存分配(如果向量被怀疑为更多则更多在外部协调的迭代优化过程中使用)numba.jit( "f8[:]( f4[:], f4[:,:], ... )" )
-calling接口指令numba.jit()
可用的所有选项及其各自的默认值(可能会将版本更改为版本)(禁用GIL并更好地将目标与numba
+硬件功能保持一致将始终有助于代码的数字密集部分)@jit( signature = [ numba.float32( numba.float32, numba.int32 ), # # [_v41] @decorator with a list of calling-signatures
numba.float64( numba.float64, numba.int64 ) #
], #__________________ a list of signatures for prepared alternative code-paths, to avoid a deferred lazy-compilation if undefined
nopython = False, #__________________ forces the function to be compiled in nopython mode. If not possible, compilation will raise an error.
nogil = False, #__________________ tries to release the global interpreter lock inside the compiled function. The GIL will only be released if Numba can compile the function in nopython mode, otherwise a compilation warning will be printed.
cache = False, #__________________ enables a file-based cache to shorten compilation times when the function was already compiled in a previous invocation. The cache is maintained in the __pycache__ subdirectory of the directory containing the source file.
forceobj = False, #__________________ forces the function to be compiled in object mode. Since object mode is slower than nopython mode, this is mostly useful for testing purposes.
locals = {} #__________________ a mapping of local variable names to Numba Types.
) #____________________# [_v41] ZERO <____ TEST *ALL* CALLED sub-func()-s to @.jit() too >>>>>>>>>>>>>>>>>>>>> [DONE]
def r...(...):
...
https://stackoverflow.com/questions/-100005032
复制相似问题