本站原创文章，转载请说明来自《老饼讲解-机器学习》www.bbbdata.com
Crout分解是矩阵LU分解的一种，它将矩阵A分解为L*U的形式，且U的对角元素全为1
本文介绍Crout分解的定义、分解的流程与公式推导，并展示Crout分解的具体代码实现
通过本文，可以快速了解Crout分解是什么，Crout分解的方法，以及如何用代码实现Cholesky分解
本节介绍Crout分解是什么，以及Crout分解的分解方法
LU分解与Crout分解
LU分解是什么
LU分解是指将方阵A分解为下三角矩阵L与上三角矩阵U的乘积
即将方阵A如下分解：
其中，L是下三角矩阵，U是上三角矩阵
Crout分解是什么
Crout分解则是LU分解方法中的一种，它的特点是U的对角元素全为1
即Crout分解把方阵A分解如下：
其中
，
✍️关于LU分解的类别
LU分解从L、U的特点可以分为三种
1. Doolittle分解 ：L对角元全为1(即单位下三角矩阵)
2. Crout分解 ：U对角元全为1(即单位上三角矩阵)
3. Cholesky分解 ：U是L的转置(要求A必须是正定对称方阵)
Crout分解-方法与流程
Crout分解通过"对U逐行、对L逐列"的迭代方式来对L、U进行求解
即先算U的第1行，再算L的第1列
然后算U的第2行，再算L的第2列....如此类推
L的计算公式
在已算出L的前(ｉ-1)列，U的(ｉ-1)行时
L的第ｉ列计算公式如下：
特别地，L的第1列元素
U的计算公式
在已算出L的前ｉ列，U的前(ｉ-1)行时
U的第ｉ行计算公式如下：
特别地，U的第1行元素
本节介绍Crout分解方法中的公式是如何推导出来的
L的求解公式-推导
L第ｉ列的求解，可由下式得到：
A的第i列为 L乘以U 的i列
表示为连加的形式
U是上三角矩阵，k>i时为0
独立抽出第i项
这是因为Uii=1
根据上式可得到L第ｉ列的求解公式为：
因此，在L的前(ｉ-1)列和U的前(ｉ-1)行已求出
即可根据上式求出L第ｉ列的所有元素，
👉补充：上述推导原理对于L的第一列同样适用
U的求解公式-推导
U第ｉ行的求解，可由下式得到：
A的第i行为 L的i行乘以U
表示为连加的形式
这是因为L是下三角矩阵，k>i时为0
独立抽出第i项
根据上式可得到U第ｉ行的求解公式为：
因此，在L的前ｉ列和U的前(ｉ-1)行已求出
即可根据上式求出U第i行的所有元素
👉补充：上述推导原理对于U的第一行同样适用
本节展示Crout矩阵分解的具体实现代码
Crout矩阵分解-代码实现
按上文的Crout分解流程，使用python编写Crout分解的程序
具体代码如下：
# -*- coding: utf-8 -*-
"""
本代码用于展示如何实现Crout矩阵分解(LU分解的一种)
本代码来自《老饼讲解-机器学习》www.bbbdata.com
"""
import numpy as np
# Crout矩阵分解函数
def CroutDecompose(A):
size = A.shape[0] # 方阵的大小
L = np.zeros((size,size)) # 初始化L全为0
U = np.identity(size) # 初始化U为单位矩阵
L[:,0] = A[:,0].copy() # 计算L的第一列
U[0,:] = A[0,:]/L[0,0] # 计算U的第一行
for i in range(1,size): # 逐列、行计算L和U
L[:,i] = A[:,i] - L[:,:i]@U[:i,i] # 计算L的第i列
U[i,:] = (A[i,:] - L[i,:i]@U[:i,:])/L[i,i] # 计算U的第i行
# 为避免产生数值问题,将L上三角元素、U的下三角元素置0
L=np.tril(L) # L只取下三角
U=np.triu(U) # U只取上三角
return L,U
# 测试样例
if __name__ == "__main__":
A = np.array([[1.,2.,5,8],[3.,5.,4,2],[6.,4,3,1],[2.,3,5,5]]).T # 生成需要分解的矩阵A
L,U= CroutDecompose(A) # Crout矩阵分解
print('\n-----需要分解的矩阵A-----') # 打印标题
print('A=\n',A) # 打印矩阵A
print('\n--------分解的结果------') # 打印标题
print('L=\n',L) # 打印分解的L
print('\nU=\n',U) # 打印分解的U
print('\n---------结果验证--------') # 打印标题
print('L*U=\n',L@U) # 打印分解的LU乘积代码运行结果如下：
从运行结果中可以看到，从A中分解出下三角矩阵L和上三角矩阵U
将L与U相乘后，的确是等于A，说明程序进行的Crout分解是无误的
好了，以上就是矩阵LU分解中的Crout分解了~
End