Calcolo numerico/Metodi multipasso

Aumento dell'ordine di convergenza

Per aumentare l'ordine di convergenza del metodo ci sono due strade possibili:

1. metodi multipasso: usano informazioni precedenti ma rimangono lineari, ne sono esempi i metodi di Adam-Bashword (espliciti), di Adam-Multonn (impliciti), di Mistrom, i generalized Simpson (espliciti e impliciti), BDF usate nel caso dei problemi stiff;
2. metodi a un passo non lineari.

Metodi multipasso lineare

I metodi di questo tipo si definiscono come

$\sum _{j=0}^k{\alpha }_j{y}_{n+j}=h\sum _{j=0}^k{b}_j{f}_{n+j}$

dove

$f_{n+j}=f(t_{n+j},y_{n+j})$

${\alpha }_j,{\beta }_j$ sono tali che ${\alpha }_k=1$ e

$\sum |{\alpha }_0|+|{\beta }_0|\ne 0$

Osserviamo che se ${\beta }_k=0$ il metodo è esplicito perché a secondo membro non compare la soluzione da calcolare al passo corrente, mentre diversamente è implicito, e si può pensare di applicare il metodo del punto fisso

$y_{n+k}=h{\beta }_{k}f(t_{n+k},y_{n+k})+\varphi$

Dev'essere $h<1/(|{\beta }_k|\ast L)$.

Teorema 7.2

Consideriamo i seguenti polinomi

$\rho (z)=\sum _{j=0}^k{\alpha }_j{z}^j$

con $z$ variabile complessa (primo polinomio caratteristico).

$\sigma (z)=\sum _{j=0}^k{\beta }_j{z}^j$

e consideriamo il metodo

$\sum _{j=0}^k{\alpha }_j{y}_{n+j}=h{\varphi }_f(y_{n+k},y_{n+k-1},y_n,t_n,h)=h\sum _{j=0}^k{\beta }_j{f}_{n+j}$

