【问题标题】:What is the best algorithm to find a determinant of a matrix?找到矩阵行列式的最佳算法是什么?
【发布时间】:2013-05-21 01:24:25
【问题描述】:

谁能告诉我哪个是找到大小为N x N的矩阵的行列式值的最佳算法?

【问题讨论】:

  • 除了大小之外,我们对矩阵了解更多吗?稀疏吗?
  • 尽管标记了stackoverflow.com/questions/1886280/… 的答案与语言无关,所以我建议这是重复的。
  • 矩阵算法足够复杂,所以你不应该自己实现它们;使用完善的库,如 LAPACK。编写库的人已经选择了行列式的最佳实现(可能是密集矩阵的 LU 分解)。
  • numpy 使用什么算法?

标签: algorithm language-agnostic determinants


【解决方案1】:

行缩减

找到 nxn 矩阵的行列式的最简单方法(实际上也不错)是通过行归约。通过记住一些关于行列式的简单规则,我们可以解决以下形式:

det(A) = α * det(R),其中R是原矩阵的行梯形A,α是某个系数。

找到行梯形矩阵的行列式真的很容易;您只需找到对角线的乘积。求解原始矩阵 A 的行列式然后归结为计算 α,因为您找到行梯形形式 R

你需要知道的

什么是行梯形?

有关简单定义,请参阅此 [链接](http://stattrek.com/matrix-algebra/echelon-form.aspx)
**注意:** 并非所有定义都要求前导条目为 1,此算法不需要。

您可以使用基本行操作找到 R

交换行,添加另一行的倍数等。

您从行列式的行操作的属性中推导出 α

  1. 如果BA的一行乘以某个非零常数ß得到的矩阵,那么

    det(B) = ß * det(A)

    • 换句话说,您基本上可以从行中“分解”出一个常数,只需将其拉到行列式前面即可。
  2. 如果B是交换两行A得到的矩阵,那么

    det(B) = -det(A)

    • 如果交换行,请翻转标志。
  3. 如果BA中某一行的倍数与另一行相加得到的矩阵,则

    det(B) = det(A)

    • 行列式不变。

请注意,在大多数情况下,您可以仅使用规则 3 找到行列式(我相信当 A 的对角线没有零时),并且在所有情况下仅使用规则 2 和 3。规则 1 对人类有帮助在纸上做数学,尽量避免分数。

示例

(我做了一些不必要的步骤来更清楚地展示每条规则)

| 2 3 3 1 | A=| 0 4 3 -3 | | 2 -1 -1 -3 | | 0 -4 -3 2 | R2 R3, -α -> α (规则 2) | 2 3 3 1 | -| 2 -1 -1 -3 | | 0 4 3 -3 | | 0 -4 -3 2 | R2 - R1 -> R2(规则 3) | 2 3 3 1 | -| 0 -4 -4 -4 | | 0 4 3 -3 | | 0 -4 -3 2 | R2/(-4) -> R2, -4α -> α(规则 1) | 2 3 3 1 | 4| 0 1 1 1 | | 0 4 3 -3 | | 0 -4 -3 2 | R3 - 4R2 -> R3, R4 + 4R2 -> R4(规则 3,应用两次) | 2 3 3 1 | 4| 0 1 1 1 | | 0 0 -1 -7 | | 0 0 1 6 | R4 + R3 -> R3 | 2 3 3 1 | 4| 0 1 1 1 | = 4 ( 2 * 1 * -1 * -1 ) = 8 | 0 0 -1 -7 | | 0 0 0 -1 |
def echelon_form(A, size):
    for i in range(size - 1):
        for j in range(size - 1, i, -1):
            if A[j][i] == 0:
                continue
            else:
                try:
                    req_ratio = A[j][i] / A[j - 1][i]
                    # A[j] = A[j] - req_ratio*A[j-1]
                except ZeroDivisionError:
                    # A[j], A[j-1] = A[j-1], A[j]
                    for x in range(size):
                        temp = A[j][x]
                        A[j][x] = A[j-1][x]
                        A[j-1][x] = temp
                    continue
                for k in range(size):
                    A[j][k] = A[j][k] - req_ratio * A[j - 1][k]
    return A

【讨论】:

    【解决方案2】:

    Here 是一个广泛的讨论。

    有很多算法。

    一个简单的方法是使用LU decomposition。那么,由于

     det M = det LU = det L * det U
    

    并且LU 都是三角形的,行列式是LU 的对角元素的乘积。那是O(n^3)。存在更有效的算法。

    【讨论】:

      【解决方案3】:

      我对 LU 因式分解不太熟悉,但我知道为了得到 L 或 U,您需要使初始矩阵为三角形(U 为上三角或 L 为下三角)。然而,一旦你得到一些 nxn 矩阵 A 的三角形矩阵并假设你的代码使用的唯一操作是 Rb - k*Ra,你可以从 i=0 求解 det(A) = Π T(i,i)到 n (即 det(A) = T(0,0) x T(1,1) x ... x T(n,n)) 用于三角矩阵 T。检查此链接以了解我在说什么关于。 http://matrix.reshish.com/determinant.php

      【讨论】:

        【解决方案4】:

        如果您进行了初步研究,您可能会发现当 N>=4 时,矩阵行列式的计算变得相当复杂。关于算法,我会指出你Wikipedia article on Matrix determinants,特别是“算法实现”部分。

        根据我自己的经验,你可以在Alglib等现有矩阵库中轻松找到LU或QR分解算法。不过算法本身并不是很简单。

        【讨论】:

          猜你喜欢
          • 2020-10-21
          • 1970-01-01
          • 1970-01-01
          • 2015-01-16
          • 2013-09-22
          • 1970-01-01
          • 1970-01-01
          • 2015-10-20
          相关资源
          最近更新 更多