Appunti Fluidodinamica computazionale

La fattorizzazione LU

Si potrebbe pensare di usare la fattorizzazione LU in qualsiasi sistema, ma in realtà non funzionerebbe in tutti i casi. Facciamo un esempio: supponiamo di avere una matrice dei coefficienti pentadiagonale. Si potrebbe pensare di scrivere questa matrice come prodotto di una matrice triangolare bassa con le diagonali posizionate esattamente dove sono le diagonali inferiori della matrice di partenza, per una matrice triangolare alta con le diagonali posizionate esattamente dove sono nella matrice di partenza. In generale facendo il prodotto delle matrici non otteniamo la matrice dei coefficienti [A]. In casi specifici questo può essere fatto, trovando una maniera per scrivere il sistema [L][U] ma non c'è una regola generale. Facendo il prodotto delle matrice otterrei una matrice con più diagonali. Alcuni metodi di fattorizzazione sono costruiti in questa maniera, operando quella che viene chiamata approssimata ovvero ipotizzando che la matrice triangolare.bassa e alta hanno la stessa struttura della matrice di partenza. Questa è una maniera per descrivere un metodo iterativo.

Conviene comunque usare questa rappresentazione perché a ritroso posso trovare le incognite, il problema è trovare la fattorizzazione.

Metodi iterativi

In generale il metodo iterativo per la soluzione di un sistema lineare può essere scritto nel seguente modo:

Per risolvere il sistema devo individuare una procedura che a partire da una soluzione di primo tentativo, attraverso una serie di passaggi, porta alla definizione di iterazioni successive fino ad avere una soluzione che, sostituita all'interno del sistema, lo soddisfa identicamente.

Affinché si possa avere la soluzione, così come nel caso precedente, la matrice dei coefficienti [A] deve essere invertibile ovvero il suo determinante deve essere diverso da zero.

Splitting additivo

Una maniera per risolvere con un metodo iterativo è usare lo splitting additivo ovvero scrivere la matrice

[A] come differenza di due matrici: Le matrici [P] e [N] devono essere matrici dello stesso ordine rispetto ad [A], inoltre la matrice di precondizionamento [P] deve essere una matrice non singolare, ovvero deve essere invertibile.

Posso usare la (9.3) come formula ricorsiva per la valutazione dell'incognita.

Per adoperare lo splitting additivo ci sono vari criteri:

• l'inversione della matrice [P] deve essere più semplice dell'inversione della matrice [A] altrimenti il costo computazionale di una operazione diventerebbe confrontabile con quella del metodo diretto e non avrebbe più senso usare questo metodo;
• La matrice [P] deve essere sufficientemente simile alla matrice dei coefficienti [A] in modo che ogni singola iterazione riduca in maniera efficace l'errore.

L'errore viene definito nella seguente maniera: matrice di iterazione [M].

Posso comprendere le proprietà del metodo iterativo attraverso una.

Per capire questa matrice.

