Esame Calcolo numerico I

Teoria calcolo numerico

Lorenzo Quaranta
VR455934
18 maggio 2023

Condizionamento in generale

Per certi problemi è possibile che un piccolo errore nell’input risulta in un errore molto più grande nell’output. Quando questo succede si dice che il problema è mal condizionato. Un’altra caratteristica di un problema mal condizionato è che l’errore nella risposta è molto più grande del residuo. Quando gli errori non sono molto amplificati, cioè un piccolo errore nell’input risulta in un errore nell’output solo leggermente più grande, si dice che il problema è ben condizionato.

-r := f (^a) - f (a) Residuo dove ^a è il candidato per la soluzione, mentre a è la soluzione esatta; la quantità r dovrebbe essere 0 se ^a fosse la soluzione.

Dal Teorema della Media abbiamo che -f (^a) - f (a) = f (c) - ^a per un punto c fra a e ^a. Risulta che 1 - - ^a a = f (^a) - f (a) / f (c)

Ovvero errore = fattore di espansione * residuo. Se il fattore di espansione dell’errore è grande vuol dire che un piccolo residuo può risultare in un errore molto più grande. In questo caso il problema sarà mal condizionato.

Supponiamo che vogliamo calcolare la funzione y = F (x) per un certo input x = a (e output y) ma facciamo un piccolo errore nell’input cioè usiamo a + h invece di a. Secondo il Teorema della Media - F (a + h) - F (a) = F (c) (a + h - a) dove c è fra a e a+h. Se il fattore di espansione dell’errore F (c) è grande significa che un piccolo errore nell’input può essere amplificato in un errore molto più grande nell’output.

Metodo di bisezione

Il Teorema di Bolzano garantisce l’esistenza di uno zero: se f(x) è continua sull’intervallo [a, b] e f (a) e f (b) hanno segni diversi, ovvero f (a)f (b) < 0, allora esiste almeno un punto z (a, b) tale che f (z) = 0.

Input: la funzione f (x), un intervallo [a, b] tale che f (a)f (b) < 0 e una tolleranza tol per decidere quando fermarci. Siano a = a, b = b e n = 1.

• Mentre a > tol calcoliamo c = (a + b) / 2
• Se f (a_n)f (c_n) < 0 allora a_n+1 = a_n e b_n+1 = c_n.
• Se invece f (c_n)f (b_n) < 0 allora a_n+1 = c_n e b_n+1 = b_n.
• n = n + 1

Potremo ritornare c := (a + b) / 2 per la nostra stima. Avremo che z ∈ (a_n, b_n) e quindi l’errore sarà dimezzato per ogni iterazione. Grazie alla stima dell’errore possiamo analizzare quante iterazioni siano necessarie per avere un errore sotto una certa tolleranza. Dato ε l’errore che vogliamo raggiungere, cerchiamo n tale che (b - a) / (2ⁿ) ≤ ε ⇔ n + 1 ≥ (log(b - a) - log(ε)) / log(2).

Vantaggi:

• È molto semplice
• Dipende solo dal poter calcolare il segno
• Converge sempre se f (a)f (b) < 0

Svantaggi:

• Non trova zeri doppi
• Non trova zeri complessi
• Se confrontato ad altri metodi, è lento

Iterazioni di punto fisso

x = z è uno zero di f (x) se e solo se z ha la proprietà g(z) = z. Si dice che x = z è un punto fisso di g(x). Geometricamente, un punto fisso è l’intersezione di y = g(x) e y = x. Una volta messa nella forma di punto fisso possiamo provare a trovare uno zero: scegliamo x₀ e creiamo la successione x_n+1 := g(x_n) per n = 0, 1, 2...

Se l’iterazione converge, che tipo di punto è questo limite? Sia p = lim x_n che assumiamo esista. Allora se g(x) è continua p = lim x_n = lim x_n+1 = lim g(x_n) = g(lim x_n) = g(p). In altre parole, se un’iterazione di punto fisso converge, converge a un punto fisso della funzione di integrazione g(x). (NON converge sempre).

