2 回答
TA贡献1842条经验 获得超12个赞
您的代码有几个问题。
触发您观察到的错误的那个是从返回值scipy.linalg.bicgstab()atuple而不是您期望的 NumPy 数组派生的。
另一个问题是您尝试访问(nx, ny)索引i, j范围分别为0到ny和的形状对象nx。因此,除非您有nx == ny上面的代码,否则在某些时候会超出数组边界。
最后,所有这些都是通过显式循环完成的。然而,NumPy 提供了更好的工具来获得你想要的东西,特别是np.reshape().
例如:
import numpy as np
nx = 800
ny = 1200
y = np.random.randint(0, 100, nx * ny)
def reshape_loops(y):
solution = np.zeros((nx, ny))
for i in range(0, nx):
for j in range(0, ny):
solution[i, j] = y[i + nx * j]
return solution
def reshape_np(y):
return np.reshape(y.copy(), (nx, ny), order='F')
print(np.allclose(reshape_loops(y), reshape_np(y)))
# True
%timeit reshape_loops(y)
# 1 loop, best of 3: 319 ms per loop
%timeit reshape_np(y)
# 1000 loops, best of 3: 1.25 ms per loop
矢量化方法快了约 250 倍。
TA贡献1828条经验 获得超4个赞
我在官方文档中四处寻找。事实证明,所有稀疏迭代求解器都会返回两件事:解和收敛信息。如果直接写成y = sp.linalg.bicgstab(A, b)
,y 就变成了一个形状为 (2,) 的元组,其中第一个元素是解,第二个元素是收敛信息。我通过这样做来修复它y, exit_code = sp.linalg.bicgstab(A, b)
。现在它工作正常
添加回答
举报