Metodo de las potencias para valores propios, Algoritmos para métodos de la Potencia para hallar el valor propio eigenvalue dominante

 %Algoritmos para métodos de la Potencia para hallar el valor propio eigenvalue dominante

%Autor: Alexander Arias

%Fecha: 2023-04-18

%https://www.geocities.ws/hectormora/mns.pdf

%pag 320

disp("Método de Potencias");

clear all;

clc;

lambda=[]

A=[-1 -2 -3;-4 -5 -6;-7 -8 -8]

x=[1;1;1]

%--------------------

disp("Método de Potencias sin normalizar");

printf(" i        x1         x2         x3      lambda\n")

x=[1;1;1];

for i=1:12

    z=A*x;

    lambda=z(1)/x(1);

    x=z;

    printf("%2d  %-+2.3e  %-+3.3e  %-+3.3e  %-+3.3f\n",i,x(1),x(2),x(3),lambda)

end %for i

%////////////////////////////////////////

%////////////////////////////////////////

disp("-----------------");

disp("Método de Potencias  normalizando con máximo");

printf(" i      x1      x2    x3    lambda\n");

x=[1;1;1];

for i=1:12

    z=A*x;

    [maximo,p]=max(abs(z));

    m=z(p);

    lambda=z(1)/x(1);

    x=z/m;

    printf("%2d  %-+2.3f  %-+3.3f  %-+3.3f  %-+3.3f\n",i,x(1),x(2),x(3),lambda);

end %for i

%////////////////////////////////////////

%////////////////////////////////////////

disp("-----------------");

disp("Método de Potencias  normalizando con norm2");

printf(" i      x1      x2    x3    lambda\n");

x=[1;1;1];

for i=1:12

    z=A*x;

    m=norm(z);

    lambda=z(1)/x(1);

    x=z/m;

    printf("%2d  %-+2.3f  %-+3.3f  %-+3.3f  %-+3.3f\n",i,x(1),x(2),x(3),lambda);

end %for i

%////////////////////////////////////////

%////////////////////////////////////////

disp("-----------------");

disp("Método de Potencias norm2 y lambda=x'*z, multiplicando vectores");

printf(" i      x1      x2    x3    lambda\n");

x=[1;1;1];

for i=1:12

    z=A*x;

    x=z/norm(z);

    lambda=x'*z;

    printf("%2d  %-+2.3f  %-+3.3f  %-+3.3f  %-+3.3f\n",i,x(1),x(2),x(3),lambda);

end %for i

%-------------------

disp("-----------------");

disp(" Se observa que la convergencia se logra en los tres componentes del vector z./x")

z./x

%-------------------

disp("-----------------");

%Ejemplo de deflacion:

%http://riul.unanleon.edu.ni:8080/jspui/bitstream/123456789/5523/1/217575.pdf

disp("Deflacion");

v=z/norm(z)

v'*x

B=zeros(2,2)

B(1,:)=A(2,2:end)-v(2)/v(1)*A(1,2:end)

B(2,:)=A(3,2:end)-v(3)/v(1)*A(1,2:end)

disp('Valores propios de B')

eig(B)

disp('Valores propios de A')

eig(A)

disp("-----------------");

disp("Valores propios para un sistema 2x2 det(lambda*I-B)=0");

a=1

b=-(B(2,2)+B(1,1))

c=B(1,1)*B(2,2)-B(1,2)*B(2,1)

p=[a b c]

roots(p)

eig(B)