Iterazioni con il metodo delle potenze in Matlab

[Alessio_Ale]

Ciao a tutti, non riesco a risolvere l'esercizio di cui vi metto il testo.




parto dicendo che so come scrivere una fuction per usare il metodo delle potenze inverse

Codice:

function [lambda_p,w,m] = potenze_inverse(A,p,z,tol,m_max)

n=size(A);
w=z/norm(z);
lambda_p(1)=p;
[L,U,P]=lu(A-p*eye(n));

for m=1:m_max
    y=L\(P*w);
    z=U\y;
    lambda_p(m+1)=p+1/(w'*z);
    w=z/norm(z);
    if abs(lambda_p(m+1)-lambda_p(m)) <= tol*abs(lambda_p(m+1))
        break
    end
end
end


o per lo meno così mi è stato spiegato durante un' esercitazione in aula; però non riesco a trovare il modo di risolvere questo esercizio.
Tra l'altro, la domanda fa parte di un quiz a crocette a cui rispondere in pochi minuti quindi una soluzione così articolata mi pare faccia perdere troppo tempo.
Mi date una mano a risolverlo e magari a trovare un metodo più immediato?

Grazie a tutti!

[Alessio_Ale]

Nessuno che riesca ad aiutarmi?

[feddy]

Ciao Alessio_Ale.

Cosa è che ti manca per risolvere l'esercizio?

[Alessio_Ale]

In realtà tutto e niente
mi spiego meglio, quello che ho postato è il codice copiato dalla teoria, se ci metto dentro i miei dati però non ottengo il risultato che dovrei ottenere.
Poi come dicevo, dovendo dare una risposta in breve tempo nell' ambito di un test mi pare una soluzione troppo elaborata da riscrivere a memoria in un paio di minuti.
Quindi mi manca che non riesco a farlo funzionare e non riesco a semplificarlo.

[feddy]

Non ho controllato il codice. Ma assumendo sia corretto, dovrebbe darti una delle risposte elencate. Intendi dire che dopo aver calcolato l'errore relativo trovi un risultato che non è in elenco, giusto?

[Alessio_Ale]

Non ci arrivo nemmeno a scrivere l'errore, in realtà non ho idea di come calcolarlo l'errore.

[Alessio_Ale]

Ti copio function e script così come li ho scritti io:

La function

Codice:

function [lambda_p,w,m] = pot_inverse(A,p,z,tol,m_max)
n=size(A);
w=z/norm(z);
lambda_p(1)=p;
[L,U,P]=lu(A-p*eye(n));

for m=1:m_max
    y=L\(P*w);
    z=U\y;
    lambda_p(m+1)=p+1/(w'*z);
    w=z/norm(z);
    if abs(lambda_p(m+1)-lambda_p(m)) <= tol*abs(lambda_p(m+1))
        break
    end
end
end

e lo script

Codice:

m_max=6;
tol=1.0e-10;
z=ones(8,1);
p=1.5;
A=hilb(8);
[lambda_p,w1,m1] = pot_inverse(A,p,z,tol,m_max)

plot(1:m1+1,lambda_p,'b+-','linewidth',2)
lambda_max_p=lambda_p(end)
autovalori_A_eig=eig(A)

% ipotesi di calcolo dell' errore
err=(norm(autovalori_max_p)-norm(lambda_A_eig))/norm(autovalori_max_p)

[feddy]

Quanto fa

Codice:

abs(lambda_p(end)-autovalori_A_eig(end))/autovalori_A_eig(end)

?

[feddy]

L'errore lo devi calcolare sulla quantità che hai approssimato, cioè l'autovalore $\lambda_p$ più vicino a $p$. Nel tuo codice stai calcolando la norma di un vettore di cui in realtà ti interessa una sola componente, l'ultima. Cioè l'output del metodo delle potenze inverse. Altra cosa, nell'errore relativo devi dividere per il valore esatto, non per quello approssimato. Nota che

Codice:

autovalori_max_p

non è nemmeno definito nel tuo listato. Non è un problema, ma in genere è buona norma fornire un codice che compili/esegue senza errori