- Studiamo l’errore e_n := x_n - p dove x_n = p è il punto fisso.

• e_n+1 = x_n+1 - p = g(x_n) - p = g(x_n) - g(p)
• = g'(c_n)(x_n - p) Teorema della Media
• = g'(c_n) × e_n

In altre parole, l’errore è amplificato se |g'(c_n)| > 1 e diminuito se |g'(c_n)| < 1. c è fra x e p e allora, nella zona di p g'(c) ≈ g'(p). Infatti se lim x_n = p allora lim g'(c_n) = g'(p) e tutto dipende da |g'(p)|.

Se |g'(p)| < 1 diciamo che x = p è un punto di attrazione (convergenza), in caso contrario avremo un punto di repulsione (divergenza).

È possibile che |g'(p)| = 1, in questo caso potrebbe convergere o divergere; è richiesta un’analisi dettagliata caso per caso.

Confronto punto fisso e bisezione

Per bisezione una soprastima per l’errore è dimezzato ad ogni iterazione, mentre per punto fisso l’errore è moltiplicato del fattore |g'(z)|. Se quest’ultima è minore di 1/2 vuol dire che punto fisso è più veloce di bisezione, in caso contrario sarà meno veloce. C’è un metodo che è generalmente più veloce di tutti e due.

Metodo delle tangenti - Metodo di Newton

L’idea è molto geometrica, cominciamo con una stima iniziale dello zero x₀ e troviamo un miglioramento x₁ sostituendo la funzione dalla sua tangente a x = x₀. La formula è:

y = f(x₀) + f'(x₀)(x - x₀)

Il miglioramento x₁ è lo zero della tangente, da cui:

• 0 = f(x₀) + f'(x₀)(x - x₀)
• -f(x₀) = f'(x₀)(x - x₀)
• x = x₀ - f(x₀) / f'(x₀)

In generale avremmo quindi:

• x_n+1 = x_n - f(x_n) / f'(x_n) per n = 0, 1, 2...

Per capire la convergenza basta calcolare g'(z) dove x = z è uno zero di f(x), cioè, f(z) = 0. Assumiamo per ora che f''(z) = 0 non potendo dividere per 0.

Calcoliamo per g(x) = x - f(x) / f'(x):

