0

我需要实现 LU 分解，然后将其与np.linalg.solvenumpy 中的函数进行比较。

代码中的函数（见下文）运行没有任何问题，但是当我使用它来求解矩阵时，我不断收到错误消息：

IndexError: list index out of range

在线上： 

L[i][j] = (A2[i][j] - s2) / U[j][j]

这是整个代码：

 def matrixMul(A, B):
        TB = zip(*B)
        return [[sum(ea*eb for ea,eb in zip(a,b)) for b in TB] for a in A]

    def pivotize(m):
        #Creates the pivoting matrix for m.
        n = len(m)
        ID = [[float(i == j) for i in range(n)] for j in range(n)]
        for j in range(n):
            row = max(range(j, n), key=lambda i: abs(m[i][j]))
            if j != row:
                ID[j], ID[row] = ID[row], ID[j]
        return ID

    def lu(A):
        #Decomposes a nxn matrix A by PA=LU and returns L, U and P.
        n = len(A)
        L = [[0.0] * n for i in range(n)]
        U = [[0.0] * n for i in range(n)]
        P = pivotize(A)
        A2 = matrixMul(P, A)
        for j in range(n):
            L[j][j] = 1.0
            for i in range(j+1):
                s1 = sum(U[k][j] * L[i][k] for k in range(i))
                U[i][j] = A2[i][j] - s1
            for i in range(j, n):
                s2 = sum(U[k][j] * L[i][k] for k in range(j))
                L[i][j] = (A2[i][j] - s2) / U[j][j]
        return (L)

    A = np.array([[1,1,3],[5,3,1],[2,3,1]])
    b = np.array([2,3,-1])
    print('LU factorization: ', lu(A))

    A = np.array([[1,1,3],[5,3,1],[2,3,1]])
    b = np.array([2,3,-1])
    print('Internal solver : ', np.linalg.solve(A,b))

有任何想法吗？谢谢！

4

1 回答 1

0

你的matrixMul方法不对。尝试这个：

matrixMul([[1, 0], [0, 1]], [[1, 0], [0, 1]])

这是将两个单位矩阵相乘，应该返回[[1, 0], [0, 1]]. 当我运行它时，它返回[[1, 0], []]. 这意味着A2，您的内部lu只有一行，其余为空。i == 1因此，当和时的索引错误j == 0。

这个失败的原因是那TB是一个zip对象。这些只能迭代一次，消耗迭代器。我实际上并不认为您根本需要该TB对象，只需遍历B.

def matrixMul(A, B):
    return [[sum(ea*eb for ea,eb in zip(a,b)) for b in B] for a in A]

这将返回所需的输出：

>>> matrixMul([[1, 0], [0, 1]], [[1, 0], [0, 1]])
>>> [[1, 0], [0, 1]]

编辑：

顺便说一句，您的解决方案还有其他问题。你的和 NumPy 的版本仍然不匹配。但是这里的解决方案将修复您的索引错误。

于 2017-01-21T21:42:01.330 回答