﻿ 130213033张艳作业七 - 计算数学达人 - 专，学者，数值代数，微分方程数值解

### 130213033张艳作业七

function [A,b,x]=matrix(n)

format long

x=rand(n,1)

N=-1*ones(n,n);

A=tril(N);

for i=1:n

A(i,i)=1;

end

A(:,n)=1;

b=A*x

function [L,U,q]=guass0(A)

% Gaussion Elimination with complete pivoting

[m,n]=size(A);

p=[];

tol=1e-16;

if m==n

for k=1:n-1

[r,c]=max(abs(A(k:n,k)));

q(k)=c+k-1;

A([k,q(k)],:)=A([q(k),k],:); %交换k,p

if abs(A(k,k))>tol

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

'A(k,k)=0!'

break;

end

end

else

'm~=n';

end

L=eye(n)+tril(A,-1);

U=triu(A);

%test_guass0

function x1=test_gauss0(A,b)

[L,U,q]=guass0（A);

n=size(A,1);

% the interchanges of b by p

for k=1:length(q)

b([k,q(k)])=b([q(k),k]);

end

%Ly=b

b(1)=b(1)/L(1,1);

for j=2:n

b(j)=(b(j)-L(j,1:j-1)*b(1:j-1))/L(j,j);

end

%Ux=y

b(n)=b(n)/U(n,n);

for i=n-1:-1:1

b(i)=(b(i)-U(i,i+1:n)*b(i+1:n))/U(i,i);

end

x1=b;%计算解

clc

clear all

format long

nn=[5:30];

error=[];

terror=[];

error1=[];

terror1=[];

for n=nn

[A,b,x]=matrix(n);

normb=norm(b,inf);

x1=test_gausslu_pp(A,b);

r=A*x-A*x1;

normr=norm(r,inf);

condA=cond(A,inf);

error=condA*(normr/normb);

error1=[error1;error];

terror=norm((abs(x1-x)),inf)/norm(x,inf);

terror1=[terror1;terror];

end

error1

terror1

err=abs(error1-terror1)