Python numpy 如何优雅的进行矩阵的并行计算

2021-04-19 17:13:00 +08:00
 sunhk25

有一列表对其自身的矩阵求和。 现在是通过双循环来计算的,numpy 有什么优雅的写法吗?

import numpy as np
def how(x, y):
    return x + y
arr = [1, 2, 3, 4, 5, 6]
match_arr = np.zeros((len(arr), len(arr)))
for i in range(0, len(arr)):
    for j in range(0, len(arr)):
        if j >= i:
            X = how(arr[i], arr[j])
            match_arr[i, j] = X
            match_arr[j, i] = X

print(match_arr)
[[ 2.  3.  4.  5.  6.  7.]
 [ 3.  4.  5.  6.  7.  8.]
 [ 4.  5.  6.  7.  8.  9.]
 [ 5.  6.  7.  8.  9. 10.]
 [ 6.  7.  8.  9. 10. 11.]
 [ 7.  8.  9. 10. 11. 12.]]

2583 次点击
所在节点    Python
18 条回复
hsfzxjy
2021-04-19 17:24:35 +08:00
arr = np.array([...])
match_arr = arr[None] + arr[:, None]
geelaw
2021-04-19 17:28:07 +08:00
arr = np.array(arr)
match_arr = np.reshape(arr, (-1, 1)) + np.reshape(arr, (1, -1))

大概是这么个意思,参考 broadcast 的概念。

如果每个元素的计算过程不能用 numpy 内置的运算表达则无法实现,因为 GIL 的存在。
princelai
2021-04-19 18:05:47 +08:00
```python
a,b = np.meshgrid(arr,arr)
match_arr = a+b
```
不一定高效,但是简单
princelai
2021-04-19 18:11:08 +08:00
还有一种方法,原理一样
```python
match_arr = np.mgrid[1:7,1:7].sum(axis=0)
```
sunhk25
2021-04-19 18:17:56 +08:00
@geelaw 那就是说用自定义的 how 函数来循环计算时还是没有优化方法呗
sunhk25
2021-04-19 18:19:49 +08:00
@princelai sum 是我自定义的一个函数
princelai
2021-04-19 18:25:03 +08:00
@sunhk25 #6 你是想说 how 是你自定义的函数?你不是简单的相加是吗?那上 numba,循环放到 numba 里很快,比 numpy 还快。或者你都有两个传播好的 array 了,你改一下 how 函数不就完了
sunhk25
2021-04-19 18:32:56 +08:00
@princelai 对 是这个意思。我研究一下 numba,谢谢
nikan999
2021-04-19 18:49:26 +08:00
先用 numba 如果还想快 就上进程
hsfzxjy
2021-04-19 18:50:36 +08:00
除了 numba,cython 也可以试试。门槛有点高,但是性能优化的上限也高
kickcellardoor
2021-04-19 19:58:50 +08:00
numba,数据量够大甚至可以 PyTorch, GPU 上来并行
Harry1993
2021-04-20 06:22:33 +08:00
樓上說 PyTorch,那我來說 TensorFlow 吧
necomancer
2021-04-20 23:41:49 +08:00
In [4]: from numba import guvectorize, float64, jit
In [5]: @jit(nopython=True)
...: def how(x, y):
...: return x + y
In [6]: @guvectorize([(float64[:], float64[:,:])], '(n)->(n, n)', nopython=True)
...: def f(arr, ret):
...: for i in range(arr.shape[0]):
...: for j in range(arr.shape[0]):
...: if j >= i:
...: tmp = how(arr[i], arr[j])
...: ret[i, j] = tmp
...: ret[j, i] = tmp
In [11]: arr = [np.arange(3), np.arange(10, 13)]

In [12]: f(arr)
Out[12]:
array([[[ 0., 1., 2.],
[ 1., 2., 3.],
[ 2., 3., 4.]],

[[20., 21., 22.],
[21., 22., 23.],
[22., 23., 24.]]])
In [13]: arr = np.arange(3)

In [14]: f(arr)
Out[14]:
array([[0., 1., 2.],
[1., 2., 3.],
[2., 3., 4.]])
necomancer
2021-04-20 23:42:20 +08:00
想要速度一定要用 nopython=True,但是代码得注意一定不能有 object
necomancer
2021-04-20 23:43:18 +08:00
这狗 shit 的排版……
necomancer
2021-04-20 23:47:29 +08:00
另,guvectorize 可以 target='cpu', 'gpu', 'parallel'
necomancer
2021-04-21 00:07:49 +08:00
我测试了一下,有个简单一点的方法,但是会慢一些:
arr = np.arange(3)
def how(x, y):
....if x < y:
........return x + y
....return x * y
np.frompyfunc(how, 2, 1)(arr[:,None], arr)

frompyfunc 会返回一个 ufunc,从而让 numpy 可以 broadcast 自定义的函数。但是效率似乎没有 numba 的 vectorize/guvectorize 高,尤其是 numba 可以 target='gpu'或者'parallel'
sunhk25
2021-04-21 07:20:06 +08:00
@necomancer target 参数学习了

这是一个专为移动设备优化的页面(即为了让你能够在 Google 搜索结果里秒开这个页面),如果你希望参与 V2EX 社区的讨论,你可以继续到 V2EX 上打开本讨论主题的完整版本。

https://www.v2ex.com/t/771716

V2EX 是创意工作者们的社区,是一个分享自己正在做的有趣事物、交流想法,可以遇见新朋友甚至新机会的地方。

V2EX is a community of developers, designers and creative people.

© 2021 V2EX