Metodi numerici per catene di Markov - laboratorio
Docente: Beatrice Meini

Lezione 2: Metodi iterativi per catene di Markov finite

La function [x,iter,vres]=potenze(P,eps,maxiter) calcola il vettore invariante di P mediante il metodo delle potenze:

  1. Costruire la matrice P con la function creaptrid.m, scegliendo ad esempio n=20, a=0.8, c=0.2 e applicare il metodo delle potenze con ad esempio eps=1.e-13, maxiter=100. Tracciare in scala semilogaritmica il vettore dei residui, con il comando semilogy( [1:iter+1], vres)

  2. Calcolare il secondo autovalore di modulo massimo di P: es=sort(abs(eig(P)),'descend'); r =es(2). Sulla stessa figura tracciare il grafico, in scala semilogaritmica, della funzione r^i, per i=1:iter, e dei residui: semilogy([1:iter+1], vres, [1:iter], r.^[1:iter]). Cosa si osserva?

  3. Calcolare l'errore assoluto e relativo, componente per componente, dell'approssimazione x ottenuta con il metodo delle potenze (ricordare che la function creap.m calcola il vettore invariante esatto). Cosa si osserva? Provare ad abbassare il valore eps.

  4. Ripetere gli stessi esperimenti applicando la permutazione della lezione precedente alla matrice P

La function [x,iter,vres]=rs(P,eps,maxiter) calcola il vettore invariante di P mediante un regular splitting, da scegliere (s)commentando parti della function:

  1. Confrontare la velocita' di convergenza del metodo di Gauss-Seidel con il metodo delle potenze: [xgs,itergs,vresgs]=rs(P,1.e-13,100); semilogy([1:length(vres)],vres,[1:length(vresgs)],vresgs)

  2. Modificare la function in modo che applichi il metodo di Jacobi. Cosa si osserva?

Fare gli stessi esperimenti con altre matrici stocastiche, ad esempio:

  1. P = [ 0.5 0.5 0 0; 0 0.8 0.2 0; 0 0 0.6 0.4; 0.9 0 0 0.1];

  2. P = [0 1 0 0 0; 0 0 0 0.5 0.5; 1 0 0 0 0; 0 0 0 0 1; 0 0 1 0 0];

Interpretare i risultati sperimentali calcolando gli autovalori della matrice di iterazione H.

Nel caso in cui un metodo iterativo non converge, usare la matrice di iterazione H1=(1-a)*eye(n)+a*H dove a e' compreso tra 0 e 1, e osservare cosa succede, per diversi valori di a.