numpy复发 [英] Recurrence with numpy

问题描述

是否有任何方法可以在不使用numpy中使用for的情况下进行复发?

Is there any way to make recurrence without using for's in numpy?

使用np.add和out关键字可以解决dtype="int"

Using np.add with out keyword do the trick with dtype="int"

import numpy as np
N = 100

fib = np.zeros(N, dtype=np.int)
fib[:2] = 1.
np.add(fib[:-2], fib[1:-1], out=fib[2:])

print(fib[:10])

输出:[ 1 1 2 3 5 8 13 21 34 55]

但是，如果dtype更改为np.float

However, if dtype is changed to np.float

import numpy as np
N = 100

fib = np.zeros(N, dtype=np.float)
fib[:2] = 1.
np.add(fib[:-2], fib[1:-1], out=fib[2:])

print(fib[:10])

输出:[ 1. 1. 2. 1. 3. 1. 4. 1. 5. 1.]

有人可以告诉我为什么吗?或其他任何递归方法?

Could someone tell me why? or any other way to make recursion?

推荐答案

以下是使用scipy过滤器的相当有效的方法:

Here is a reasonably efficient method using a scipy filter:

import numpy as np
from scipy import signal

def scipy_fib(n):
    x = np.zeros(n, dtype=np.uint8)
    x[0] = 1

    sos = signal.tf2sos([1], [1, -1, -1])
    return signal.sosfilt(sos, x)

(请注意，我们不能使用lfilter，因为就信号处理而言，这是一个不稳定的过滤器，它会使Python解释器崩溃，因此我们必须转换为二阶节形式.)

(Note that we can't use lfilter because in terms of signal processing this is an unstable filter which makes the Python interpreter crash, so we have to convert to second-order sections form.)

过滤方法的问题在于我不确定是否可以将其推广到解决ODE的实际问题中.

The problem with the filtering approach is that I'm not sure if it can be generalized to the actual problem of solving an ODE.

另一种解决方案是简单地用Python编写循环并使用 JIT编译器进行加速:

Another solution is to simply write the loop in Python and accelerate it with a JIT compiler:

import numba as nb    

@nb.jit(nopython=True)
def numba_fib(n):
    y = np.empty(n)
    y[:2] = 1

    for i in range(2, n):
        y[i] = y[i-1] + y[i-2]

    return y

JIT方法的优点是您可以实现各种花哨的东西，但是如果坚持使用简单的循环和分支而不调用任何(非JITted)Python函数，则效果最佳.

The advantage of the JIT approach is that you can implement all kinds of fancy stuff, but it works best if you stick to simple loops and branches without calling any (non-JITted) Python functions.

计时结果:

# Baseline without JIT:
%timeit numba_fib(10000)  # 100 loops, best of 3: 5.47 ms per loop

%timeit scipy_fib(10000)  # 1000 loops, best of 3: 719 µs per loop
%timeit numba_fib(10000)  # 10000 loops, best of 3: 33.8 µs per loop

# For fun, this is how Paul Panzer's solve_banded approach compares:
%timeit banded_fib(10000)  # 1000 loops, best of 3: 1.33 ms per loop

scipy过滤器方法比纯Python循环快，但比JITted循环慢.我猜想filter函数是相对通用的，它可以完成我们不需要的东西，而JITted循环可以编译成很小的循环.

The scipy filter approach is faster than the pure Python loop but slower than the JITted loop. I guess the filter function is relatively generic and does stuff we don't need here, while the JITted loop compiles down to a very small loop.

这篇关于numpy复发的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！