用 broadcasting/vectorized 解决方案用内部函数调用替换 for 循环
Replacing for loops with function call inside with broadcasting/vectorized solution
问题：
当使用广播而不是广播标量来匹配数组时，矢量化函数会由于某种原因将数组缩小为标量。
MWE：
下面是一个MWE。它包含一个双 for 循环。我在编写不使用 for 循环而是使用 broadcasting/vectorized numpy.
的更快代码时遇到了问题
import numpy as np
def OneD(x, y, z):
ret = np.exp(x)**(y+1) / (z+1)
return ret
def ThreeD(a,b,c):
value = OneD(a[0],b[0], c)
value *= OneD(a[1],b[1], c)
value *= OneD(a[2],b[2], c)
return value
M_1 = M_2 = [[0,0,0],[0,0,1], [1,1,1], [1,0,2]]
scales0 = scales1 = [1.1, 2.2, 3.3, 4.4]
cc0 = cc1 = 1.77
results = np.zeros((4,4))
for s0, n0, in enumerate(M_1):
for s1, n1, in enumerate(M_2):
v = ThreeD(n0, n1, s1)
v *= cc0 * cc1 * scales0[s0] * scales1[s1]
results[s0, s1] += v
虽然我想删除两个 for 循环，但为了简单起见，我首先尝试摆脱内部循环。不过，请随意回答并删除两者。
尝试失败：
这是我改变循环的方式
rr = [0,1,2,3]
myfun = np.vectorize(ThreeD)
for s0, n0, in enumerate(M_1):
#for s1, n1, in enumerate(M_2):
v = myfun(n0, M_2, rr)
v *= cc0 * cc1 * scales0[s0] * scales1[rr]
results[s0, rr] += v
错误信息：
Traceback (most recent call last):
File "main.py", line 36, in <module>
v = myfun(n0, M_2, rr)
File "/usr/lib/python3/dist-packages/numpy/lib/function_base.py", line 1573, in __call__
return self._vectorize_call(func=func, args=vargs)
File "/usr/lib/python3/dist-packages/numpy/lib/function_base.py", line 1633, in _vectorize_call
ufunc, otypes = self._get_ufunc_and_otypes(func=func, args=args)
File "/usr/lib/python3/dist-packages/numpy/lib/function_base.py", line 1597, in _get_ufunc_and_otypes
outputs = func(*inputs)
File "main.py", line 18, in ThreeD
value = OneD(a[0],b[0], c)
IndexError: invalid index to scalar variable.
我还需要向量化 OneD 函数吗？我希望通过矢量化 ThreeD 函数，它可以进行适当的簿记。
在您的循环中，n0 和 n1 是嵌套 M_ 列表的元素，每个列表有 3 个元素。
In [78]: ThreeD(np.arange(3),np.arange(3),3)
Out[78]: 46.577468547527005
OneD 使用数组，所以可以获得完整的 n lists/arrays:
In [79]: OneD(np.arange(3), np.arange(3),3)
Out[79]: array([ 0.25 , 1.84726402, 100.85719837])
In [80]: np.prod(_)
Out[80]: 46.577468547527005
并且产品匹配 ThreeD。
仅查看双循环的 ThreeD 部分：
In [81]: for s0, n0, in enumerate(M_1):
...: for s1, n1, in enumerate(M_2):
...: print(n0,n1,s1, ThreeD(n0, n1, s1))
...:
[0, 0, 0] [0, 0, 0] 0 1.0
[0, 0, 0] [0, 0, 1] 1 0.125
[0, 0, 0] [1, 1, 1] 2 0.037037037037037035
[0, 0, 0] [1, 0, 2] 3 0.015625
[0, 0, 1] [0, 0, 0] 0 2.718281828459045
...
[1, 0, 2] [1, 0, 2] 3 46.577468547527005
从您的列表中创建数组：
In [82]: M1 = np.array(M_1); M2 = np.array(M_2)
In [83]: M1.shape
Out[83]: (4, 3)
我用这个广播调用复制了那些 ThreeD 结果：
In [87]: np.prod(OneD(M1[:,None,:], M2[None,:,:], np.arange(4)[None,:,None]), axis=2)
Out[87]:
array([[1.00000000e+00, 1.25000000e-01, 3.70370370e-02, 1.56250000e-02],
[2.71828183e+00, 9.23632012e-01, 2.73668744e-01, 3.13836514e-01],
[2.00855369e+01, 6.82476875e+00, 1.49418072e+01, 6.30357490e+00],
[2.00855369e+01, 1.85516449e+01, 1.49418072e+01, 4.65774685e+01]])
我正在将 (4,1,3)、(1,4,3) 和 (1,4,1) 数组传递给 OneD。结果是 (4,4,3)，然后我在最后一个轴上相乘得到 (4,4).
剩下的计算是：
In [88]: (cc0*cc1*np.array(scales0)[:,None]*np.array(scales1)[None,:])
Out[88]:
array([[ 3.790809, 7.581618, 11.372427, 15.163236],
[ 7.581618, 15.163236, 22.744854, 30.326472],
[11.372427, 22.744854, 34.117281, 45.489708],
[15.163236, 30.326472, 45.489708, 60.652944]])
In [89]: _87*_88 # multiplying these two 4x4 arrays
Out[89]:
array([[3.79080900e+00, 9.47702250e-01, 4.21201000e-01, 2.36925563e-01],
[2.06089744e+01, 1.40052502e+01, 6.22455564e+00, 9.51755427e+00],
[2.28421302e+02, 1.55228369e+02, 5.09773834e+02, 2.86747781e+02],
[3.04561737e+02, 5.62605939e+02, 6.79698445e+02, 2.82506059e+03]])
匹配`结果：
In [90]: results
Out[90]:
array([[3.79080900e+00, 9.47702250e-01, 4.21201000e-01, 2.36925563e-01],
[2.06089744e+01, 1.40052502e+01, 6.22455564e+00, 9.51755427e+00],
[2.28421302e+02, 1.55228369e+02, 5.09773834e+02, 2.86747781e+02],
[3.04561737e+02, 5.62605939e+02, 6.79698445e+02, 2.82506059e+03]])
问题：
当使用广播而不是广播标量来匹配数组时，矢量化函数会由于某种原因将数组缩小为标量。
MWE：
下面是一个MWE。它包含一个双 for 循环。我在编写不使用 for 循环而是使用 broadcasting/vectorized numpy.
的更快代码时遇到了问题import numpy as np
def OneD(x, y, z):
ret = np.exp(x)**(y+1) / (z+1)
return ret
def ThreeD(a,b,c):
value = OneD(a[0],b[0], c)
value *= OneD(a[1],b[1], c)
value *= OneD(a[2],b[2], c)
return value
M_1 = M_2 = [[0,0,0],[0,0,1], [1,1,1], [1,0,2]]
scales0 = scales1 = [1.1, 2.2, 3.3, 4.4]
cc0 = cc1 = 1.77
results = np.zeros((4,4))
for s0, n0, in enumerate(M_1):
for s1, n1, in enumerate(M_2):
v = ThreeD(n0, n1, s1)
v *= cc0 * cc1 * scales0[s0] * scales1[s1]
results[s0, s1] += v
虽然我想删除两个 for 循环，但为了简单起见，我首先尝试摆脱内部循环。不过，请随意回答并删除两者。
尝试失败：
这是我改变循环的方式
rr = [0,1,2,3]
myfun = np.vectorize(ThreeD)
for s0, n0, in enumerate(M_1):
#for s1, n1, in enumerate(M_2):
v = myfun(n0, M_2, rr)
v *= cc0 * cc1 * scales0[s0] * scales1[rr]
results[s0, rr] += v
错误信息：
Traceback (most recent call last):
File "main.py", line 36, in <module>
v = myfun(n0, M_2, rr)
File "/usr/lib/python3/dist-packages/numpy/lib/function_base.py", line 1573, in __call__
return self._vectorize_call(func=func, args=vargs)
File "/usr/lib/python3/dist-packages/numpy/lib/function_base.py", line 1633, in _vectorize_call
ufunc, otypes = self._get_ufunc_and_otypes(func=func, args=args)
File "/usr/lib/python3/dist-packages/numpy/lib/function_base.py", line 1597, in _get_ufunc_and_otypes
outputs = func(*inputs)
File "main.py", line 18, in ThreeD
value = OneD(a[0],b[0], c)
IndexError: invalid index to scalar variable.
我还需要向量化 OneD 函数吗？我希望通过矢量化 ThreeD 函数，它可以进行适当的簿记。
在您的循环中，n0 和 n1 是嵌套 M_ 列表的元素，每个列表有 3 个元素。
In [78]: ThreeD(np.arange(3),np.arange(3),3)
Out[78]: 46.577468547527005
OneD 使用数组，所以可以获得完整的 n lists/arrays:
In [79]: OneD(np.arange(3), np.arange(3),3)
Out[79]: array([ 0.25 , 1.84726402, 100.85719837])
In [80]: np.prod(_)
Out[80]: 46.577468547527005
并且产品匹配 ThreeD。
仅查看双循环的 ThreeD 部分：
In [81]: for s0, n0, in enumerate(M_1):
...: for s1, n1, in enumerate(M_2):
...: print(n0,n1,s1, ThreeD(n0, n1, s1))
...:
[0, 0, 0] [0, 0, 0] 0 1.0
[0, 0, 0] [0, 0, 1] 1 0.125
[0, 0, 0] [1, 1, 1] 2 0.037037037037037035
[0, 0, 0] [1, 0, 2] 3 0.015625
[0, 0, 1] [0, 0, 0] 0 2.718281828459045
...
[1, 0, 2] [1, 0, 2] 3 46.577468547527005
从您的列表中创建数组：
In [82]: M1 = np.array(M_1); M2 = np.array(M_2)
In [83]: M1.shape
Out[83]: (4, 3)
我用这个广播调用复制了那些 ThreeD 结果：
In [87]: np.prod(OneD(M1[:,None,:], M2[None,:,:], np.arange(4)[None,:,None]), axis=2)
Out[87]:
array([[1.00000000e+00, 1.25000000e-01, 3.70370370e-02, 1.56250000e-02],
[2.71828183e+00, 9.23632012e-01, 2.73668744e-01, 3.13836514e-01],
[2.00855369e+01, 6.82476875e+00, 1.49418072e+01, 6.30357490e+00],
[2.00855369e+01, 1.85516449e+01, 1.49418072e+01, 4.65774685e+01]])
我正在将 (4,1,3)、(1,4,3) 和 (1,4,1) 数组传递给 OneD。结果是 (4,4,3)，然后我在最后一个轴上相乘得到 (4,4).
剩下的计算是：
In [88]: (cc0*cc1*np.array(scales0)[:,None]*np.array(scales1)[None,:])
Out[88]:
array([[ 3.790809, 7.581618, 11.372427, 15.163236],
[ 7.581618, 15.163236, 22.744854, 30.326472],
[11.372427, 22.744854, 34.117281, 45.489708],
[15.163236, 30.326472, 45.489708, 60.652944]])
In [89]: _87*_88 # multiplying these two 4x4 arrays
Out[89]:
array([[3.79080900e+00, 9.47702250e-01, 4.21201000e-01, 2.36925563e-01],
[2.06089744e+01, 1.40052502e+01, 6.22455564e+00, 9.51755427e+00],
[2.28421302e+02, 1.55228369e+02, 5.09773834e+02, 2.86747781e+02],
[3.04561737e+02, 5.62605939e+02, 6.79698445e+02, 2.82506059e+03]])
匹配`结果：
In [90]: results
Out[90]:
array([[3.79080900e+00, 9.47702250e-01, 4.21201000e-01, 2.36925563e-01],
[2.06089744e+01, 1.40052502e+01, 6.22455564e+00, 9.51755427e+00],
[2.28421302e+02, 1.55228369e+02, 5.09773834e+02, 2.86747781e+02],
[3.04561737e+02, 5.62605939e+02, 6.79698445e+02, 2.82506059e+03]])