标签: matrix-inverse

我有这个问题,这需要解决X在AX=B.A是15000 x 15000的顺序,是稀疏和对称的.B是15000 X 7500并且不稀疏.解决X的最快方法是什么？

我可以想到两种方式.

1. 最简单的方式, X = A\B

2. 使用for循环,

   invA = A\speye(size(A))
   for i = 1:size(B,2)
       X(:,i) = invA*B(:,i);
   end
   

   Run Code Online (Sandbox Code Playgroud)

有没有比上面两个更好的方法？如果没有,我提到的两者之间哪一个最好？

我的矩阵形状为（4000，4000），我想求逆。（我对倒置矩阵的直觉分解成这么大的矩阵。）

起始矩阵的大小为e-10，具有以下值：print matrix提供输出

[[  2.19885119e-10   2.16462810e-10   2.13062782e-10 ...,  -2.16462810e-10
   -2.19885119e-10  -2.16462810e-10]
 [  2.16462810e-10   2.19885119e-10   2.16462810e-10 ...,  -2.13062782e-10
   -2.16462810e-10  -2.19885119e-10]
 [  2.13062782e-10   2.16462810e-10   2.19885119e-10 ...,  -2.16462810e-10
   -2.13062782e-10  -2.16462810e-10]
 ..., 
 [ -2.16462810e-10  -2.13062782e-10  -2.16462810e-10 ...,   2.19885119e-10
    2.16462810e-10   2.13062782e-10]
 [ -2.19885119e-10  -2.16462810e-10  -2.13062782e-10 ...,   2.16462810e-10
    2.19885119e-10   2.16462810e-10]
 [ -2.16462810e-10  -2.19885119e-10  -2.16462810e-10 ...,   2.13062782e-10
    2.16462810e-10   2.19885119e-10]]

然后，我使用NumPy的numpy.linalg.inv（）反转矩阵。

import numpy as np
new_matrix = np.linalg.inv(matrix)
print new_matrix

这是我得到的输出：

