Doolittle方法（即LU分解）及Python实现

目录

一、LU分解原理

二、LU分解过程

三、Python实现完整代码

四、矩阵的一些补充概念

---

一、LU分解原理

1.

Gauss elimination takes more computational time for higher-order matrices. To reduce time consumption a matrix can be decomposed using LU factoring methods.

A system of linear equations can be written in a matrix form Ax = b where |A| ≠0

即系数矩阵A需要是可逆矩阵。

2.

Using LU factorization A can be decomposed as A = LU where L is lower triangular matrix and U is an upper triangular matrix.

由可逆矩阵的性质，若A是可逆矩阵，L、U均为方阵，则L、U也是可逆矩阵。

这里L为下三角矩阵，U为上三角矩阵。

另外，并不是每个矩阵都有LU分解。但是这些矩阵可借由排列其各行顺序来解决，所以最终会得到一个PLU 分解。

3.

Hence we have

If we let

then we have

二、LU分解过程

1.

Decompose matrix A as where

2.

Determine all values of unknown entries in L and U

这一步解并没有它看起来那么麻烦。

Alternatively we can also find L and U by using elementary row operation (addition of 2 rows)

我们也可以用仅进行列变换的方法来得到L，用仅进行行变换的方法来得到U

3.

Solve

using forward substitution

using backward substitution

三、Python实现完整代码

#LU分解求Ax=b
#过程：A=LU,LZ=b,Ux=Z

import numpy as np
import pandas as pd
from scipy import linalg

np.random.seed(2)
def LU_decomposition(A,b):
    n=len(A[0])
    L = np.zeros([n,n])
    U = np.zeros([n, n])
    for i in range(n):
        L[i][i]=1
        if i==0:
            U[0][0] = A[0][0]
            for j in range(1,n):
                U[0][j]=A[0][j]
                L[j][0]=A[j][0]/U[0][0]
        else:
                for j in range(i, n):#U
                    temp=0
                    for k in range(0, i):
                        temp = temp+L[i][k] * U[k][j]
                    U[i][j]=A[i][j]-temp
                for j in range(i+1, n):#L
                    temp = 0
                    for k in range(0, i ):
                        temp = temp + L[j][k] * U[k][i]
                    L[j][i] = (A[j][i] - temp)/U[i][i]
    
    Z=linalg.solve(L, b)
    x=linalg.solve(U, Z)
    
    print("L=",L)
    print("U=",U)
    print("Z=",Z)
    print("Exact solution: x=",x)
    
    return

#定义A,b的值
A=[[-3,-5,10,13],[2,0,-4,6],[8,-2,1,-3],[5,-1,1,15]]
b=[-14,8,7,0]

LU_decomposition(A,b)

四、矩阵的一些补充概念

pivot:
Find the entry in the left column with the largest absolute value. This entry is called the pivot.

即首元，取当前列的最大元素

naive Gaussian elimination:
meaning no row exchanges allowed

gaussian elimination with partial pivoting

gaussian elimination with full pivoting

gaussian elimination with scaled partial pivoting

Crout method

Cholesky method