数值线性代数第二版徐树方高立张平文上机习题第二章实验报告.docx
《数值线性代数第二版徐树方高立张平文上机习题第二章实验报告.docx》由会员分享，可在线阅读，更多相关《数值线性代数第二版徐树方高立张平文上机习题第二章实验报告.docx（15页珍藏版）》请在冰豆网上搜索。
数值线性代数第二版徐树方高立张平文上机习题第二章实验报告
（1）估计5到20阶Hilbert矩阵的
范数条件数
（2）设
，先随机地选取
，并计算出
；然后再用列主元Gauss消去法求解该方程组，假定计算解为
。
试对n从5到30估计计算解
的精度，并且与真实相对误差作比较。
解
（1）分析：
利用
使
从5循环到20，利用
函数得到Hilbert矩阵
；先将算法编制成通用的子程序，利用算法编成的子程序
，对
求解，得到
的一个估计值
；再利用
得到
；则条件数
。
另，矩阵
的
范数条件数可由
直接算出，两者可进行比较。
程序为
1算法编成的子程序
functionv=opt（B）
k=1;
n=length（B）;
x=1./n*ones（n,1）;
whilek==1
w=B*x;
v=sign（w）;
z=B'*v;
ifnorm（z,inf）<=z'*x
v=norm（w,1）;
k=0;
else
x=zeros（n,1）;
[s,t]=max（abs（z））;
x（t）=1;
k=1;
end
end
end
2问题
（1）求解ex2_1
forn=5:
20
A=hilb（n）;
B=inv（A.'）;
v=opt（B）;
K1=v*norm（A,inf）;
K2=cond（A,inf）;
disp（['n=',num2str（n）]）
disp（['估计条件数为',num2str（K1）]）
disp（['实际条件数为',num2str（K2）]）
end
计算结果为
n=5
估计条件数为943656
实际条件数为943656
n=6
估计条件数为.0028
实际条件数为.0028
n=7
估计条件数为
实际条件数为
n=8
估计条件数为
实际条件数为
n=9
估计条件数为.422
实际条件数为.422
n=10
估计条件数为
实际条件数为
n=11
估计条件数为344
实际条件数为344
Warning:
Matrixisclosetosingularorbadlyscaled.Resultsmaybe
inaccurate.RCOND=.
>Inex2_1at3
Warning:
Matrixisclosetosingularorbadlyscaled.Resultsmaybe
inaccurate.RCOND=.
>Incondat47
Inex2_1at6
n=12
估计条件数为+16
实际条件数为+16
Warning:
Matrixisclosetosingularorbadlyscaled.Resultsmaybe
inaccurate.RCOND=.
>Inex2_1at3
Warning:
Matrixisclosetosingularorbadlyscaled.Resultsmaybe
inaccurate.RCOND=.
>Incondat47
Inex2_1at6
n=13
估计条件数为+18
实际条件数为+18
Warning:
Matrixisclosetosingularorbadlyscaled.Resultsmaybe
inaccurate.RCOND=.
>Inex2_1at3
Warning:
Matrixisclosetosingularorbadlyscaled.Resultsmaybe
inaccurate.RCOND=.
>Incondat47
Inex2_1at6
n=14
估计条件数为+17
实际条件数为+17
Warning:
Matrixisclosetosingularorbadlyscaled.Resultsmaybe
inaccurate.RCOND=.
>Inex2_1at3
Warning:
Matrixisclosetosingularorbadlyscaled.Resultsmaybe
inaccurate.RCOND=.
>Incondat47
Inex2_1at6
n=15
估计条件数为+17
实际条件数为+17
Warning:
Matrixisclosetosingularorbadlyscaled.Resultsmaybe
inaccurate.RCOND=.
>Inex2_1at3
Warning:
Matrixisclosetosingularorbadlyscaled.Resultsmaybe
inaccurate.RCOND=.
>Incondat47
Inex2_1at6
n=16
估计条件数为+17
实际条件数为+17
Warning:
Matrixisclosetosingularorbadlyscaled.Resultsmaybe
inaccurate.RCOND=.
>Inex2_1at3
Warning:
Matrixisclosetosingularorbadlyscaled.Resultsmaybe
inaccurate.RCOND=.
>Incondat47
Inex2_1at6
n=17
估计条件数为+18
实际条件数为+18
Warning:
Matrixisclosetosingularorbadlyscaled.Resultsmaybe
inaccurate.RCOND=.
>Inex2_1at3
Warning:
Matrixisclosetosingularorbadlyscaled.Resultsmaybe
inaccurate.RCOND=.
>Incondat47
Inex2_1at6
n=18
估计条件数为+18
实际条件数为+18
Warning:
Matrixisclosetosingularorbadlyscaled.Resultsmaybe
inaccurate.RCOND=.
>Inex2_1at3
Warning:
Matrixisclosetosingularorbadlyscaled.Resultsmaybe
inaccurate.RCOND=.
>Incondat47
Inex2_1at6
n=19
估计条件数为+18
实际条件数为+18
Warning:
Matrixisclosetosingularorbadlyscaled.Resultsmaybe
inaccurate.RCOND=.
>Inex2_1at3
Warning:
Matrixisclosetosingularorbadlyscaled.Resultsmaybe
inaccurate.RCOND=.
>Incondat47
Inex2_1at6
n=20
估计条件数为+18
实际条件数为+18
结果分析
随着矩阵阶数增加，估计值误差开始出现，
时估计条件数与实际值存在误差；且条件数很大，Hilbert矩阵为病态的。
解
（2）分析：
先根据题目要求，利用
和
使
从5循环到30，作出
和随机的
，并计算出
；然后再利用第一章习题中得到的
和
用列主元Gauss消去法求解该方程组，假定计算解为
，得
，利用第
（1）问所得函数
计算
的一个估计值，利用
计算
的无穷范数，则
的相对误差估计为
，真实相对误差为
。
程序为
1列主元Gauss消去法求解该方程组的程序为
的
分解：
function[L,U,P]=GaussCol（A）
n=length（A）;
fork=1:
n-1
[s,t]=max（abs（A（k:
n,k）））;
p=t+k-1;
temp=A（k,1:
n）;
A（k,1:
n）=A（p,1:
n）;
A（p,1:
n）=temp;
u（k）=p;
ifA（k,k）~=0
A（k+1:
n,k）=A（k+1:
n,k）/A（k,k）;
A（k+1:
n,k+1:
n）=A（k+1:
n,k+1:
n）-A（k+1:
n,k）*A（k,k+1:
n）;
else
break;
end
end
L=tril（A）;
U=triu（A）;
L=L-diag（diag（L））+diag（ones（1,n））;
P=eye（n）;
fori=1:
n-1
temp=P（i,:
）;
P（i,:
）=P（u（i）,:
）;
P（u（i）,:
）=temp;
end
end
高斯消去法解线性方程组
functionx=Gauss（A,b,L,U,P）
ifnargin<5
P=eye（length（A））;
end
n=length（A）;
b=P*b;
forj=1:
n-1
b（j）=b（j）/L（j,j）;
b（j+1:
n）=b（j+1:
n）-b（j）*L（j+1:
n,j）;
end
b（n）=b（n）/L（n,n）;
y=b;
forj=n:
-1:
2
y（j）=y（j）/U（j,j）;
y（1:
j-1）=y（1:
j-1）-y（j）*U（1:
j-1,j）;
end
y
（1）=y
（1）/U（1,1）;
x=y;
end
2问题
（2）求解ex2_2
forn=5:
30
A=2*eye（n）+tril（-1*ones（n））;A（1:
n-1,n）=ones（n-1,1）;
x=100*rand（n,1）;
b=A*x;
[L,U,P]=GaussCol（A）;x1=Gauss（A,b,L,U,P）;
r=b-A*x1;
p1=norm（r,inf）*opt（inv（A.'））*norm（A,inf）/norm（b,inf）;
p2=norm（x-x1,inf）/norm（x,inf）;
disp（['n=',num2str（n）]）
disp（['估计相对误差为',num2str（p1）]）
disp（['实际相对误差为',num2str（p2）]）
y1（n-4）=p1;y2（n-4）=p2;
end
plot（5:
30,y1,5:
30,y2）
legend（'估计相对误差','实际相对误差'）
计算结果为
n=5
估计相对误差为
实际相对误差为
n=6
估计相对误差为
实际相对误差为
n=7
估计相对误差为
实际相对误差为
n=8
估计相对误差为
实际相对误差为
n=9
估计相对误差为
实际相对误差为
n=10
估计相对误差为
实际相对误差为
n=11
估计相对误差为
实际相对误差为
n=12
估计相对误差为
实际相对误差为
n=13
估计相对误差为
实际相对误差为
n=14
估计相对误差为
实际相对误差为
n=15
估计相对误差为
实际相对误差为
n=16
估计相对误差为
实际相对误差为
n=17
估计相对误差为
实际相对误差为
n=18
估计相对误差为
实际相对误差为
n=19
估计相对误差为
实际相对误差为
n=20
估计相对误差为
实际相对误差为
n=21
估计相对误差为
实际相对误差为
n=22
估计相对误差为
实际相对误差为
n=23
估计相对误差为
实际相对误差为
n=24
估计相对误差为
实际相对误差为
n=25
估计相对误差为
实际相对误差为
n=26
估计相对误差为
实际相对误差为
n=27
估计相对误差为
实际相对误差为
n=28
估计相对误差为
实际相对误差为
n=29
估计相对误差为
实际相对误差为
n=30
估计相对误差为
实际相对误差为
结果分析
n较小时估计的较好，随着n的增大估计值误差增大
升级会员 