%LU分解,Doolittle分解
%L为单位下三角阵,U为上三角阵
clear;
clc;
A=[2,2,3;4,7,7;-2,4,5];%系数矩阵
b=[1;2;3];%n维向量
n=length(b);%方程个数n
x=zeros(n,1);%n行一列的零向量
A(2:n,1)=A(2:n,1)./A(1,1);
for i=2:n-1
A(i,i)=A(i,i)-sum(A(i,1:i-1)'.*A(1:i-1,i));
for j=i+1:n
A(i,j)=A(i,j)-sum(A(i,1:i-1)'.*A(1:i-1,j));
A(j,i)=(A(j,i)-sum(A(j,1:i-1)'.*A(1:i-1,i)))/A(i,i);
end
end
A(n,n)=A(n,n)-sum(A(n,1:n-1)'.*A(1:n-1,n));
A
U=A;L=A;
for i=1:n
L(i,i)=1;
end
for i=1:n-1
for j=i+1:n
L(i,j)=0;
end
end
L %单位下三角阵
for i=2:n
for j=1:i-1
U(i,j)=0;
end
end
U %上三角阵???????
%Ly=b,向后迭代;
y=zeros(n,1);%n行1列的零向量
y(1)=b(1);
for i=2:n
y(i)=b(i)-sum(L(i,1:i-1)'.*y(1:i-1));
end
y
%输出y
%Ux=y,向前迭代
x(n)=y(n)/U(n,n);
for i=n-1:-1:1
x(i)=(y(i)-sum(U(i,i+1:n)'.*x(i+1)))/U(i,i);
end
x
%输出x