Il metodo è consistente (cioè l'errore di troncamento locale tende a 0 per $h\to 0$), se valgono le seguenti condizioni:

$\sum _{j=0}^k{\alpha }_j=0,⟶\rho (1)=0$

$\sum _{j=0}^{k}j{\alpha }_{j}f(t_n,y(t_n))={\varphi }_f(t_n,t_n,t_n$

($h$ tende a 0 quindi tutti gli istanti tendono $t_n$)

con ${\varphi }_f=\sum _j{\beta }_j{f}_{n+j}$.

Dimostrazione

Affermiamo che

$y(t_{n+j})-y(t_n)=jh{y}^{\prime }({\xi }_j)=jhf({\xi }_j)$

con ${\xi }_j\in (t_n,t_{n+j})$ con $j=0,\dots ,k$ per il teorema del valor medio.

Consideriamo l'errore di troncamento:

$\tau =1/h[\sum _{j=0}^k{\alpha }_{j}y(t_{n+j})-h{\varphi }_f]$

e in base a quanto appena ricavato per $y_{n+j}$:

$=1/h[\sum _{j=0}^k{\alpha }_j\ast (y(t_n)+jh{y}^{\prime }({\xi }_j))]-{\varphi }_f$

e spezzando le somme

$=1/h\ast \sum _{j=0}^k{\alpha }_j\ast y(t_n)+\sum _{j=0}^{k}j{\alpha }_j{y}^{\prime }({\xi }_j)-{\varphi }_f$

Per $h\to 0$ $y^{\prime }({\xi }_j)\to y^{\prime }(t_n)$ e ${\varphi }_f(y^{\prime }(t),y^{″}(t),\dots )={\varphi }_f(t_n,t_n,t_n)$

Segue quindi che si ha consistenza se

$\sum _{j=1}^k{\alpha }_j=0$

e se

$\sum _{j=0}^{k}j{\alpha }_j{y}^{\prime }(t_n)-{\varphi }_f=0$

Nella pratica tenendo conto che $y^{\prime }(t_n)=f(t_n,y_n)$ si impongono le condizioni di consistenza:

$\sum _j{\alpha }_j=0$

$\sum _j[j{\alpha }_j-{\beta }_j]=0$

$\sum _j[j^2/2{\alpha }_j-j{\beta }_j]=0$

$\sum _j[j^3/6{\alpha }_j-j^2/2{\beta }_j]=0\dots$

Relazione tra convergenza e consistenza

Se ho convergenza, allora ho consistenza, ma non vale viceversa. Per avere convergenza, oltre alla consistenza dobbiamo richiedere la 0-stabilità.

Si dice che il problema è ben posto o totalmente stabile se per qualunque coppia di perturbazioni $(\delta (t),\delta )$ e $({\delta }^{\ast }(t),{\delta }^{\ast })$ tali che

$|\delta (t)-{\delta }^{\ast }(t)|<|\delta -{\delta }^{\ast }|$

segue che, posto $z(t),z^{\ast }(t)$ le soluzioni dei corrispondenti problemi di Cauchy, esiste una costante $s$ positiva tale che per ogni $t\in I_0$,

$|z(t)-z^{\ast }(t)|\le s\epsilon$

Definizione corrispondente nel discreto: consideriamo un sistema di differenze perturbato.

$\{\begin{array}{c}{z}_{n+1}=z_n+h{\varphi }_f(t_n,y_n,h)+{\delta }_{n+1})\\ z_0=y_0+{\delta }_0\end{array}$

Il metodo è 0-stabile se, date le successioni delle perturbazioni $\{{\delta }_n{\}}_{n=0}^N$, e $\{{\delta }_n^{\ast }{\}}_{n=0}^N$ che corrispondono a soluzioni $z_n$ per $n=0,\dots ,N$ e $z_n^{\ast }$ per $n=0,\dots ,N$, se esiste una costante $s$ positiva ed esiste $h_0>0$ tale che scelto $0<H<H_0$, segue che la norma della differenza tra le due soluzioni discrete è minore di $s\epsilon$, dove supponiamo che $|{\delta }_n-{\delta }_n^{\ast }|<\epsilon$ (la differenza rispetto alla definizione precedente è che deve esistere $H_0$ tale per cui si abbia la conservazione dell'entità della perturbazione).

Stabilità come condizione necessaria

Questo esempio mostra che è necessario richiedere la stabilità.

Esempio 7.1

Consideriamo il problema di Cauchy

$\{\begin{array}{c}{y}^{\prime }=0\\ y(t_0)=0\end{array}$

Applicando il metodo multipasso:

$\sum _{j=0}^k{\alpha }_j{y}_{n+j}=0$

perché $f$ (e quindi il secondo membro) è identicamente nullo.

$(y_1,\dots ,y_{k-1})$ sono molto piccoli ma non esattamente nulli.

Questa è un'equazione alle differenze omogenea. Il polinomio

$\rho (z)=\sum _{j=0}^k{\alpha }_j{z}^j$

ha radici $z_1,\dots ,z_s$ distinte rispettivamente con molteplicità ${\mu }_1,\dots ,{\mu }_s$, allora la soluzione dell'equazione alle differenze è

$y_n=\sum _{i=0}^s\sum _{j=0}^{{\mu }_{i-1}}{c}_{ij}{z}_i^n\ast n^j$

(simile alla scrittura delle soluzioni nelle equazioni differenziali nel continuo).

Se $|z_j|>1$, il metodo diverge, e lo stesso avviene nel caso in cui il modulo di $z$ è uguale a 1 e ha molteplicità maggiore di 1.

Root condition

Teorema 7.3

Il metodo numerico

$\sum _{j=0}^s{\alpha }_j{y}_{n+j}=h\sum _j{\beta }_j{f}_{n+j}$

è 0-stabile se e solo se vale la root condition.

Teorema 7.4

Il metodo

$\sum _{j=0}^s{\alpha }_j{y}_{n+j}=h\sum _j{\beta }_j{f}_{n+j}$

è convergente se e solo se è consistente e 0-stabile.

Teorema 7.5 (prima barriera di Dahlquist)

Nessun metodo lineare multipasso a $k$ passi che sia 0-stabile può avere ordine di convergenza maggiore di $k+1$ se $k$ è dispari, e di $k+2$ se $k$ è pari.

Espressione dei polinomi caratteristici nei metodi multipasso

Nei metodi di Adam-Bashword e Adam-Multonn il polinomio caratteristico è della forma

$\rho (z)=z^{k-1}\ast (z-1)$

e $z=1$ è la radice principale necessaria per la consistenza, e c'è poi la radice 0 con molteplicità $k-1$.

Nel Simpson generalizzato

$\rho (z)=z^k-z^{k-2}=z^{k-2}\ast (z^2-1)$

dove le radici sono $z_1=1,z_2=-1$ con molteplicità 1, e $z=0$ con molteplicità $k-2$.

Le BDF hanno il secondo polinomio caratteristico definito come

$\sigma (z)={\beta }^k\ast z^k$

Metodi di Adam-Bashword

I metodi di Adam-Bashword derivano da formule di quadratura. Consideriamo l'equazione integrale

$y(t_{n+1})-y(t_n)={\int }_{t_n}^{t_{n+1}}f(t,y(t))dt$

Consideriamo il polinomio interpolante $f$ con nodi equispaziati: in particolare, per i metodi espliciti consideriamo il polinomio interpolante di grado $p$ nei nodi $(t_n,\dots ,t_{n-p})$, per gli Adam-Multonn consideriamo il polinomio di grado $p+1$ nei nodi $t_{n+1},\dots t_{n-p}$ e otteniamo la formula

$y_{n+1}=y_n+h\sum _{j=-1}^p{\beta }_j{f}_{n-j}$

Se ${\beta }_{-1}=0$ il metodo è esplicito, altrimenti è implicito.

$\text{Metodi di Adam-Bashword:}{y}_{n+1}=y_n+h{f}_n\text{1 nodo, Eulero esplicito}$

$y_{n+1}=y_n+h/2(3f_n-f_{n-1})\text{2 nodi, trapezi esplicito}$

$y_{n+1}=y_n+h/12(23f_n-16f_{n-1}+5f_{n-2})\text{3 nodi}$

Adam-Multon

$y_{n+1}=y_n+h{f}_{n+1}\text{1 nodo, Eulero implicito}$

$y_{n+1}=y_n+h/2(f_{n+1}+f_n)\text{2 nodi, trapezi implicito}$

$y_{n+1}=y_n+h/12(5f_{n+1}+8f_n-f_{n-1})\text{3 nodi}$

In conclusione vale il seguente teorema.

Teorema 7.6

I metodi espliciti di Adam-Bashword hanno ordine $q$ mentre quelli di Adam-Multonn hanno ordine $p=q+1<math>aeccezionedel{l}^{\prime }Euleropercui<math>q=1$ e $p=1$./

Tecnica del predictor-corrector

Consideriamo un metodo multipasso implicito, da risolvere con un metodo di punto fisso. Poniamo

$z=y_{n+1}$

e consideriamo la relazione

$z=y_n+h{\beta }_{-1}f(t_{n+1},z)+h\sum _{j=0}^p{\beta }_j{f}_{n-j}$

dove il secondo addendo è noto e ne conosciamo tutte le quantità. $z_{k+1}=g(z_k)$ con $g=y_n+h{\beta }_{-1}f(t_{n+1},z)+h\sum _{j=0}^p{\beta }_j{f}_{n-j}$, e si dà come innesco $y_n$ soluzione all'istante precedente. Si impone $g^{\prime }<1$ quindi

$|g^{\prime }(z)|=h{\beta }_{-1}\ast |\frac{\partial f}{\partial z}|\le h{\beta }_{-1}L$

e supponendo di avere le condizioni del teorema di esistenza e unicità

$h{\beta }_{-1}L<1$

cioè sul passo $h$ otteniamo la condizione $h<|{\beta }_{-1}|\ast L$ e la costante di Lipschitz condiziona la convergenza del metodo di punto fisso.

Il metodo è costoso, perché a ogni iterazione del metodo di punto fisso si richiede una valutazione della funzione $f$. Più l'innesco fornito è buono, più velocemente viene soddisfatto il criterio d'arresto, e meno oneroso è il metodo dal punto di vista computazionale.

Metodo predictor-corrector

Il metodo predictor-corrector consiste dei seguenti passi:

1. predictor ($𝒫$): $y_{n+1}(x_0)$ è ottenuto applicando uno schema multipasso esplicito, cioè

$y_{n+1}(0)=y_n+h\sum _{j=0}^{\overline{p}}{\beta }_j{f}_{n-j}$

2. evaluation ($ℰ$):

$f_{n+1}=f(t_{n+1},y_{n1+}(0))$

(ora un'espressione approssimativa di $f_{n+1}$ è nota)

3. corrector ($𝒞$):

$y_{n+1}=y_n+h{\beta }_{-1}{f}_{n+1}+h\sum _{j=0}^p{\beta }_j{f}_{n-j}$

I passi $ℰ,𝒞$ possono essere iterati nel seguente modo:

$ℰf_{n+1}(k)=f(t_{n+1},y_{n+1}(k))$

$𝒞y_{n+1}(k)=y_n+h{\beta }_{-1}{f}_{n+1}(k)$

dove $f_{n+1}(k)$ è stato calcolato prima.

Si può quindi considerare anche il metodo:

$𝒫(ℰ𝒞)m$

in cui valutazione e correzione vengono alternate $m$ volte.

Il passo finale è

$f_{n+1}^m=f(t_{n+1},y_{n+1}^m)$

cioè il metodo si prepara per un'eventuale iterazione. Se $𝒫$ ha ordine $p-1$ i dati di innesco aggiuntivi sono calcolabili con un errore di ordine di $o(h^p)$ e se il corrector è di ordine $p$, allora $(𝒫𝒞)$ ha ordine $p$.

Esempio 7.2 (esempio di metodo predictor corrector)

Usiamo come predictor l'Eulero esplicito

$y_{n+1}(0)=y_n+h{f}_n$

e come corrector il metodo dei trapezi:

$y_{n+1}=y_n+h/2(f_n+f(t_{n+1},y_{n+1}(0)))$

Questo metodo può anche essere considerato come metodo a un passo.