2

在我的程序代码中,我有numpy值数组和numpy索引数组。两者都是在程序初始化期间预先分配和预定义的。
程序的每一部分都有一个values用于计算的数组和三个索引数组idx_from_exch,idx_valuesidx_to_exch。有一个全局值数组来交换几个部分的值:exch_arr.
大多数时候索引数组有 2 到 5 个索引,很少(很可能永远不会)需要更多索引。dtype=np.int32,shape并且值在整个程序运行期间是恒定的。因此我ndarray.flags.writeable=False在初始化后设置,但这是可选的。索引数组的索引值idx_values并按idx_to_exch数字顺序排序,idx_source 可以排序,但没有办法定义它。与一个值数组/部分对应的所有索引数组都具有相同的shape. 数组通常有 50 到 1000 个元素
。并且在整个程序运行期间是恒定的,数组的值在每次迭代中都会发生变化。 以下是示例数组:valuesexch_arrshapedtype=np.float64

import numpy as np
import numba as nb

values = np.random.rand(100) * 100  # just some random numbers
exch_arr = np.random.rand(60) * 3  # just some random numbers
idx_values = np.array((0, 4, 55, -1), dtype=np.int32)  # sorted but varying steps
idx_to_exch = np.array((7, 8, 9, 10), dtype=np.int32)  # sorted and constant steps!
idx_from_exch = np.array((19, 4, 7, 43), dtype=np.int32)  # not sorted and varying steps

示例索引操作如下所示:

values[idx_values] = exch_arr[idx_from_exch]  # get values from exchange array
values *= 1.1  # some inplace array operations, this is just a dummy for more complex things
exch_arr[idx_to_exch] = values[idx_values]  # pass some values back to exchange array

由于这些操作每次迭代都会应用一次,持续数百万次迭代,因此速度至关重要。在上一个问题中,我一直在研究提高索引速度的许多不同方法,但是考虑到我的应用程序(尤其是通过使用常量索引数组进行索引并将它们传递给另一个索引数组来获取值),我忘记了足够具体。
到目前为止,最好的方法似乎是花哨的索引。我目前也在尝试numba guvectorize,但似乎不值得付出努力,因为我的阵列很小。 memoryviews会很好,但由于索引数组不一定有一致的步骤,我知道没有办法使用memoryviews.

那么有没有更快的方法来进行重复索引?为每个索引操作预定义内存地址数组的某种方式,dtype并且shape始终保持不变?ndarray.__array_interface__给了我一个内存地址,但我无法将它用于索引。我想到了类似的东西:

stride_exch = exch_arr.strides[0]
mem_address = exch_arr.__array_interface__['data'][0]
idx_to_exch = idx_to_exch * stride_exch + mem_address

这可行吗?
我也一直在考虑strides直接使用 with as_strided,但据我所知,只允许一致的步幅,我的问题需要不一致strides

任何帮助表示赞赏!提前致谢!


编辑:
我刚刚纠正了我的示例计算中的一个巨大错误!
该操作values = values * 1.1 改变了数组的内存地址。我在程序代码中的所有操作都布置为不更改数组的内存地址,因为许多其他操作都依赖于使用 memoryviews。因此,我用正确的就地操作替换了虚拟操作:values *= 1.1

4

1 回答 1

0

使用 numpy 布尔数组绕过昂贵的花式索引的一种解决方案是使用 numba 并跳过 numpy 布尔数组中的 False 值。

示例实现:

@numba.guvectorize(['float64[:], float64[:,:], float64[:]'], '(n),(m,n)->(m)', nopython=True, target="cpu")
def test_func(arr1, arr2, inds, res):
    for i in range(arr1.shape[0]):
        if not inds[i]:
            continue
        for j in range(arr2.shape[0]):
            res[j, i] = arr1[i] + arr2[j, i]

当然,使用 numpy 数据类型(较小的字节大小将运行得更快)和目标是"cpu"or "parallel"

于 2022-01-01T17:42:07.213 回答