数值分析上机作业总.docx

数值分析上机实验
一、解线性方程组直接法（教材49页14题）
追赶法程序如下：
functionx=followup(A,b)
n=rank(A);
for(i=1:n)
if(A(i,i)==0)
disp('Error:对角有元素为0');
return;
end
end;
d=ones(n,1);
a=ones(n-1,1);
c=ones(n-1);
for(i=1:n-1)
a(i,1)=A(i+1,i);
c(i,1)=A(i,i+1);
d(i,1)=A(i,i);
end
d(n,1)=A(n,n);
for(i=2:n)
d(i,1)=d(i,1)-(a(i-1,1)/d(i-1,1))*c(i-1,1);
b(i,1)=b(i,1)-(a(i-1,1)/d(i-1,1))*b(i-1,1);
end
x(n,1)=b(n,1)/d(n,1);
for(i=(n-1):-1:1)
x(i,1)=(b(i,1)-c(i,1)*x(i+1,1))/d(i,1);
end
主程序如下：
functionzhunganfa
A=[2-2000000;-25-200000;0-25-20000;00-25-2000;000-25-200;0000-25-20;00000-25-2;000000-25];
b=[220/27;0;0;0;0;0;0;0];
x=followup(A,b)
计算结果：

x=
8.1478
4.0737
2.0365
1.0175
0.5073
0.2506
0.1194
0.0477
二、解线性方程组直接法（教材49页15题）
程序如下：
functiontiaojianshu(n)
A=zeros(n);
forj=1:1:n
fori=1:1:n
A(i,j)=(1+0.1*i)^(j-1);
end
end
c=cond(A)
d=rcond(A)
当n=5时
c=
5.3615e+005
d=
9.4327e-007
当n=10时
c=
8.6823e+011
d=
5.0894e-013
当n=20时
c=
3.4205e+022
d=
8.1226e-024
备注：对于病态矩阵A来说，d为接近0的数；对于非病态矩阵A来说，d为接近1的数。
三、解线性方程组的迭代法（教材74页14题）
（1）用Jacobi迭代法求：
Jacobi迭代法程序如下：
function[x,n]=jacobi(A,b,x0,eps,varargin)
ifnargin==3
eps=1.0e-6;
M=200;
elseifnargin<3
error
return
elseifnargin==5
M=varargin{1};
end
D=diag(diag(A));
L=-tril(A,-1);
U=-triu(A,1);
B=D\(L+U);
f=D\b;
x=B*x0+f;
n=1;
whilenorm(x-x0)>=eps
x0=x;
x=B*x0+f;
n=n+1;
if(n>=M)
disp('Warning:迭代次数太多，可能不收敛');
return;
end
end
本题主程序如下：
functionyakebidiedai
A=[101234;19-12-3;2-173-5;32312-1;4-3-5-115];
b=[12;-27;14;-17;12];
x0=[0;0;0;0;0];
[x,n]=jacobi(A,b,x0)
计算结果：
x=
1.0000
-2.0000
3.0000
-2.0000
1.0000
n=

67
经过67次迭代，得到最终结果
（2）用Gauss-Seidel迭代法求：
Gauss-Seidel迭代法程序如下：
function[x,n]=gauseidel(A,b,x0,eps,M)
ifnargin==3
eps=1.0e-6;
M=200;
elseifnargin==4
M=200;
elseifnargin<3
error
return;
end
D=diag(diag(A));
L=-tril(A,-1);
U=-triu(A,1);
G=(D-L)\U;
f=(D-L)\b;
x=G*x0+f;
n=1;
whilenorm(x-x0)>=eps
x0=x;
x=G*x0+f;
n=n+1;
if(n>=M)
disp('Warning:迭代次数太多，可能不收敛');
return;
end
end
本题主程序如下：
functiongaosidiedai
A=[101234;19-12-3;2-173-5;32312-1;4-3-5-115];
b=[12;-27;14;-17;12];
x0=[0;0;0;0;0];
[x,n