首页 文章

以递归方式计算矩阵（nxn）的行列式

提问于
浏览 次
4

我将编写一些代码来计算方阵（nxn）的行列式，使用拉普拉斯算法（意义递归算法）编写Wikipedia's Laplace Expansion .

我已经有了类 Matrix ，其中包括 init ， setitem ， getitem ， repr 以及计算行列式所需的所有东西（包括 minor(i,j) ） .

所以我尝试了下面的代码：

def determinant(self,i=0)  # i can be any of the matrix's rows
    assert isinstance(self,Matrix)
    n,m = self.dim()    # Q.dim() returns the size of the matrix Q
    assert n == m
    if (n,m) == (1,1):
        return self[0,0]
    det = 0
    for j in range(n):
        det += ((-1)**(i+j))*(self[i,j])*((self.minor(i,j)).determinant())
    return det

正如预期的那样，在每次递归调用中， self 都变成了一个合适的小调 . 但是当从递归调用返回时，它没有原始矩阵 . 这会在 for 循环中导致麻烦（当函数到达 (n,m)==(1,1) 时，返回矩阵的这一个值，但在 for 循环中， self 仍然是1x1矩阵 - 为什么？）

python recursion

4 回答

  • -2

    我发布了这段代码，因为我无法在互联网上使用它，如何使用标准库来解决n * n行列式 . 目的是与那些觉得有用的人分享 . 我开始计算与a（0，i）相关的子矩阵Ai . 我使用递归行列式使它变短 .

    def submatrix(M, c):
    B = [[1] * len(M) for i in range(len(M))]


    for l in range(len(M)):
        for k in range(len(M)):
            B[l][k] = M[l][k]

    B.pop(0)

    for i in range(len(B)):
        B[i].pop(c)
    return B


def det(M):
    X = 0
    if len(M) != len(M[0]):
        print('matrice non carrée')
    else:
        if len(M) <= 2:
            return M[0][0] * M[1][1] - M[0][1] * M[1][0]
        else:
            for i in range(len(M)):
                X = X + ((-1) ** (i)) * M[0][i] * det(submatrix(M, i))
    return X

    对不起之前没有评论:)如果你需要任何进一步的解释，请不要犹豫 .

    回复于
  • -2

    你确定 minor 返回一个新对象而不是对原始矩阵对象的引用吗？我使用了您的确切行列式方法并为您的类实现了 minor 方法，它对我来说很好 .

    下面是矩阵类的快速/脏实现，因为我不太重要，因为我们正在处理决定因素 . 注意 det 方法，与你的方法相同，并且 minor 方法（其余的方法都是为了方便实现和测试）：

    class matrix:
    def __init__(self, n):
        self.data = [0.0 for i in range(n*n)]
        self.dim = n
    @classmethod
    def rand(self, n):
        import random
        a = matrix(n)
        for i in range(n):
            for j in range(n):
                a[i,j] = random.random()
        return a
    @classmethod
    def eye(self, n):
        a = matrix(n)
        for i in range(n):
            a[i,i] = 1.0
        return a        
    def __repr__(self):
        n = self.dim
        for i in range(n):
            print str(self.data[i*n: i*n+n])
        return ''    
    def __getitem__(self,(i,j)):
        assert i < self.dim and j < self.dim
        return self.data[self.dim*i + j]
    def __setitem__(self, (i, j), val):
        assert i < self.dim and j < self.dim
        self.data[self.dim*i + j] = float(val)
    #
    def minor(self, i,j):
        n = self.dim
        assert i < n and j < n
        a = matrix(self.dim-1)
        for k in range(n):
            for l in range(n):
                if k == i or l == j: continue
                if k < i:
                    K = k
                else:
                    K = k-1
                if l < j:
                    L = l
                else:
                    L = l-1
                a[K,L] = self[k,l]
        return a
    def det(self, i=0):
        n = self.dim    
        if n == 1:
            return self[0,0]
        d = 0
        for j in range(n):
            d += ((-1)**(i+j))*(self[i,j])*((self.minor(i,j)).det())
        return d
    def __mul__(self, v):
        n = self.dim
        a = matrix(n)
        for i in range(n):
            for j in range(n):
                a[i,j] = v * self[i,j]
        return a
    __rmul__ = __mul__

    现在进行测试

    import numpy as np
a = matrix(3)
# same matrix from the Wikipedia page
a[0,0] = 1
a[0,1] = 2
a[0,2] = 3
a[1,0] = 4
a[1,1] = 5
a[1,2] = 6
a[2,0] = 7
a[2,1] = 8
a[2,2] = 9
a.det()   # returns 0.0
# trying with numpy the same matrix
A = np.array(a.data).reshape([3,3])
print np.linalg.det(A)  # returns -9.51619735393e-16

    在numpy情况下的残差是因为它通过（高斯）消除方法而不是拉普拉斯展开计算行列式 . 您还可以比较随机矩阵的结果，看看您的行列式函数和numpy 's doesn' t之间的差异超出 float 精度：

    import numpy as np
a = 10*matrix.rand(4)
A = np.array( a.data ).reshape([4,4])
print (np.linalg.det(A) - a.det())/a.det() # varies between zero and 1e-14
    回复于
  • 1

    这是python 3中的函数 .

    注意：我使用一维列表来容纳矩阵，而size数组是方形数组中的行数或列数 . 它使用递归算法来查找行列式 .

    def solve(matrix,size):
    c = []
    d = 0
    print_matrix(matrix,size)
    if size == 0:
        for i in range(len(matrix)):
            d = d + matrix[i]
        return d
    elif len(matrix) == 4:
        c = (matrix[0] * matrix[3]) - (matrix[1] * matrix[2])
        print(c)
        return c
    else:
        for j in range(size):
            new_matrix = []
            for i in range(size*size):
                if i % size != j and i > = size:
                    new_matrix.append(matrix[i])

            c.append(solve(new_matrix,size-1) * matrix[j] * ((-1)**(j+2)))

        d = solve(c,0)
        return d
    回复于
  • 0

    嘿，我已经使用递归函数在MATLAB中编写了代码 . 这可能对您有所帮助 .

    function value = determinant(A)
% Calculates determinant of a square matrix A.
% This is a recursive function. Not suitable for large matrices.

[rows, columns] = size(A);
if rows ~= columns
    error('input matrix is not a square matrix.')
end

value = 0;
if rows == 2
    for i = 1:rows
        value = A(1,1)*A(2,2) - A(1,2)*A(2,1);
    end
else
    for i = 1:rows
        columnIndices = [1:i-1 i+1:rows];
        value = value + (-1)^(i+1)*A(1,i)*... 
            determinant(A(2:rows, columnIndices));
    end
end
end
    回复于

相关问题