Exercício - Método de Gauss com substituições ascendentes

Escreva um procedimento que, dada uma matriz A e um vector b, determine a solução

do sistema Ax=b aplicando o método de Gauss

> restart;
with(linalg):

Warning, new definition for norm

Warning, new definition for trace

> gauss:=proc(M::matrix,v::vector)
local n,i,k,lik,x,s,j,A,b;
n:=rowdim(M);
A:=evalm(M);
b:=evalm(v);
x:=vector(n);

### ETAPA 1: triangularização da matriz ###

for k from 1 to n-1 do
for i from k+1 to n do

lik:=evalf(A[i,k]/A[k,k]);
A:=addrow(A,k,i,-lik);
b[i]:=b[i]-b[k]*lik;

od;
od;
print(A); print(b);

### ETAPA 2: resolução pelo método das substituições ascendentes ###

x[n]:=evalf(b[n]/A[n,n]);

for i from n-1 to 1 by -1 do

s:=0;
for j from i+1 to n do

s:=evalf(s+A[i,j]*x[j]);
od;
x[i]:=evalf((b[i]-s)/A[i,i]);
od;

RETURN(evalm(x));

end:

> ####################################################################

> Gauss:=proc(M::matrix,v::vector)
local b,n,i,k,lik,x,s,j,A;
n:=rowdim(M);
A:=evalm(M);
b:=evalm(v);
x:=vector(n);
for k from 1 to n-1 do
for i from k+1 to n do
lik:=A[i,k]/A[k,k];
A:=addrow(A,k,i,-lik);
b[i]:=b[i]-b[k]*lik;
od;
od;
x[n]:=b[n]/A[n,n];
for i from n-1 to 1 by -1 do
s:=0;
for j from i+1 to n do
s:=s+A[i,j]*x[j];
od;
x[i]:=(b[i]-s)/A[i,i];
od;
RETURN(evalm(x));
end:

> A:=matrix(3,3,[4,-9,2,2,-4,4,-1,2,2]);
b:=vector(3,[2,3,1]);

> gauss(A,b);

> evalm(A);evalm(b);

> Gauss(A,b);

> evalm(A);evalm(b);

>