% Test 1, calculo numerico (521230), S2-2011
% Tema 7: Sistemas de ecuaciones lineales
% Pregunta 2

% A5 y T5

[A,T] = matrizT(5);

% vector parte derecha

b = ones(20,1);

% solucion al sistema con metodo directo

x = T\b;

% descomposicion LU de A
[L,U,P] = lu(A);

% Inicializando u1,u2,u3,u4
u1 = zeros(5,1);
u2 = zeros(5,1);
u3 = zeros(5,1);
u4 = zeros(5,1);

% Calculando y1,y2,y,y4
y1 = b(1:5);
y2 = b(6:10);
y3 = b(11:15);
y4 = b(16:20);

% Jacobi por bloques
for i=1:20
    u1n = A\(u2 + y1);
    u2n = A\(u1 + u3 + y2);
    u3n = A\(u2 + u4 + y3);
    u4n = A\(u3 + y4);
    u1 = u1n;
    u2 = u2n;
    u3 = u3n;
    u4 = u4n;
end

% solucion del sistema Tx = b mediante Jacobi por bloques
xj = [u1;u2;u3;u4];

% norma diferencia entre soluciones
norm(x-xj)

% mejorando iteracion anterior
% hay al menos 2 formas de hacerlo
% una es reusando la descomposicion LU de A calculada antes (en
% comentarios)
% la otra es, resolviendo los 4 sistemas de cada iteracion simultaneamente
% la primera es la mas eficiente de las 2, pero se acepta cualquiera de las
% dos
for i=1:20
%     u1n = U\(L\(P*(u2 + y1)));
%     u2n = U\(L\(P*(u1 + u3 + y2)));
%     u3n = U\(L\(P*(u2 + u4 + y3)));
%     u4n = U\(L\(P*(u3 + y4)));
    B = [u2+y1 u1+u3+y2 u2+u4+y3 u3+y4];
    z = A\B;
    u1n = z(:,1);
    u2n = z(:,2);
    u3n = z(:,3);
    u4n = z(:,4);
    u1 = u1n;
    u2 = u2n;
    u3 = u3n;
    u4 = u4n;
end

xjm = [u1;u2;u3;u4];

norm(xjm-x)