From a28aad0a4d9ab2da602113a02256b82ec00aac45 Mon Sep 17 00:00:00 2001 From: Leonardo Robol Date: Thu, 22 Oct 2009 10:58:09 +0200 Subject: [PATCH] Completata la lezione del 21 ottobre di Gemignani. --- capitolo2.tex | 82 ++++++++++++++++++++++++++++++++++++++++++++++++++++++----- 1 file changed, 75 insertions(+), 7 deletions(-) diff --git a/capitolo2.tex b/capitolo2.tex index 80d24ff..23c05e2 100644 --- a/capitolo2.tex +++ b/capitolo2.tex @@ -423,7 +423,8 @@ e calcoliamo $B = \herm Q A Q$ otteniamo \[ B = \herm Q A Q = \hat D + \beta_{\frac{n}{2}}z\herm z \] -con $z = \left[ \begin{array}{c} q_1 \\ q_2 \end{array}\right]$ dove $q_i$ è l'ultima colonna di $\herm Q_i$. +con $z = \left[ \begin{array}{c} q_1 \\ q_2 \end{array}\right]$ dove $q_i$ è l'ultima colonna di $\herm Q_i$ +per $i = 1,2$. \subsection{La ricerca degli autovalori} Abbiamo quindi ricondotto il problema a calcolare gli autovalori di una matrice $D + w\herm w$ dove $D$ @@ -438,10 +439,10 @@ del polinomio caratteristico vale solo per $\lambda \in \C \setminus (\spe{D_1} Notiamo inoltre che $( I + \theta(\lambda I - D)^{-1}z\herm z)$ è una matrice elementare e quindi ha determinante\footnote{% Ricordiamo che in genearale una matrice elementare $I - \sigma u \herm v$ ha determinante $1 - \sigma \herm v u$} $1 - \theta (\lambda I - D)^{-1} \herm zz$; sviluppando ulteriormente si ottiene -\[ +\begin{equation} \label{eqn:polmindivetimp} p(\lambda) = \prod_{j=1}^{n} (\lambda - \hat d_j) ( 1 + \theta (\lambda I - D)^{-1} \herm zz ) = \prod_{j=1}^{n} (\lambda - \hat d_j) ( 1 + \theta \sum_{k=1}^{n} \frac{z_k^2}{\lambda - \hat d_k} ) -\] +\end{equation} ed estendendo per continuità a tutto $\C$ si ottiene infine \[ p(\lambda) = \prod_{j=1}^n (\lambda - \hat d_j) + \theta\sum_{k=1}^n z_k \prod_{\substack{s=1 \\ s\neq k}}^n (\lambda - \hat d_s) @@ -458,9 +459,76 @@ Assumiamo dunque che \item $z_j \neq 0$ per ogni $j = 1 \ldots n$; \item $\hat d_k \neq \hat d_j$ per ogni $k \neq j$, ovvero che tutti gli autovalori siano distinti; \end{itemize} -A questo punto dalle relazioni precedenti ricaviamo che -\[ +A questo punto dalla relazione~\ref{eqn:polmindivetimp} ricaviamo che +\begin{equation} \label{eqn:eqsecdivetimp} p(\lambda) = 0 \iff 1 + \theta \sum_{j=1}^n \frac{z_j^2}{\lambda - \hat d_j} = 0 -\] +\end{equation} + e quindi gli autovalori delle matrici di ordine $\frac{n}{2}$ separano quelli della matrice grande, e si -può applicare (ad esempio) il metodo di bisezione. +può applicare (ad esempio) il metodo di bisezione, o in generale altri procedimenti di iterazione +funzionale. + +\begin{os} + Notiamo che per applicare la bisezione nella (\ref{eqn:eqsecdivetimp}) abbiamo bisogno di conoscere, oltre + ai $\hat d_j$ anche gli $z_j$ e ricordando che $z$ è costruito a partire dalle matrici $Q_1$ e $Q_2$, è chiaro + che dobbiamo conoscere anche queste due. Queste erano le matrici che diagonalizzavano $T_1$ e $T_2$ e si + ottenevano quindi \textbf{conoscendo gli autovettori} di $T_1$ e $T_2$. Dobbiamo allora, per poter procedere + con il metodo, calcolarci tutti gli autovettori. +\end{os} + +\subsection{Autovettori e problemi di stabilità} +Il fatto di dover calcolare gli autovettori ha un lato positivo ed uno negativo, precisamente +\begin{enumerate}[(a)] + \item Il calcolo effettivo risulta inaspettatamente semplice, tanto che, trovati gli autovalori saremo + capaci di trovare ogni autovettore con costo $O(n)$ e quindi di ottenerli tutti con costo $O(n^2)$; + \item Il problema del calcolo degli autovettori non è in generale ben condizionato e non si riescono a + dare condizioni perchè lo sia. Dato che influenzerà anche il calcolo degli autovalori nei passaggi successivi, + non possiamo garantire l'accuratezza del risultato finale; +\end{enumerate} + +L'idea per trovare gli autovettori è cercare quelli di $\hat D + \beta_{\frac n 2}z\herm z$ e poi una +volta ottenuta la matrice $U$ con questi, moltiplicarla per la matrice a blocchi con $Q_1$ e $Q_2$ sulla +diagonale (che d'ora in poi chiameremo $\hat Q$) per trovare gli autovettori della matrice $T$. \\ +Osserviamo che se $v$ è autovettore per la matrice $D + \beta_{\frac n 2}z\herm z$ allora +\[ + (\hat D + \theta z \herm z - \lambda I)v = 0 +\] +e quindi +\[ + (\hat D - \lambda I)v = -\theta z \herm z v +\] +Se fosse $\herm z v = 0$ avrei che $(\hat D - \lambda I)v = 0$ ma questo non è possibile perché $v$ +è non nullo e $\hat D - \lambda I$ è invertibile (come assunto precedentemente). Quindi deve essere +$\herm z v \neq 0$ e quindi dato che $v$ è determinato a meno di uno scalare possiamo assumere +$\herm z v = 1$. In definitiva si ottiene +\[ + v = - \theta z (\hat D - \lambda I)^{-1} +\] +che si può scrivere anche più esplicitamente come +\begin{equation} \label{eqn:calcoloautovettdivetimp} + v_j = \frac{- \theta z}{d_j - \lambda} +\end{equation} +Riassumendo, vogliamo trovare la matrice $U$ con i $v_j$ sulla diagonale e questo si può fare con +la Formula~\ref{eqn:calcoloautovettdivetimp} con costo $O(n^2)$. Ottenuta la matrice $U$ avremo +la matrice per diagonalizzare $T$ calcolando $\hat Q U$, e potremmo quindi procedere con il metodo. + +Sorgono però dei problemi computazionali non banali, più precisamente +\begin{enumerate}[(i)] + \item La matrice $U$ deve essere ortogonale ma questo non è garantito dall'algoritmo. Se ad un passo + abbiamo perdita di ortogonalità questo poi causerà un accumulo di errore in tutti i passi successivi + portando ad una imprecisione del calcolo. Sembra che in una implementazione del 2006 si sia riusciti + a risolvere questo problema. %% TODO: inserire una citazione + \item Come già sottolineato prima il calcolo degli autovalori potrebbe essere mal condizionato dal principio + portando ad un accumulo di errore notevole ad ogni passo. + \item Apparentemente il costo computazionale del calcolo di $\hat Q U$ non è trascurabile in quanto + un prodotto matrice per matrice in generale costa $O(n^3)$. Fortunatamente $U$ non è una matrice qualsiasi: + matrici del tipo $u_{ij} = \frac{r_i s_j}{x_i - y_j}$ sono dette matrici \emph{Cauchy-like} ed esistono + algoritmi per calcolare il prodotto fra una Cauchy-like ed una matrice generica con costo $O(n^2\log{n})$ +\end{enumerate} + +Come ultima nota, ricordiamo che questo algoritmo è applicabile solo alle matrici tridiagonali hermitiante +(e, previa tridagonalizzazione, a qualsiasi matrice hermitiana) ma non alle matrici in forma di Hessenberg +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 -- 2.1.4