【问题标题】:Python: How to avoid the loops in building this matrix and to calculate the determinant of it faster?Python:如何避免构建这个矩阵的循环并更快地计算它的行​​列式?
【发布时间】:2026-01-29 17:15:01
【问题描述】:

在 Python 中构建矩阵时遇到了一些问题。 每个元素都有一个循环,each element A_{ij}的形式如图,这里x是q个元素的数组(在下面的代码中用xi表示)。

我已经尝试了以下代码,但是花费了太多时间。我认为这是因为循环的数量,所以我正在考虑将其视为两个矩阵的乘积,但由于 lambda 有两个维度,因此它不起作用。

既然这些代码会以函数的形式出现,而且会被多次使用,有没有什么办法可以让它跑的更快呢?非常感谢!!

def lambdak(i,j,alpha,rho):
    return math.pi * alpha**2 * rho * math.exp(-math.pi**2 * alpha**2 *(i**2 + j**2))
def phik(i,j,x,alpha,rho):
    return cmath.exp(2 * math.pi * 1j * (i*x[0] + j*x[1]))
alpha = 0.5
rho = 50
num = 30
x = np.random.uniform(-0.5,0.5,num)
y = np.random.uniform(-0.5,0.5,num)
xi = np.zeros((num,3))
for i in range(num):
    xi[i] = np.array([x[i], y[i], 0])
q = len(xi)
A = [[np.sum(list(map(lambda j:
                     np.sum(list(map(lambda i:
                                    lambdak(i,j,alpha,rho)/(1-lambdak(i,j,alpha,rho))* phik(i,j,xi[x]-xi[y],alpha,rho),
                                    range(-N,N+1)))),
                     range(-N,N+1)))) for x in range(q)] for y in range(q)]
a = np.linalg.inv(A)

【问题讨论】:

  • 只要看看你的代码,我就可以给你一些建议。 1) 您可以将计算 lambdak(i,j,alpha,rho) 移动到一个单独的函数中,并将其存储在一个二维数组中。您不必为每个 q 重新计算它。 2)此代码也可以并行化,您可以独立计算每个值,但 python 有 GIL 限制。基本上这意味着,即使你使用 python 实现多线程,你也不会看到显着的加速。但是有一些微妙的优化,比如缓存,绝对可以让你的代码更快。
  • @skippy 非常感谢您的回答!我会尝试第一点!
  • @skippy 但是对于您的第二点,您能进一步解释一下吗?我研究过“缓存”,但由于我对它不太熟悉,所以我完全迷失了 XO
  • 没问题。所以我试图用 (2) 得到的是,也许 python 可能不是并行化代码的正确工具,它基本上嵌套了 for 循环,但每次迭代都不依赖于前一个或下一个。你可能想看看 go-lang 或 c++,它们为你提供更好的并发模型,不像 python 有 GIL (wiki.python.org/moin/GlobalInterpreterLock)。此外,您的上述代码可以推送到 GPU 以显着提高速度。 Torch 是一个基于 Lua 的脚本,可以在 GPU (torch.ch) 上运行,但 GPU 可能是答案。

标签: python loops matrix time coding-efficiency


【解决方案1】:

正如您所怀疑的,性能不佳的原因是您的循环。我假设您使用的是numpy,因为您在循环中调用了np.sum。诀窍是把循环翻过来,然后将更大的结构(即矩阵)传递给numpy 函数。

这样做可以显着提高性能。上面给出的代码,可以修改为:

import numpy as np

def lambdak(i,j,alpha,rho):
    return np.pi * alpha**2 * rho * np.exp(-math.pi**2 * alpha**2 *(i**2 + j**2))

def phik(i,j,x,alpha,rho):
    return np.exp(2 * np.pi * 1j * (i*x[:, :, 0] + j*x[:, :, 1]))

alpha = 0.5
rho = 50
num = 30
x = np.random.uniform(-0.5,0.5,num)
y = np.random.uniform(-0.5,0.5,num)
xi = np.zeros((num,3))
for i in range(num):
    xi[i] = np.array([x[i], y[i], 0])

X = np.arange(num).reshape(1,num)
Y = np.arange(num).reshape(num,1)

xi_diff = xi[X] - xi[Y]

N = 30

A = np.sum(map(lambda j:
                np.sum(map(lambda i:
                        lambdak(i,j,alpha,rho)/(1-lambdak(i,j,alpha,rho))* phik(i,j,xi_diff,alpha,rho),
                        range(-N,N+1)), 0),
                     range(-N,N+1)), 0)

a = np.linalg.inv(A)

这里将外部循环转换为矩阵,将其作为一个整体传递给numpy 函数。另外我预先计算了 xi_diff,因为每次调用都会传递整个结构(即使phik() 只使用了一部分)。

这会显着加快速度。但是,计算中的数值稳定性可能会受到影响,当我比较两种方法的输出时,它们的差异约为 0.1%。不过,希望这没问题。

【讨论】:

    最近更新 更多