Lezione di Gemignani del 23 ottobre e modifiche nell'introduzione

Leonardo Robol [2009-10-23 13:49]
Lezione di Gemignani del 23 ottobre e modifiche nell'introduzione
Filename
capitolo2.tex
introduzione.tex
diff --git a/capitolo2.tex b/capitolo2.tex
index 23c05e2..d31a7a3 100644
--- a/capitolo2.tex
+++ b/capitolo2.tex
@@ -531,4 +531,126 @@ Come ultima nota, ricordiamo che questo algoritmo è applicabile solo alle matri
 che non presentano purtroppo alcuna proprietà di separazione\footnote{ovvero quella che ci permette di trovare
 i limite su cui applicare la bisezione o qualsiasi altra iterazione funzionale}. \\
 Sono stati fatti molti tentativi a questo proposito, ma non esiste attualmente nessuna implementazione
-robusta del metodo Divide et Impera sulle matrici di Hessenberg.
\ No newline at end of file
+robusta del metodo Divide et Impera sulle matrici di Hessenberg.
+
+\section{Il metodo delle potenze}
+Per ultimo analizzeremo il metodo delle potenze che ci permette di trovare gli autovettori di una matrice
+conoscendo gli autovalori. In generale è particolarmente comodo per calcolare l'autovalore dominante.
+Consideriamo dunque questo caso per primo
+
+\subsection{Il metodo diretto}
+Supponiamo $A \in \mat{\C}{n}$ diagonalizzabile, per semplicità. Elimineremo questa ipotesi in seguito.
+Consideriamo gli autovalori ordinati in modo che
+\[
+ |\lambda_1| \geq |\lambda_2| \geq \ldots \geq |\lambda_n|
+\]
+e aggiungiamo l'ipotesi che $\lambda_1$ sia semplice e che sia in modulo strettamente maggiore degli altri.
+Il metodo si basa su questa semplice osservazione (la cui paternità viene, ancora una volta, attribuita a Gauss):
+Se prendiamo un vettore $y_0 \in \C^n$ e la nostra scelta non è incredibilmente sfortunata\footnote{%
+dove con incredibilmente sfortunata intendiamo che il prodotto scalare tra $y_0$ e l'autovettore relativo
+a $\lambda_1$ è esattamente uguale a $0$} possiamo considerare la successione
+\begin{equation}
+ \left\{ \begin{array}{ll}
+          y_0 & = y_0 \\
+          y_k & = Ay_{k-1}
+         \end{array} \right.
+\end{equation}
+Se scriviamo $y_0$ nella base di autovettori per $A$ abbiamo
+\begin{equation}
+ y_0 = \sum_{i=1}^{n} \alpha_i v_i
+\end{equation}
+e quindi osservando che $y_k = Ay_{k-1} = A^{k}y_0$ dopo $k$ iterazioni otteniamo
+\begin{equation} \label{eqn:metpotite}
+ y_k = A^ky_0 = \sum_{i=1}^n \alpha_i \lambda_i^k v_i =
+ \lambda_i^k \cdot \Big( \alpha_1 v_1 + \underbrace{\sum_{i=2}^n \alpha_i (\frac{\lambda_i}{\lambda_1})^k v_i}_{%
+\text{tende a}\ 0\ \text{se}\ k \to \infty} \Big)
+\end{equation}
+e ricordando l'ipotesi fatta precedentemente, ovvero che $|\lambda_1| > |\lambda_i|$ per ogni $i$ abbiamo
+che la seconda parte dell'uguaglianza va a $0$ per $k \to \infty$ e quindi la successione di vettori
+tende ad un multiplo di $v_1$, ovvero l'autovettore cercato.
+
+Purtroppo questa implementazione del metodo non è realistica, perché si incorrerebbe sicuramente nel caso di
+overflow o underflow in poche iterazioni, e probabilmente prima di avere una buona approssimazione dell'autovettore.
+Si ricorre quindi al semplice stratagemma di riscalare ad ogni iterazione il vettore ottenuto, ad esempio
+dividendolo per la sua norma infinito. In questo modo otteniamo una successione di vettori unitari che convergeranno
+all' autovettore unitario. \\
+Possiamo ora emprimere la successione in questo modo
+\[
+ \left\{ \begin{array}{ll}
+          y_0 & = y_0 \\
+	  y_k & = \displaystyle \frac{Ay_{k-1}}{||Ay_{k-1}||_\infty}
+         \end{array} \right.
+\]
+\begin{os}
+ Questo metodo mi permette anche di valutare e raffinare la stima dell'autovalore dominante. Se considero
+ infatti il rapporto fra le $j$-esime componenti di due vettori consecutivi, ovvero $\frac{y_k,j}{y_{k-1,j}}$
+ non è difficile verificare tramite la (\ref{eqn:metpotite}) che
+ converge all'autovalore cercato.
+\end{os}
+
+Esistono delle varianti del metodo delle potenze che sfruttano le idee qui esposte per risolvere problemi simili
+e migliorare le condizioni di convergenza.
+
+\subsection{Il metodo inverso con shift}
+Osserviamo come prima cosa che con pochi cambiamenti si potrebbe calcolare l'autovalore di modulo più piccolo
+considerando, invece che la matrice $A$, la matrice $A^{-1}$. \\
+Potremmo riscrivere l'iterazione in questo modo
+\begin{equation}
+ \left\{ \begin{array}{ll}
+          y_0 & = y_0 \\
+          y_k & = \beta_{k} \cdot A^{-1}y_{k-1}
+         \end{array} \right.
+\end{equation}
+dove $\beta_{k} = ||y_{k}||_\infty$\footnote{adotteremo questa notazione per tutta la sezione, d'ora in poi}.
+Questo procedimento ci imporrebbe però di calcolare l'inversa di $A$, un procedimento che in generale si tende ad
+evitare per problemi di instabilità. Possiamo quindi rileggere l'espressione precedente come risoluzione di
+un sistema lineare, ovvero
+\begin{equation}
+ \left\{ \begin{array}{lll}
+          y_0 &=& y_0 \\
+          Az_k &=& y_{k-1} \\
+          y_k &=& \beta_k \cdot z_k
+         \end{array} \right.
+\end{equation}
+A questo punto possiamo scegliere varie tecniche per risolvere il sistema lineare; sottolineamo che un metodo diretto
+che utilizza una fattorizzazione in questo caso è molto conveniente rispetto ad uno iterativo. Se calcoliamo
+la fattorizzazione di $A$, infatti, possiamo poi usarla per risolvere tutti i sistemi lineari della successione con
+un costo basso.
+
+In particolare, ricordiamo che per una tridiagonale hermitiana e per una matrice in forma di Hessenberg superiore
+il costo di una fattorizzazione QR sono rispettivamente $O(n)$ e $O(n^2)$, e sono anche il costo della risoluzione
+del sistema lineare. Si verifica sperimentalmente che, partendo da una buona approssimazione degli autovalori,
+il metodo converge in genere in 3-4 passi, indipendemente dalla dimensione di $A$. Si conclude quindi che il calcolo
+dell'autovettore dominante risulta essere $O(n)$ per le tridiagonali e $O(n^2)$ per le matrici in forma di
+Hessenberg.
+
+Ci poniamo ora il problema di come calcolare gli autovettori che stanno nella parte centrale dello spettro, ed è qui
+che interviene lo \emph{shifting}. Non è difficile osservare che gli autovalori della matrice $A - \gamma I$ sono gli
+stessi di $A$ ``shiftati'' di $\gamma$, ovvero sono
+\[
+ \spe{A - \gamma I} = \{ \lambda - \gamma \ | \ \lambda \in \spe{A} \ \}
+\]
+Supponiamo ora di avere una buona approssimazione $\hat \lambda_j$ dell'autovalore $\lambda_j$, dove con ``buona
+approssimazione'' intendiamo tale che
+\[
+ |\hat \lambda_j - \lambda_j| < |\hat \lambda_j - \lambda_k| \quad \forall k \in \spe{A}\setminus\{\lambda_j\}
+\]
+In questo caso $\hat \lambda_j - \lambda$ è l'autovalore di modulo più piccolo della matrice $A - \hat \lambda_j I$
+e il suo autovettore è quello relativo all'autovalore $\lambda_j$ della matrice $A$. Possiamo quindi applicare
+il metodo delle potenze inverso a $A - \lambda_j I$ e ottenere l'autovettore cercato. Questa operazione prende
+il nome di \emph{metodo delle potenze inverso con shift}.
+
+\begin{os}
+ In questo modo si può applicare un metodo qualsiasi per ottenere un'approssimazione degli autovalori (ad esempio
+ il metodo QR) e poi calcolare tutti gli autovettori della matrice con il metodo delle potenze inverso con shift.
+ Se la matrice è tridiagonale hermitiana si possono ottnere tutti con costo $O(n^2)$, mentre se è in forma di
+ Hessenberg superiore il costo diventa $O(n^3)$ per le considerazioni fatte in precedenza.
+\end{os}
+
+\begin{os}
+ Un'altra cosa che si può notare è che il metodo fornisce anche un possibile raffinamente per l'autovalore
+ $\lambda_j$ tramite il rapporto fra le componenti di $y_k$ e di $y_{k-1}$. Si potrebbe essere tentati
+ di sostituire questo eventuale raffinamento al posto dello shifting ma si verifica che questa non è necessariamente
+ una buona pratica. Esistono alcuni casi in cui può inficiare la convergenza del metodo.
+\end{os}
+
diff --git a/introduzione.tex b/introduzione.tex
index 8df7a72..d74cd15 100644
--- a/introduzione.tex
+++ b/introduzione.tex
@@ -27,6 +27,10 @@ Questi appunti si dividono in due parti
  Durante il corso del prof. Bini abbiamo usato il linguaggio Fortrano 95 e quindi gli esempi saranno
  esposti con quel linguaggio. Sto comunque valutando l'opportunità di esporre anche qualche algoritmo in
  qualche linguaggio diverso.
+ \item[Parte III -- Appendice] Tutte le cose che in realtà non fanno strettamente parte del corso ma
+ che potrebbero essere utili per comprendere meglio quello che si trova scritto in questi appunti. Ci
+ si trova ad esempio, una breve introduzione al linguaggio Fortran che è usato nella parte II per
+ implementare gli algoritmi.
 \end{description}

 \section*{Notazione}
ViewGit