combiniamo le equazioni (9.3) e (9.4) ottenendo: La modalità con cui definisco lo splitting additivo definisce come evolve l'errore da una iterazione all'altra: Affinché il metodo sia convergente l'errore deve tendere a zero per iterazioni successive ovvero: Posso scrivere in maniera del tutto equivalente: Abbiamo visto un caso analogo quando abbiamo fatto l'analisi di stabilità, per questo possiamo dire che il raggio spettrale che affinché si abbia convergenza, il raggio spettrale della matrice [M], ovvero il massimo dei suoi autovalori deve essere minore di 1. Questa condizione non è sempre facile valutarla, ma possiamo usare anche condizioni equivalenti che riguardano una opportuna norma della matrice [M]. Di seguito sono riportate alcune condizioni equivalenti: Metodo di Jacobi (per punti) In questo metodo si utilizza uno splitting particolare. La matrice [A] viene divisa in più parti: Nel sistema originario consideravo la matrice b costituita soltantodai termine della diagonale, cioè valutavo all'interazione k+1 solo le incognite corrispondenti ai valori sulla diagonale della matrice [A]. Se ci vogliamo riferire al caso generico possiamo dire che: Per le nostre applicazioni, supponendo di voler risolvere un problema di Laplace 2D utilizzando una formula a 5 punti posso scrivere la seguente discretizzazione, già vista sopra: Possiamo applicare il metodo di Jacobi alla formula valutando l'incognita u in (i,j) all'interazione k+1, mentre le altre incognite devono essere valutate all'interazione precedente k. Possiamo anche vedere la scrittura del codice di questo metodo, che è molto semplice: dopo aver fornito la soluzione iniziale u0, devo memorizzare l'interazione k, valutare u(i,j) utilizzando l'espressione sopra, riversare su u(k) la soluzione ottenuta, fare un check sulla convergenza (calcolando la forma dell'equazione discreta in ciascun punto e calcolare il residuo) e rifare.le stesse operazioni sarebbero state eseguite, ma in un ordine diverso. Il metodo di Gauss-Seidel per punti permette quindi di ottimizzare il calcolo delle incognite, sfruttando i valori già calcolati nel processo iterativo.le variabili appena valutate saranno diverse e dipenderanno dall'ordine con cui sto aggiornando il processo iterativo nel vettore delle incognite. In questo senso ha un significato scrivere la formula di Gauss-Seidell in maniera esplicita: Questo metodo è sicuramente più efficiente del precedente perché non devo aspettare il completamento dell'aggiornamento delle componenti del vettore per poter utilizzare i valori all'interazione k+1. Dal punto di vista implementativo questo comporta che non ho bisogno di memorizzare il vettore delle incognite all'interazione k e differenziarlo da quello all'interazione k+1, posso semplicemente aggiornarli: non ho quindi la necessità di memorizzare due vettori, me ne basta uno. Possiamo vedere il codice:

La verifica sulla convergenza in questo caso è relativamente semplice: i metodi (sia Gauss-Seindel che dominanza diagonale: Jacobi per punti) convergono se la matrice è a

Sicuri che entrambi i metodi sono convergenti, ma poiché Gauss-Seindel utilizza subito i valori aggiornati, questo metodo ha bisogno di meno iterazioni per giungere a convergenza. Poiché abbiamo utilizzato parzialmente la matrice dei coefficienti, per questo metodo lo splitting additivo è cambiato: Metodo di sovrarilassamento metodo.

I metodi precedentemente visti possono essere modificati utilizzando una tecnica nota come di sovrarilassamento. Supponiamo di utilizzare un metodo iterativo qualsiasi. A valle otteniamo la soluzione all'interazione k+1: partendo dalla soluzione all'interazione k sono arrivato alla soluzione all'interazione successiva attraverso un delta. La correzione Δx posso pensare di aumentarla per raggiungere più velocemente la soluzione stessa, questo è il senso di rilassamento. Indichiamo con x la soluzione del generico metodo iterativo, posso pensare di scrivere la soluzione alla nuova iterazione come un coefficiente.

ω che moltiplica la soluzione fornita dal metodo iterativo più (1-ω) per la soluzione all'interazione k:

Se ω>1 sto correggendo l'interazione precedente attraverso una correzione che è una correzione maggiorata, nell'ottica di arrivare più velocemente alla soluzione stessa. In questo caso si parla di sovrarilassamento.

Il coefficiente di sovrarilassamento ω lo posso utilizzare grande purché preservi la convergenza del metodo. Devo scegliere ω in modo da preservare la convergenza, quindi preservare la dominanza diagonale.

Il coefficiente può essere anche minore di 1 se osservo che il metodo non sia convergente. In questo caso si parla di sottorilassamento.

Con la condizione di dominanza diagonale si può verificare che il coefficiente ω deve essere compreso tra 0 e 2 se utilizziamo il metodo di Gauss-Seindel.

Può essere necessario usare il sottorilassamento quando arrivo alla soluzione del

Partendo da un sistema non lineare che ho opportunamente linearizzato e che risolvo attraverso questa tecnica. La linearizzazione introduce una ulteriore variabile, in termini di convergenza, rispetto a quello analizzato sopra.

Quando il processo iterativo parte da un sistema non lineare, arriva ad un sistema lineare attraverso una linearizzazione, e questo sistema lineare lo risolvo con un metodo di rilassamento. È possibile che abbia bisogno di un sotto rilassamento per evitare che il metodo diverga.

Per il metodo di Jacobi invece non posso sovrarilassare.

Metodi iterativi a blocchi:

I metodi iterativi visti finora non sono i più efficienti, quelli più efficienti sono i metodi a blocchi.

Questo metodo non considera singolarmente le equazioni, ma, nell'ottica di incrementare l'efficienza del metodo e quindi ridurre il numero di iterazioni, le raggruppa e risolve una serie di sistemi lineari di dimensioni più piccole.

Raggruppando le

equazioni e diminuendo le dimensioni del problema, sarà più facile invertire le matrici dei coefficienti. Avrò bisogno di più operazioni per completare un'interazione, perché sarò costretto a invertire sistemi più piccoli, però avrò un vantaggio nel numero di iterazioni.

Dobbiamo capire come raggruppare in maniera intelligente le incognite per ottenere l'obiettivo.

Nel caso delle griglie cartesiane, il raggruppamento delle incognite è semplice, posso infatti pensare di raggruppare le incognite per righe (ad eccezione dalla prima e ultima dove sono applicate le condizioni al contorno).

Posso pensare di scrivere le equazioni discrete e contemporaneamente risolvere le incognite che sono su una stessa riga. Raggruppando in questa maniera le incognite e adottando una tecnica di Jacobi, l'equazione di Laplace può essere scritta nella seguente maniera:

Raggruppando per righe, le incognite

All'interazione k+1 sono quindi tutte quelle con j fisso. Il vantaggio di aver raggruppato le incognite in questa maniera.