3

我正在尝试有效地解决 Python 3.6 中的以下线性系统：

b = A x

其中 A 是 N x N Toeplitz 矩阵（即，从左到右的对角线是常数），x, b 是 N x N 矩阵。在 Scipy 中实现为 scipy.linalg.solve_toeplitz 的 Durbin-Levinson 应该利用 A 的结构有效地解决上述系统，并且当 x,b 是向量时确实可以做到这一点。

但是，当 x, b 是矩阵时，它比 scipy.linalg.solve 和 numpy.linalg.solve 慢得多（请参阅我的测试脚本的输出）

执行时间在我的应用程序中至关重要。我有哪些选择可以加快流程并利用 Toeplitz 结构？

scipy.linalg.solve_toeplitz 速度慢的一个可能解释是，它对矩阵输入使用慢速 for 循环来求解 x、b 为向量的单个线性系统（一次求解每一列）。

"""This script serves to benchmark several methods for solving a
linear equation involving a Toeplitz matrix and determine which one is
faster.

We consider the equation: b = A x, where A is Toeplitz, b and x are matrices and we wish to
solve for x.
"""
import numpy as np
import numpy.linalg
import scipy.linalg
import time

N = 500
np.random.seed(1)

# Construct random vectors/matrices
x = np.random.rand(N, N)
c,r = np.random.rand(2, N)
A = scipy.linalg.toeplitz(c, r)
b = A.dot(x)

# Validate various solutions
x_sol = []
x_sol.append(('numpy.linalg.solve(A, b)', numpy.linalg.solve(A, b)))
x_sol.append(('scipy.linalg.solve(A, b)', scipy.linalg.solve(A, b)))
x_sol.append(('scipy.linalg.solve_toeplitz((c, r), b)', scipy.linalg.solve_toeplitz((c, r), b)))
for solution in x_sol:
    print("Method: {} - error: {}".format(solution[0], numpy.linalg.norm(solution[1] - x)))

# Time the solutions
x_time = []
for solution in x_sol:
    t_start = time.time()
    eval(solution[0])
    t_end = time.time()
    x_time.append(t_end - t_start)
print("Timings:")
print(x_time)

脚本输出：

Method: numpy.linalg.solve(A, b) - error: 4.8794411236474704e-11
Method: scipy.linalg.solve(A, b) - error: 4.824494916488265e-11
Method: scipy.linalg.solve_toeplitz((c, r), b) - error: 2.7607739766053664e-08
Timings:
[0.08000469207763672, 0.03300189971923828, 0.6740384101867676]

文档：https ://docs.scipy.org/doc/scipy-0.17.0/reference/generated/scipy.linalg.solve_toeplitz.html

4

0 回答 0