?0.02x1?x2?4x3?3x4?x5?11??x?x?2x?x?3x?1412345???4x1?2x2?3x3?3x4?x5?4??3x?x?3x?2x?4x?1612345???x1?3x2?x3?4x4?4x5?18
解:
Matlab程序如下:
A=[2 -1 4 -3 1;-1 1 2 1 3;4 2 3 3 -1;-3 1 3 2 4;1 3 -1 4 4] b=[11 14 4 16 18] function x=zuoye6(A,b) [n,v]=size(b);
D=[A,b;eye(n),zeros(n,v)],[s1,m]=size(D); for k=1:(n-1)
s=abs(A(k,k));p=k;q=k; for i=k:n for j=k:n if abs(A(i,j))>s
s=abs(A(i,j));p=i;q=j; end end end
if p>k t=D(k,:); D(k,:)=D(p,:); D(p,:)=t; end if q>k t1=D(:,k); D(:,k)=D(:,q); D(:,q)=t1; end h=D(k+1:n,k)/D(k,k);
D(k+1:n,k+1:m)=D(k+1:n,k+1:m)-h*D(k,k+1:m); D(k+1:n,k)=zeros(n-k,1);
6
end for k=n:-1:1
D(k,k:m)=D(k,k:m)/D(k,k);
for r=1:k-1 D(r,:)=D(r,:)-D(r,k)*D(k,:); end end
P=D(n+1:2*n,1:n);Q=D(1:n,n+1:m);x=P*Q 在Matlab窗口中执行:
A=[0.02 -1 4 -3 1;-1 1 2 1 3;4 2 3 3 -1;-3 1 3 2 4;1 3 -1 4 4]; b=[11 14 4 16 18]'; zuoye6(A,b) 运行结果如下: x =
2.94117647058824 -3.82352941176472 1.00000000000000 0.94117647058824 5.94117647058824
0?2?10??12?10??0?12?1?0?12??0?000?17. 用追赶法解线性方程组 ?
7
0??x1???x?0??2???x3??0????x4??1?2???x5??????????1??0??0?0??0 ?要求给出程序和运行结果. 解:
?1?1??00??2?100100000??2??0??0?130???100?00???2A???12?100???2??2?0?12?10???03100?????00?00?12?1?LU???????000?12?0?310????0?4??0???0??000?451???????00{Ly?b于是有求解Ax=b即为求解Ux?y,式中b=(1 0 0 0 0)T
??10000??1???1???1000????1???2??0?23100??y1????1???2???y2??0??1????y??3???0??00?310????y??40??3?1???4?????y5??????0????4??000?4??据
?51??1?? y=??5??
?5??2?1000?????1???6??03?1?????2?100????2???4??2?3?1???003?10???x???1??1?2?????x?2?0005?1????3?1???1??4??6??x3????x???4???3?4据??0000?5?????1???x5??=??5?? x=?1??6??
Matlab程序如下:
8
?43?10???05?4?1?006??5?? %定义zuoye7.m文件 function x=zuoye7(a,b,c,d) a1=[0;a];
n=length(b);q=zeros(n,1);p=zeros(n,1); %LU分解
q(1)=b(1);for k=2:n,p(k)=a1(k)/q(k-1); q(k)=b(k)-p(k)*c(k-1); end %解Ly=d
y=zeros(n,1);y(1)=d(1);for k=2:n, y(k)=d(k)-p(k)*y(k-1);end %解Ux=y
x=zeros(n,1); x(n)=y(n)/q(n);for k=n-1:-1:1,x(k)=(y(k)-c(k)*x(k+1))/q(k);end x
在Matlab窗口中执行: a=[-1 -1 -1 -1]'; b=[2 2 2 2 2]'; c=[-1 -1 -1 -1]'; d=[1 0 0 0 0]'; x=zuoye7(a,b,c,d) 运行结果如下: x = 0.83333333333333 0.66666666666667 0.50000000000000 0.33333333333333 0.16666666666667
9
9. 分别用Jacobi迭代法和Gauss-Seidel迭代法求解程组(编写程序)
?10x1+3x2?x3?14????2x1-10x2?3x3?-5? ?x?3x?10x?14?23?1?T取初值X(0)=(0,0,0),精确到小数后面四位。
解:
(1) Jacobi迭代法的Matlab程序:
% x0为初始向量,ep为精度,N为最大次数,x是近似解向量 A=[10 3 1;2 -10 3;1 3 10];
b=[14 -5 14];n=length(b);N=500;ep=1e-6;x0=zero(n,1); n=length(b);x0=zeros(n,1);x=zeros(n,1);k=0 while K x(i)=(b(i)-A(I,[1:i-1,i+1:n])*x0(1:i-1,i+1:n))/A(i,i); end if norm(x-x0,inf) if k==N,Warning(‘已达到迭代次数上限’);end disp([‘k=’,num2str(k)]),x 运行结果如下: k=15 x= 1.000000327906423 10