作者：凯鲁嘎吉 - 博客园
http://www.cnblogs.com/kailugaji/
三、实验程序
五、解答（按如下顺序提交电子版）
1.(程序)
(1)LU分解源程序：
function [l,u]=lu12(a,n)
for k=1:n-1
for i=k+1:n
a(i,k)=a(i,k)/a(k,k);
for j=k+1:n
a(i,j)=a(i,j)-a(i,k)*a(k,j);
end
end
end
l=eye(n);
u=zeros(n,n);
for k=1:n
for i=k:n
u(k,i)=a(k,i);
end
end
for k=1:n
for j=1:k-1
l(k,j)=a(k,j);
end
end
(2)直接三角分解法源程序：
function [a,l,u,y,x]=direct_triangle(a,b,n)
%a为N*N矩阵，b为n*1列向量
for k=1:n-1
for i=k+1:n
a(i,k)=a(i,k)/a(k,k);
for j=k+1:n
a(i,j)=a(i,j)-a(i,k)*a(k,j);
end
end
end
l=eye(n);
u=zeros(n,n);
for k=1:n
for i=k:n
u(k,i)=a(k,i);
end
end
for k=1:n
for j=1:k-1
l(k,j)=a(k,j);
end
end
y=ones(n,1);
x=ones(n,1);
y(1,1)=b(1,1);
for i=2:n
s=0;
for k=1:i-1
s=s+l(i,k)*y(k,1);
end
y(i,1)=b(i,1)-s;
end x(n,1)=y(n,1)/u(n,n);
for j=n-1:-1:1
s1=0;
for k1=j+1:n
s1=s1+u(j,k1)*x(k1,1);
end
x(j,1)=(y(j,1)-s1)/u(j,j);
end
2.(运算结果)
（1）求一个4阶矩阵的LU分解。
>> a=[10,7,8,7;7,5,6,5;8,6,10,9;7,5,9,10];
>> [l,u]=lu12(a,4) l = 1.0000 0 0 0
0.7000 1.0000 0 0
0.8000 4.0000 1.0000 0
0.7000 1.0000 1.5000 1.0000 u = 10.0000 7.0000 8.0000 7.0000
0 0.1000 0.4000 0.1000
0 0 2.0000 3.0000
0 0 0 0.5000
>> a=[10 7 8 7;7 5 6 5;8 6 10 9;7 5 9 10];b=[32 23 33 31]';
>> [a,l,u,y,x]=direct_triangle(a,b,4) a = 10.0000 7.0000 8.0000 7.0000
0.7000 0.1000 0.4000 0.1000
0.8000 4.0000 2.0000 3.0000
0.7000 1.0000 1.5000 0.5000 l = 1.0000 0 0 0
0.7000 1.0000 0 0
0.8000 4.0000 1.0000 0
0.7000 1.0000 1.5000 1.0000 u = 10.0000 7.0000 8.0000 7.0000
0 0.1000 0.4000 0.1000
0 0 2.0000 3.0000
0 0 0 0.5000 y = 32.0000
0.6000
5.0000
0.5000 x = 1.0000
1.0000
1.0000
1.0000
比如，希尔伯特矩阵就是一个病态矩阵，在方程组问题求解之前，可以先判断其条件数是否较大。
源程序：hilbert.m:
function [A,cond1]=hilbert(k)
format rat
A=zeros(k,k);
for m=1:k
for n=1:k
A(m,n)=1/(m+n-1);
end
end
cond1=cond(A,inf);
运行结果：
>> [A,cond1]=hilbert(3)
A =
1 1/2 1/3
1/2 1/3 1/4
1/3 1/4 1/5
cond1 =
748
>> [A,cond1]=hilbert(4)
A =
1 1/2 1/3 1/4
1/2 1/3 1/4 1/5
1/3 1/4 1/5 1/6
1/4 1/5 1/6 1/7
cond1 =
28375
>> [A,cond1]=hilbert(5)
A =
1 1/2 1/3 1/4 1/5
1/2 1/3 1/4 1/5 1/6
1/3 1/4 1/5 1/6 1/7
1/4 1/5 1/6 1/7 1/8
1/5 1/6 1/7 1/8 1/9
cond1 =
943656
从结果可见希尔伯特矩阵是一个病态矩阵，用一般的直接法和迭代法会有较大的误差，甚至严重失真。