首页 > 解决方案 > 构建一个numpy数组作为前一个元素的函数

问题描述

我想创建一个 numpy 数组,其中第一个元素是一个已定义的常量,并且每个下一个元素都通过以下方式定义为前一个元素的函数:

import numpy as np

def build_array_recursively(length, V_0, function):
    returnList = np.empty(length)
    returnList[0] = V_0
    for i in range(1,length):
        returnList[i] = function(returnList[i-1])
    return returnList

d_t = 0.05
print(build_array_recursively(20, 0.3, lambda x: x-x*d_t+x*x/2*d_t*d_t-x*x*x/6*d_t*d_t*d_t))

上面的打印方法输出

[0.3 0.28511194 0.27095747 0.25750095 0.24470843 0.23254756 0.22098752
 0.20999896 0.19955394 0.18962586 0.18018937 0.17122037 0.16269589
 0.15459409 0.14689418 0.13957638 0.13262186 0.1260127 0.11973187 0.11376316]

在没有 for 循环的情况下,有没有一种快速的方法在 numpy 中执行此操作?如果是这样,有没有办法在当前元素之前处理两个元素,例如可以类似地构造斐波那契数组吗?我在这里发现了一个类似的问题

是否可以对 NumPy 数组的递归计算进行矢量化,其中每个元素都依赖于前一个元素?

但一般没有回答。在我的示例中,差分方程很难手动求解。

标签: pythonnumpyscipy

解决方案


使用函数使代码更清晰(而不是内联 lambda):

def fn(x):
    return x-x*d_t+x*x/2*d_t*d_t-x*x*x/6*d_t*d_t*d_t

还有一个结合build_array_recursively和元素的函数method2

def foo1(length, V_0, function):
    returnList = np.empty(length)
    returnList[0] = x = V_0
    for i in range(1,length):
        returnList[i] = x = function(x)
    return returnList

In [887]: timeit build_array_recursively(20,0.3, fn);                                                                     
61.4 µs ± 63 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)
In [888]: timeit method2(20,0.3, fn);                                                                  
16.9 µs ± 103 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
In [889]: timeit foo1(20,0.3, fn);                                                                     
13 µs ± 29.7 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)

和 中的主要节省时间的方法是将最后一个值从一次迭代转移到下一次迭代,而不是使用 进行method2索引。foo2xreturnList[i-1]

分配给预先分配的数组或列表追加的累积方法不太重要。性能通常是相似的。

这里的计算很简单,你在循环中所做的细节会对整个时间产生很大的影响。

所有这些都是循环。有些ufunc有一个reduce(and accumulate) 方法,可以将函数重复应用于输入数组的元素。 np.sum, np.cumsum, 等利用这个。但是你不能用一个通用的 Python 函数来做到这一点。

您必须使用某种编译工具numba来更快地执行这种循环。


推荐阅读