Presentazione Multigrid

Metodi Multigrid e analisi di metodi iterativi (Jacobi e Gauss-Seidel)

Margherita Maria M63001118

Martina Russo M63001128

Michelle Pepe M63001196

Indice

• Generalità
• Metodi iterativi
• Teorema della convergenza
• Jacobi e Gauss-Seidel: confronto
• Errore e residuo
• Elementi Multigrid
• Multigrid ed errore
• Analisi Multigrid
• Aliasing
• Applicazioni del metodo Multigrid
• Modi di Fourier
• Primo schema
• Schemi Multigrid
• Secondo schema
• Terzo schema
• Full grid

Metodi iterativi

• Un metodo iterativo si distingue da metodi esatti, in quanto non arriva a risolvere perfettamente un problema, ma si basa sulla convergenza verso la soluzione di quest’ultimo.

Forma del problema lineare

Solitamente un problema lineare si presenta nella seguente forma:

Ax=b

Ax=b

Da questa, possiamo ricavare la soluzione come:

x=(I-A)x+b

x=(I-A)x+b

Attraverso l’utilizzo dei metodi iterativi, è possibile esprimere la soluzione dello stesso problema in relazione al passo di iterazione k-esimo:

x_k+1=(I-A)x_k+b

x_k+1=(I-A)x_k+b

Precondizionamento

Ovviamente in una situazione del genere nulla si può dire sulla convergenza del metodo!

Tuttavia sappiamo che per ogni successiva iterazione, una nuova configurazione x_k+1 viene calcolata a partire da x_k.

Consideriamo il problema scritto a questo punto sotto questa forma:

Px_k+1=(P - A)x_k+b

Px_k+1=(P - A)x_k+b

La matrice P introdotta risulta essere molto utile in un metodo iterativo e prende il nome di matrice di precondizionamento.

Tra le possibili scelte, abbiamo: matrice di Jacobi, matrice di Gauss-Seidel.

Errore nei metodi iterativi

Un metodo iterativo converge se l’errore passo dopo passo diminuisce.

e_k=x_k-x

e_k=x_k-x

Sfruttando la matrice di precondizionamento, introdotta in precedenza, possiamo descrivere l’errore al passo successivo come segue:

Pe_k+1=Px_k+1 - Px = (P-A)x_k + b - (P-A)x - b = (P-A)e_k

Moltiplicando ambo i membri per la matrice P^-1 otteniamo:

e_k+1 = (I - P^-1A)e_k

Tale matrice la definiamo matrice M: l’obiettivo è far annullare la matrice per abbattere l’errore.

e_k+1 = Me_k

Metodo iterativo: teorema della convergenza dell’errore

Per poter annullare l’errore, faremo riferimento agli autovalori di M.

Qualora l’errore e_k fosse un autovettore dello spazio definito da M, varrebbe la relazione:

Me_k = λ_ke_k

Altrimenti l’errore potrebbe essere comunque scritto come combinazione lineare degli autovettori di M che costituiscono una base dello spazio di M.

Me_k = c₁λ₁v₁ + c₂λ₂v₂ + … + c_nλ_nv_n

Me_k = c₁λ₁v₁ + c₂λ₂v₂ + … + c_nλ_nv_n

Possiamo notare come l’abbattimento dell’errore dipenda dal valore degli autovalori di M.

Nell’ipotesi in cui valga il teorema della convergenza tale per cui |λ_j(M)|< 1 per ogni j, la velocità di convergenza è tanto più rapida, tanto minore è il più grande raggio spettrale degli autovalori di M; tale valore prende il nome di raggio spettrale.

ρ(M) = max {|λ_j(M)|}

Metodi iterativi: Jacobi e Gauss-Seidel

Tra i metodi iterativi vogliamo porre attenzione sui metodi di Jacobi e Gauss-Seidel.

Jacobi e Gauss-Seidel generano una successione di valori che converge alla soluzione del sistema:

==

Partiamo dal sistema lineare Ax=b con P non singolare.

A=P-N

A=P-N

Ax=b Px=Nx+b

Ax=b Px=Nx+b

In generale, con un metodo iterativo, abbiamo che:

�+1 � �^{-1 -1}� � � �=��=� +� +� con k ≥ 0

b

Jacobi

Jacobi A=L+D+U

A=L+D+U

P=D

P=D

A= N=-(L+U)=D-A

N=-(L+U)=D-A

con k ≥ 0

Jacobi

^{-1 -1 -} 1( )=� ( )=�� �+� �− � �=−� −�^{-1 -1 -} 1��

Con matrice diagonale con elementi 1/�

+1 � ^-1� � �=� +��

+1 � ^-1��

In forma scalare

Jacobi - analisi

Il metodo di Jacobi richiede di tenere in memoria almeno due vettori (2n), ma i calcoli possono essere svolti in parallelo sulle componenti.

Per ogni iterazione il calcolo di ognuna delle n componenti richiede n-1 moltiplicazioni e n-1 somme.

Di conseguenza ricaviamo una complessità computazionale:

2� (� )

Esempio

Presentiamo un esempio relativo alla risoluzione con metodo iterativo di un sistema del tipo Ax=b, con matrice A tridiagonale.

• Questa matrice presenta autovalori del tipo: λ_j(A) = 2 – 2cos().
• Siamo tuttavia interessati agli autovalori della matrice: M = (I-P^-1A).
• P matrice di Jacobi: P = diag(A) = 2I.

Possiamo quindi calcolare gli autovalori di M, come:

λ_j(M) = 1 – λ_j(A) = cos().

Dove con j = 1 otteniamo il raggio spettrale, cioè il più grande tra gli autovalori.

Esempio: autovettore associato

L’autovettore associato a j = 1 è una funzione trigonometrica la cui frequenza è molto bassa.

Tale funzione trigonometrica è: sin( kπ N + 1 ) k = 1, 2, . . . , N.

Dunque all’autovalore dal valore più alto è associato un autovettore derivato dal campionamento di un segnale a bassa frequenza.

Dato che l’autovalore nell’esempio appena presentato è associato ad un angolo prossimo a zero, vale la approssimazione di Taylor:

cos() () 2

Questo risultato non deve sorprendere; i metodi iterativi riescono con buona efficienza ad eliminare componenti in alta frequenza, ma non quelle in bassa frequenza.

Esempio in Matlab

Nel seguente codice Matlab è stato implementato l’esempio appena descritto per mostrare come avviene la convergenza verso la soluzione del sistema lineare attraverso Jacobi:

cos() () 2

Esempio