• g(x) = x - f(x) / f'(x)
• g'(x) = 1 - f'(x)f''(x) / (f'(x))² = 0

Quindi il fattore di espansione dell’errore g'(c_n) ≈ g'(z) = 0 è molto piccolo e il più vicino è s a z. Per punto fisso abbiamo usato il Teorema della Media; notiamo che è un caso speciale di uno sviluppo di Taylor con errore:

• f(x) = f(a) + f'(a)(x - a) + f''(a)(x - a)²/2 + ...

Usando quello di grado 1 con a = z:

• e_n+1 = x_n+1 - z = g(x_n) - z = g(x_n) - g(z)
• = g'(z)(x_n - z) + g''(c_n)(x_n - z)²/2
• = g'(z)(x_n - z) + (x_n - z)²/2 g''(c_n)
• = e_n²/2 g''(c_n)

Nel caso che lim x_n = z, lim c_n = z. Quindi il comportamento dell’errore è come e_n+1 = Ce_n² convergenza quadratica. Un comportamento come quello di punto fisso è detto convergenza lineare, ovvero |C|e_n+1 = Ce_n con |C| < 1.

Attenzione:

• Se |C| > 1 si ha la divergenza anziché la convergenza
• Se |e_n| > 1 l’errore crescerà

Se invece f'(z) = 0 si può fattorizzare f(x) = (x - z)h(x).

Se anche f''(z) = 0 vuol dire che x = z è almeno uno zero doppio e f(x) = (x - z)²h(x).

In generale x = z è uno zero di ordine r se:

• f(x) = (x - z)^rh(x) con h(z) ≠ 0

Per la convergenza calcoliamo g'(z). Per r ≥ 1, |g'(z)| < 1 ⇔ r = 1. Quindi, per uno zero con molteplicità 2 c’è solo convergenza lineare. Ricordiamo che gli zeri doppi sono anche mal condizionati e Newton è meno veloce; se sappiamo la molteplicità di r allora la media delle r stime è una stima buona dello zero. Ho bisogno che x₀ sia vicino sufficientemente allo zero, cioè che x₀ sia nella zona di attrazione dello zero. Z = {x ∈ C : lim x_n = z} zona di attrazione.

C’è un algoritmo speciale per il caso p(x) = Σ_k=0ⁿa_kx^k. Per un polinomio l’iterazione è sempre:

• x_n+1 = x_n - p(x_n) / p'(x_n)

C’è un algoritmo efficiente per valutare p(x) a un punto generale x = z (Newton-Horner).

Condizionamento degli zeri di un polinomio

Per il polinomio p(x) = Σ_k=0ⁿa_kx^k, l’input è il vettore (a₀ ... a_n) e l’output uno degli zeri z = z_j.

Per analizzare il condizionamento di uno zero vogliamo calcolare il suo fattore di espansione dell’errore per un errore nel coefficiente dell’input a_k. Per gradi superiori al 5 non esiste una formula esplicita, possiamo calcolarne la derivata implicita:

∂z/∂a_k = -z^k/p'(z)

che è il fattore di espansione dell’errore nel coefficiente di input a_k. Quando uno di questi fattori è grande vuol dire che lo zero è mal condizionato.

Un difetto del Metodo di Newton è che è necessario sapere anche la derivata della funzione, che potrebbe non essere possibile calcolare, o comunque molto difficile. Il Metodo delle Secanti sostituisce la tangente dalla secante fra gli ultimi due punti.

Metodo delle secanti

Per trovare il nuovo punto x_n+1 usiamo i due precedenti x_n e x_n-1. La secante è la retta fra i due punti (x_n, f(x_n)) e (x_n-1, f(x_n-1)):

y = f(x_n-1) + m(x - x_n-1) dove m è il coefficiente angolare della retta.

m = (f(x_n) - f(x_n-1)) / (x_n - x_n-1)

Troviamo lo zero di questa retta:

x_n+1 = x_n - f(x_n) / m

che è proprio come Newton, ma con la derivata f'(x_n) sostituita dal coefficiente angolare m della secante.

Ordine di convergenza

Per punto fisso l’errore ha il comportamento:

|e_n+1| ≈ C |e_n|

Mentre per Newton:

|e_n+1| ≈ C |e_n|²

In generale, se l’errore ha il comportamento:

|e_n+1| ≈ C |e_n|^p

si dice che la convergenza è di ordine p.

Norme di matrici

Consideriamo il sistema lineare Ax = b dove A ∈ R^n×n e b ∈ Rⁿ. Sia ^x candidato per la soluzione. Calcoliamo il residuo r := b - A^x. Notiamo che, dato il candidato, il residuo si può sempre calcolare. Usiamo il residuo per capire quanto grande può essere l’errore e = x - ^x dove x è la soluzione esatta.

r = b - A^x = Ax - A^x = A(x - ^x)

e = x - ^x = A^-1r

Essendo un vettore, usiamo le norme:

• ||x||₂ norma euclidea
• ||x||₁ := Σ|x_i| norma uno
• ||x||_∞ := max |x_i| norma infinito

Tutte e tre fanno parte di una famiglia generale di norme:

||x||_p := (Σ|x_i|^p)^1/p dove p ∈ R, p ≥ 1.

Da notare che p = ∞ da ||x||_∞ = max |x_i| = lim ||x||_p quando p → ∞.

Quindi abbiamo che:

||x||₁ ≤ ||x||_p ≤ ||x||_∞ ≤ n^1/p ||x||_p

Visto che lim_p→∞ n^1/p = 1, allora:

||x||_∞ = max |x_i|.

Per determinare quanto più grande potrebbe essere ||y|| in confronto a ||x|| scriviamo:

||y|| / ||x||.