在python中创建对称矩阵的最快方法，元素如下

Question

我有一个数组集合K0,K1,K2....Kn定义在1D阵列上z

我想在不使用for循环的情况下以最快的方式得到下面的对称矩阵。

[np.trapz(K0*K0,z)  np.trapz(K0*K1,z)  np.trapz(K0*K2,z)  np.trapz(K0*K3,z)...]
[      .            np.trapz(K1*K1,z)  np.trapz(K1*K2,z)  np.trapz(K1*K3,z)...]
A = [      .                    .          np.trapz(K2*K2,z)  np.trapz(K2*K3,z)...]
[      .                    .             .               np.trapz(K3*K3,z)...]
[      .                    .             .                       .           ]

下面是我能处理的最快速度(对于大n来说仍然不够快…n>10000)。

我将这些K存储在称为KK的组合数组中

KK = []
for i in range(n):
KK.append(Ki)
KK = np.array(KK)

A = np.zeros((n,n))
for i in range(n):
A[i,i:] = A[i:,i] = np.trapz((KK[i]*KK[i:]),z)

有什么更快的方法?我不在乎这个解决方案有多不优雅或不像蟒蛇。我只是想加快速度。

Answer 1

你正在使用对称的性质，矩阵使它非常有效。加速的一种方法是使用Numba

import numpy as np
import numba as nb
@nb.njit(cache=True, nogil=True, parallel=True)
def fun(KK,z,n):
A = np.zeros((n,n))
for i in nb.prange(n):
A[i,i:] = A[i:,i] = np.trapz((KK[i]*KK[i:]),z)
return A

老答案

np.trapz(KK.T[:,:,None]@KK.T[:,None,:],z,axis=0) # using matrix multiplication
np.trapz(np.einsum('ik,jk->ijk',KK,KK),z,axis=2) # Using einsum