有没有办法在numpy中不使用for来进行复发?
使用np.addwithout关键字可以解决问题dtype="int"
np.add
out
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]
[ 1 1 2 3 5 8 13 21 34 55]
但是,如果dtype更改为np.float
dtype
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.]
[ 1. 1. 2. 1. 3. 1. 4. 1. 5. 1.]
有人可以告诉我为什么吗?或任何其他方式进行递归?
这是使用scipy过滤器的一种相当有效的方法:
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解释器崩溃,因此我们不能使用它,因此我们必须转换为二阶节形式。)
lfilter
过滤方法的问题在于我不确定是否可以将其推广到解决ODE的实际问题。
另一个解决方案是简单地用Python编写循环并使用JIT编译器对其进行加速:
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函数,则效果最佳。
计时结果:
# 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循环则可以编译成很小的循环。