Fernando Henrique Gomes Zucatelli 1

domingo – 17/10/2010 – 15:31

1. CÁLCULO NUMÉRICO

Algoritmo para obter inversa de A após fatoração LU. PA=LU

Algoritmo:

function [L U P Z] = Inversa_Fatoracao_LU(A)

%

% Funcao para calcular a inversa de uma matriz

% usando a fatoração LU

% A: matriz dada

% M: recebe A para resolver o sistema

% U: recebe A para se tornar Upper

% L: recebe identidade com L linhas de A e é trabalhada até ser Lower

% P: recebe identidade e sofre as permutações.

% I: recebe identidade e sofre eliminação Gaussiana para ser inv(L)

% Z: recebe L quando U é Upper para se tornar inversa

[Lm Cm] = size(A);

L = eye(Lm); %criar a identidade de tamanho Linhas de A

U = [A];

P = eye(Lm);

%pivotiamento para P e U

for i=1:Lm

[el_max pos_max] = max(abs(A(i,:)));

if (abs(el_max)~=0)

if el_max ~= U(i,1)

U_temp=U(i,:);

U(i,:)=U(pos_max,:);

U(pos_max,:)=U_temp;

P_temp=P(i,:);

P(i,:)=P(pos_max,:);

P(pos_max,:)=P_temp;

end

else

break

fprintf(1,'Matriz é singular')

return;

end

end

% Eliminação de Gauss para U = Upper e L = Lower

for i=1:Lm

for j=i+1:Lm

m=U(j,i)/U(i,i);

U(j,:) = U(j,:) - m*U(i,:);

L(j,i) = m;

end

end

% U é Upper e L é lower

% Eliminação de Gauss para obter inv(L) sem usar o comando pronto do MATLAB

I = eye(Lm);

L_Aux = [L];

for i=1:Lm

for j=i+1:Lm

m=L_Aux(j,i)/L_Aux(i,i);

L_Aux(j,:) = L_Aux(j,:) - m*L_Aux(i,:);

I(j,:) = I(j,:) - m*I(i,:);

end

end

y = I*P; % A = P'*L*U => U*inv(A) = inv(L)*P = I*P => UZ=y

Z = zeros(Lm);

for i=1:Lm

M = [U y(:,i)];

[Lm Cm] = size(M);

z1 = zeros(Lm,1); %matriz coluna

z1(Lm) = M(Lm,Cm)/M(Lm,Lm);

for j=Lm-1:-1:1

z1(j) = (M(j,Cm) - M(j,j+1:Lm) * z1(j+1:Lm))/M(j,j);

end

Z (:,i) = Z(:,i) + z1;

end

Fernando Henrique Gomes Zucatelli 2

domingo – 17/10/2010 – 15:31

Resultados

a) C =

4 -1 2 -2

3 -1 0 0

2 3 1 0

0 7 1 1

Calculada pelo MATLAB format short (nota-se que é arredondado pelos zeros após o ponto) inv(C)=

-1.0000 -1.0000 4.0000 -2.0000

-3.0000 -4.0000 12.0000 -6.0000

11.0000 14.0000 -43.0000 22.0000

10.0000 14.0000 -41.0000 21.0000

Z é a inv(D) calculada pelo programa criado

Z =

-1.0000 -1.0000 4.0000 -2.0000

-3.0000 -4.0000 12.0000 -6.0000

11.0000 14.0000 -43.0000 22.0000

10.0000 14.0000 -41.0000 21.0000

b) D =

2 1 0 0

1 0 -1 1

0 1 1 1

-1 0 0 3b

calcula pelo MATLAB inv(D) =

3 -3 -3 2

-5 6 6 -4

4 -5 -4 3

1 -1 -1 1

Z é a inv(D) calculada pelo programa criado

Z =

3 -3 -3 2

-5 6 6 -4

4 -5 -4 3

1 -1 -1 1