基于MATLAB的Cholesky分解法解线性方程组

定理

若A∈Rn∗n对称正定,则存在一个对角元为正数的下三角矩阵L∈Rn∗n,使得A=LLT成立.假设现在要求求解线性方程组AX=B,其中A为对称正定矩阵,那么X可以通过下面步骤求解

(1). 求A的Cholesky分解,得到A=LLT

(2).求解LY=B,得到Y

(3).求解LTX = Y,得到X

现在的关键问题是对A进行Cholesky分解

解释:如果矩阵X是对称正定的，则Cholesky分解将矩阵X分解成一个下三角矩阵和上三角矩阵的乘积。设上三角矩阵为R，则下三角矩阵为其转置，即X=R’R。

实现

1).自己实现:下面是自己实现的求得一个上三角举证的matlab函数

function [u] = cholesky( A )
if A~=A'
error('Invalid Argument')
end
[n,n] = size(A);
X = zeros(n,1);
Y = zeros(n,1);
B = zeros(n,n);
for k = 1:n
A(k,k) = sqrt(A(k,k))
if A(k,k)==0
'A is singular.no unique solution'
break;
end
for i = k+1:n
A(i,k) = A(i,k)/A(k,k)
end
for j=k+1:n
for i = k+1 :n
A(i,j) = A(i,j)-A(i,k)*A(i,k)
end
end
end
for x = 1:n
for y = 1:n
B(x,y) = A(x,y)'
if x<y
B(x,y) = 0;
break;
end
end
end
u = B

2).matlab自带库函数

MATLAB函数chol(X)用于对矩阵X进行Cholesky分解，其调用格式为：

R=chol(X)：产生一个上三角阵R，使R’R=X。若X为非对称正定，则输出一个出错信息。

[R,p]=chol(X)：这个命令格式将不输出出错信息。当X为对称正定的，则p=0，R与上述格式得到的结果相同；否则p为一个正整数。如果X为满秩矩阵，则R为一个阶数为q=p-1的上三角阵，且满足R’R=X(1:q,1:q)。