[[  1.95176541e+25   9.66643852e+23  -1.22660930e+25 ...,  -1.96621184e+25
   -9.41413909e+24   1.33500310e+25]
 [  2.01500967e+25   1.08946558e+24  -1.25813014e+25 ..., …

我有一个对称矩阵.现在,问题在于我需要填充这样的维度矩阵(32**3) x (32**3).我需要填充矩阵的原因是因为在我的程序中我使用它进行各种计算:我正在反转它,我将它与其他矩阵相乘......在我看来,为了执行这些各种你需要实际存储完整矩阵的计算,你不能使用例如只有一半(但我可能是错的,在这种情况下请告诉我应该怎么做).

问题是这样的矩阵对于我的计算机来说太大了,我得到以下错误:

Traceback (most recent call last):
    File "program.py", line 191, in <module>
        A = zeros((n_x*n_y*n_z. n_x*n_y*n_z), float)
MemoryError

在这里,n_x = 32.那么,我该如何解决这个问题呢？有没有办法存储这么大的矩阵,或者是一种避免存储它的聪明方法？如果我可以使用它们而不会在计算中出错,那么这两种方式对我来说都没问题.

为了完整起见,我在下面报告A矩阵的构建方式:

    n_x = n_y = n_z = 32
    L_x = L_y = L_z = n_x
    A = zeros((n_x*n_y*n_z , n_x*n_y*n_z), float)
    P_0 = 50.0
    sigma_x = sigma_y = sigma_z = 0.9
    sigma_root = np.sqrt(sigma_x**2 + sigma_y**2 + sigma_z**2)
    twosigmasquared = 2.*sigma_root**2
    for l in range(n_x*n_y*n_z):
        for m in range(n_x*n_y*n_z):
            A[l][m] …

我有一个对称正定矩阵（例如协方差矩阵），我想计算它的逆矩阵。在数学中，我知道使用 Cholesky 分解来反转矩阵会更有效，特别是当你的矩阵很大时。但我不确定“numpy.lianlg.inv()”是如何工作的。假设我有以下代码：

import numpy as np

X = np.arange(10000).reshape(100,100)
X = X + X.T - np.diag(X.diagonal()) #  symmetry 
X = np.dot(X,X.T) # positive-definite

# simple inversion:
inverse1 = np.linalg.inv(X) 

# Cholesky decomposition inversion:
c = np.linalg.inv(np.linalg.cholesky(X))
inverse2 = np.dot(c.T,c)

哪一个更有效（inverse1 或 inverse2）？如果第二个更有效，为什么 numpy.linalg.inv() 不使用它呢？

从数学上来说，通过 Cholesky 分解对正定矩阵求逆比仅使用 更快np.linalg.inv(X)。然而，当我对两者进行实验时，发现 Cholesky 分解的性能更差！

# Inversion through Cholesky\np = X.shape[0]\nIp = np.eye(p)\n%timeit scipy.linalg.cho_solve(scipy.linalg.cho_factor(X,lower=True), Ip)\n\nThe slowest run took 17.96 times longer than the fastest. This could mean that an intermediate result is being cached.\n10000 loops, best of 3: 107 \xc2\xb5s per loop\n\n\n# Simple inversion\n%timeit np.linalg.inv(X)\n\nThe slowest run took 58.81 times longer than the fastest. This could mean that an intermediate result is being cached.\n10000 loops, best of 3: 25.9 \xc2\xb5s per loop\n

后者花费的时间更短。为什么是这样？在 中R …

下面的矩阵是奇异的，AFAIK 试图反转它应该导致

numpy.linalg.linalg.LinAlgError: Singular matrix

但相反，我确实得到了一些输出矩阵。请注意，输出矩阵是一个无意义的结果，因为它有一行 0（这是不可能的，因为矩阵的逆本身应该是可逆的）！

我是否在这里遗漏了与浮点精度或伪逆计算而不是真逆计算相关的内容？

$ np.__version__ 
'1.13.1'
$ np.linalg.inv(np.array([[2,7,7],[7,7,7],[8,7,7]]))
array([[  0.00000000e+00,   0.00000000e+00,   0.00000000e+00],
       [  3.43131400e+15,  -2.05878840e+16,   1.71565700e+16],
       [ -3.43131400e+15,   2.05878840e+16,  -1.71565700e+16]])```

的Woodbury的矩阵身份指出一些矩阵的秩-K校正的逆可以通过执行秩k值校正到原始矩阵的逆矩阵来计算。

如果A是一个p × p满秩矩阵，通过UCV其中Uis p × k、Cisk × k和Vis 进行秩校正k × p，则伍德伯里恒等式为：

(A + UCV)^{-1} = A^{-1} - A^{-1} U (C^{-1} + V A^{-1} U)^{-1} V A^{-1}

关键是不是反转p × p矩阵，而是反转k × k矩阵。在许多应用中，我们可以假设k < p. A在某些情况下，反转可能很快，例如如果A是对角矩阵。

我在这里实现了这个，假设A是对角线，那C就是身份：

def woodbury(A, U, V, k):
    A_inv = np.diag(1./np.diag(A))  # Fast matrix inversion …

I\xe2\x80\x99m 需要一种有效的方法来反转 python 中的 7000x7000 空气动力学影响系数（密集）矩阵。在使用 FORTRAN 例程之前，我已经开始使用 LAPACK 中的 LU 分解例程来处理问题，我已经看到它在其他相关应用程序中的使用非常有效。不过，我读到 NumPy 和 SciPy 线性系统求解器大多基于对 C 中相同 LAPACK/BLAS 函数的直接调用，并且想知道切换到 FORTRAN 是否真的会在一定程度上减少计算时间这证明放弃一种更简单、更高级的语言是合理的。

如果有 Python 求解器可以保证该大小（1000 到 10000，平方）的矩阵具有相似的性能，那么它们是哪些？

我确实需要矩阵逆，所以切换到迭代 Ax=b 解决方案不是一个选择。

我正在研究机器学习的基础知识，当我阅读有关 Matlab 中正规方程的使用时，我有一个问题。

如果正规方程为：

和

是伪逆，所以：

并pinv(X) 计算 的伪逆X。

为什么在 Matlab/Octave 中，有些人使用：

theta = pinv(X'*X)*X'*y ...(4)

代替

theta = pinv(X)*y ...(5)

？

(4)不是多余的吗？不pinv(X)帮我们计算伪逆公式吗？

谢谢你！

我最近尝试比较不同的 python 和 C++ 矩阵库的线性代数性能，以便了解在即将进行的项目中使用哪些库。虽然线性代数运算有多种类型，但我选择主要关注矩阵求逆，因为它似乎会给出奇怪的结果。我在下面编写了以下代码进行比较，但我想我一定做错了什么。

C++代码

    #include <iostream>
    #include "eigen/Eigen/Dense"
    #include <xtensor/xarray.hpp>
    #include <xtensor/xio.hpp>
    #include <xtensor/xview.hpp>
    #include <xtensor/xrandom.hpp>
    #include <xtensor-blas/xlinalg.hpp> //-lblas -llapack for cblas,   -llapack -L OpenBLAS/OpenBLAS_Install/lib -l:libopenblas.a -pthread for openblas

    //including accurate timer
    #include <chrono>
    //including vector array
    #include <vector>

    void basicMatrixComparisonEigen(std::vector<int> dims, int numrepeats = 1000);
    void basicMatrixComparisonXtensor(std::vector<int> dims, int numrepeats = 1000);
    
    int main()
    {
      std::vector<int> sizings{1, 10, 100, 1000, 10000, 100000};
    
      basicMatrixComparisonEigen(sizings, 2);
      basicMatrixComparisonXtensor(sizings,2);
      return 0;
    }
    
    
    void basicMatrixComparisonEigen(std::vector<int> dims, int numrepeats)
    {
      std::chrono::high_resolution_clock::time_point t1;
      std::chrono::high_resolution_clock::time_point …