Metodo di bisezione e iterazioni di punto fisso

Prerequisiti: Funzioni continue su un intervalloTeoremi sulle funzioni continue su un intervallo: zeri (con la bisezione), valori intermedi (f manda intervalli in intervalli), monotonia e continuità, continuità dell'inversa, Weierstrass (su [a, b] ci sono massimo e minimo).Funzioni continue su un intervallo → (teorema degli zeri), Teoremi di Rolle e di LagrangeRolle: f continua su [a,b], derivabile in ]a,b[, f(a) = f(b) ⇒ esiste c in ]a,b[ con f'(c) = 0. Lagrange (valor medio): stesse ipotesi senza f(a) = f(b) ⇒ esiste c con f'(c) = (f(b) − f(a))/(b − a): la tangente in c è parallela alla corda. Conseguenza: su un intervallo f è costante se e solo se f' = 0; su un insieme che non è un intervallo è solo costante a tratti (arctan x + arctan(1/x)).Teoremi di Rolle e di Lagrange → (valor medio), Condizionamento di un problema e stabilità di un algoritmoL'errore finale ha tre origini (dati, troncamento, arrotondamento). Un problema è ben posto se la soluzione esiste, è unica e dipende con continuità dai dati; il numero di condizionamento K = |x f'(x)/f(x)| misura di quanto un errore relativo sui dati si amplifica nel risultato (si perdono circa log10 K cifre). La stabilità riguarda invece l'algoritmo: un algoritmo instabile amplifica gli errori anche su un problema ben condizionato (ricorrenza I_n = 1 − n I_{n−1}); un problema mal condizionato non si salva con nessun algoritmo. Uno schema numerico converge se è consistente e stabile (derivata numerica: troncamento contro arrotondamento).Condizionamento di un problema e stabilità di un algoritmo → (condizionamento di una radice).

Si cerca uno zero $\alpha$ di una funzione continua $f: [a,b] \to \mathbb{R}$, cioè $f(\alpha) = 0$. Salvo casi particolari (polinomi di grado $\le 2$) non esiste una formula risolutiva (per $x^3 - x - 2 = 0$, $e^x = 3x$, $\cos x = x$ non si trova $\alpha$ con le funzioni elementari): si costruisce una successione $x_0, x_1, x_2, \dots$ che converge ad $\alpha$ e ci si ferma quando l'errore è sotto una tolleranza.

Il metodo di bisezione

Teorema (degli zeri). Se $f$ è continua in $[a,b]$ e $f(a)f(b) < 0$, esiste $\alpha \in\, ]a,b[$ con $f(\alpha) = 0$.

L'idea: se $f$ cambia segno su $[a,b]$ la radice sta nell'intervallo; guardando il segno di $f$ nel punto medio si capisce in quale metà sta, e si ripete.

Procedura (bisezione). Dato $[a_0, b_0] = [a,b]$ con $f(a_0)f(b_0) < 0$. Per $n = 0, 1, 2, \dots$:

1. $c_n = \dfrac{a_n + b_n}{2}$ (punto medio).
2. Se $f(c_n) = 0$ ci si ferma. Se $f(a_n)f(c_n) < 0$ si pone $[a_{n+1}, b_{n+1}] = [a_n, c_n]$, altrimenti $[a_{n+1}, b_{n+1}] = [c_n, b_n]$.

Ad ogni passo $f$ cambia ancora segno nel nuovo intervallo, quindi $\alpha \in [a_{n+1}, b_{n+1}]$.

Teorema (stima dell'errore). Il punto medio $c_n$ soddisfa $$|c_n - \alpha| \le \frac{b - a}{2^{\,n+1}}.$$

Dimostrazione. A ogni passo la lunghezza dell'intervallo si dimezza: $b_n - a_n = \frac{b-a}{2^n}$. Poiché $\alpha$ e $c_n$ stanno entrambi in $[a_n, b_n]$ e $c_n$ è il punto medio, $|c_n - \alpha| \le \frac{b_n - a_n}{2}$ ∎.

Quante iterazioni? Per avere $|c_n - \alpha| \le \text{tol}$ basta $\frac{b-a}{2^{n+1}} \le \text{tol}$, cioè $$n \ge \log_2\frac{b-a}{\text{tol}} - 1.$$ È una stima a priori: si sa dall'inizio quante iterazioni servono, indipendentemente da $f$. Per $b-a = 1$ servono $n \ge 8{,}97 \to 9$ per $\text{tol} = 10^{-3}$, $n \ge 18{,}9 \to 19$ per $10^{-6}$, $n \ge 32{,}2 \to 33$ per $10^{-10}$: ogni 3,3 iterazioni si guadagna una cifra decimale.

Esempio svolto. $f(x) = x^3 - x - 2$ su $[1,2]$: $f(1) = -2 < 0$ e $f(2) = 4 > 0$, quindi c'è una radice ($\alpha = 1{,}52138\ldots$).

Come si legge una riga: a $n = 1$, $c_0 = 1{,}5$ ha $f(c_0) = 3{,}375 - 1{,}5 - 2 = -0{,}125 < 0$, come $f(a_0)$; quindi la radice è in $[c_0, b_0] = [1{,}5, 2]$. A $n=2$, $c_1 = 1{,}75$ ha $f = 5{,}359 - 1{,}75 - 2 = +1{,}609 > 0$, come $f(b)$: la radice è in $[1{,}5, 1{,}75]$.

Osservazioni:

• L'errore non decresce monotonamente (a $n=1$ è $0{,}2286$ dopo $0{,}0214$): il maggiorante $\frac{b-a}{2^{n+1}}$ sì, l'errore effettivo no. Quindi bisezione è "garantita" ma non "uniformemente migliorante".
• Il metodo è robusto: servono solo la continuità e un cambio di segno, non la derivata. Converge sempre a una radice dell'intervallo. Non si può usare per le radici di molteplicità pari (ad esempio $x^2$ in $0$), dove $f$ tocca l'asse senza cambiare segno.
• È lento: la convergenza è lineare con fattore $\frac12$ (un bit binario per passo). Per questo la si usa per localizzare la radice e poi si affina con metodi più veloci (Metodo di Newton, metodo delle secanti e ordine di convergenzaOrdine di convergenza p: |e_{k+1}| ≈ C|e_k|^p (p = 1 lineare, p = 2 quadratica: le cifre corrette raddoppiano a ogni passo). Newton x_{k+1} = x_k − f(x_k)/f'(x_k): zero della tangente, ordine 2 per radici semplici con e_{k+1} ≈ f''(α)/(2f'(α)) e_k²; per radici di molteplicità m è solo lineare con tasso (m−1)/m (si recupera con x − m f/f'). Secanti: sostituisce f' con un rapporto incrementale, ordine (1+√5)/2 ≈ 1,618, una sola valutazione di f per passo. Newton richiede un x_0 abbastanza vicino: può ciclare o divergere.Metodo di Newton, metodo delle secanti e ordine di convergenza →).
• Se nell'intervallo ci sono più radici con cambio di segno, converge a una di esse (non si sa quale).

Iterazioni di punto fisso

Si riscrive $f(x) = 0$ come $x = g(x)$. Gli zeri di $f$ sono i punti fissi di $g$: $\alpha = g(\alpha)$. Si sceglie $x_0$ e si itera $$x_{k+1} = g(x_k), \qquad k = 0, 1, 2, \dots$$ Se $g$ è continua e $x_k \to \alpha$, allora passando al limite $\alpha = g(\alpha)$: il limite è un punto fisso.

Per una data $f$ ci sono molte riscritture. Per $x^3 - x - 2 = 0$ ad esempio $x = x^3 - 2$, oppure $x = \sqrt[3]{x+2}$, oppure $x = \frac{2}{x^2 - 1}$. Non sono equivalenti dal punto di vista numerico.

Esempio: $x = \cos x$. La radice di $\cos x - x = 0$ è $\alpha = 0{,}739085133\ldots$ (il punto in cui le curve $y = \cos x$ e $y = x$ si incontrano).

Si parte da $x_0 = 1$ (in radianti) e si itera $x_{k+1} = \cos x_k$. Errore $e_k = x_k - \alpha$:

Il rapporto tra errori successivi tende a un numero costante, $-0{,}6736 = g'(\alpha) = -\sin\alpha$: l'errore si moltiplica per $g'(\alpha)$ a ogni passo. Il segno alterno dice che i valori oscillano attorno ad $\alpha$ (si vede nel grafico a "ragnatela": $g'<0$).

Condizioni di convergenza

Teorema (delle contrazioni, o del punto fisso). Sia $g: [a,b] \to [a,b]$ derivabile con $$|g'(x)| \le L < 1 \quad \forall x \in [a,b].$$ Allora: (1) $g$ ha un solo punto fisso $\alpha \in [a,b]$; (2) per ogni $x_0 \in [a,b]$ la successione $x_{k+1} = g(x_k)$ converge ad $\alpha$; (3) valgono le stime $$|x_k - \alpha| \le L^k\,|x_0 - \alpha| \quad\text{(velocità)},\qquad |x_k - \alpha| \le \frac{L^k}{1-L}\,|x_1 - x_0| \quad\text{(a priori)},\qquad |x_k - \alpha| \le \frac{L}{1-L}\,|x_k - x_{k-1}| \quad\text{(a posteriori)}.$$

Dimostrazione (idea). Per il teorema di Lagrange (Teoremi di Rolle e di LagrangeRolle: f continua su [a,b], derivabile in ]a,b[, f(a) = f(b) ⇒ esiste c in ]a,b[ con f'(c) = 0. Lagrange (valor medio): stesse ipotesi senza f(a) = f(b) ⇒ esiste c con f'(c) = (f(b) − f(a))/(b − a): la tangente in c è parallela alla corda. Conseguenza: su un intervallo f è costante se e solo se f' = 0; su un insieme che non è un intervallo è solo costante a tratti (arctan x + arctan(1/x)).Teoremi di Rolle e di Lagrange →) per ogni $x,y \in [a,b]$ esiste $\xi$ tra loro con $|g(x) - g(y)| = |g'(\xi)||x-y| \le L|x-y|$: $g$ accorcia le distanze (contrazione).

• Esistenza: $h(x) = g(x) - x$ è continua con $h(a) = g(a) - a \ge 0$ e $h(b) = g(b) - b \le 0$ (perché $g(x) \in [a,b]$), quindi per gli zeri esiste $\alpha$ con $h(\alpha) = 0$.
• Unicità: se $\alpha, \beta$ sono punti fissi, $|\alpha - \beta| = |g(\alpha) - g(\beta)| \le L|\alpha - \beta|$, quindi $(1-L)|\alpha-\beta| \le 0$ e $\alpha = \beta$.
• Convergenza e stime: $|x_k - \alpha| = |g(x_{k-1}) - g(\alpha)| \le L|x_{k-1} - \alpha| \le \dots \le L^k|x_0 - \alpha|$, che tende a $0$ perché $L<1$. Inoltre $|x_{k+1} - x_k| \le L^k|x_1 - x_0|$ e, sommando la serie geometrica, $|x_{k+m} - x_k| \le \frac{L^k}{1-L}|x_1-x_0|$; per $m \to \infty$ si ha la stima a priori ∎.

Il caso locale: se $g$ è $C^1$, $\alpha = g(\alpha)$ e $|g'(\alpha)| < 1$, per continuità $|g'| \le L<1$ in un intorno di $\alpha$ e l'iterazione converge per $x_0$ abbastanza vicino ad $\alpha$. Se $|g'(\alpha)| > 1$ l'iterazione si allontana da $\alpha$ (il punto fisso è "repulsivo").

Velocità. Da $x_{k+1} - \alpha = g(x_k) - g(\alpha) = g'(\xi_k)(x_k - \alpha)$ si ha $\frac{e_{k+1}}{e_k} \to g'(\alpha)$: la convergenza è lineare con tasso $|g'(\alpha)|$. Quanto più $|g'(\alpha)|$ è piccolo, tanto più veloce. Per guadagnare una cifra decimale servono circa $\frac{1}{-\log_{10}|g'(\alpha)|}$ iterazioni: $\approx 5{,}8$ per $|g'| = 0{,}67$ ($\cos x$), $\approx 1{,}2$ per $|g'| = 0{,}144$.

Stesso problema, riscritture diverse

Esempio: $x^3 - x - 2 = 0$, $\alpha = 1{,}5213797$.

• $g_1(x) = \sqrt[3]{x + 2}$. $g_1'(x) = \frac{1}{3}(x+2)^{-2/3}$, e $g_1'(\alpha) = \frac{1}{3\cdot 3{,}52138^{2/3}} = 0{,}1440$. Su $[1,2]$ $g_1$ è crescente con $g_1(1) = 1{,}442$ e $g_1(2) = 1{,}587$, quindi $g_1([1,2]) \subset [1,2]$ e $|g_1'| \le \frac13\cdot 3^{-2/3} = 0{,}160$. Converge, con tasso $0{,}144$:

Ogni iterazione dà quasi una cifra decimale in più.

• $g_2(x) = x^3 - 2$. $g_2'(x) = 3x^2$, $g_2'(\alpha) = 6{,}94 > 1$: il punto fisso è repulsivo. Da $x_0 = 1{,}5$ si ottiene $1{,}375,\ 0{,}5996,\ -1{,}784,\ -7{,}68\ldots$: diverge.

La scelta di $g$ è parte del progetto del metodo. Una regola: se $g'(\alpha) = 0$ la convergenza è più che lineare (è il caso del metodo di Newton, nella prossima nota).

Criteri di arresto

Si ferma l'iterazione quando una di queste grandezze è sotto la tolleranza:

• Incremento $|x_k - x_{k-1}|$. Per una contrazione, l'errore è al più $\frac{L}{1-L}$ volte l'incremento: se $L = 0{,}9$ l'errore può essere $9$ volte l'incremento. Con $L$ vicino a $1$ un incremento piccolo non garantisce un errore piccolo.
• Residuo $|f(x_k)|$. Vale $|x_k - \alpha| \approx \frac{|f(x_k)|}{|f'(\alpha)|}$: se $|f'(\alpha)|$ è piccolo (radice quasi multipla, problema mal condizionato) un residuo piccolo non significa errore piccolo. Esempio: $f(x) = (x-1)^3$, $x = 1{,}01$: $f = 10^{-6}$ ma l'errore è $10^{-2}$.
• Numero massimo di iterazioni, sempre, per evitare cicli infiniti.

Codice

import math

def bisezione(f, a, b, tol):
    n = 0
    while (b - a) / 2 > tol:          # errore di c_n <= (b-a)/2^(n+1)
        c = (a + b) / 2
        if f(a) * f(c) <= 0: b = c
        else: a = c
        n += 1
    return (a + b) / 2, n

def punto_fisso(g, x0, tol, nmax=200):
    x = x0
    for k in range(1, nmax + 1):
        xn = g(x)
        if abs(xn - x) <= tol:
            return xn, k
        x = xn
    return x, nmax

print(bisezione(lambda x: x**3 - x - 2, 1, 2, 1e-6))
print(punto_fisso(math.cos, 1.0, 1e-6))
print(punto_fisso(lambda x: (x + 2)**(1/3), 1.5, 1e-10))

Output:

(1.5213804244995117, 19)
(0.7390855263619245, 34)
(1.5213797067927244, 11)

Bisezione: $19$ passi per $10^{-6}$ (previsti $19$). Punto fisso con $\cos$: $34$ iterazioni; con $\sqrt[3]{x+2}$: $11$ iterazioni per $10^{-10}$. In MATLAB fzero(@(x) x.^3 - x - 2, [1 2]) combina bisezione e metodi più veloci.

Errori tipici

• Applicare la bisezione senza controllare $f(a)f(b) < 0$, oppure su un'intervallo con radice di molteplicità pari.
• Dimenticare che $g$ deve mandare $[a,b]$ in sé: $|g'| < 1$ da sola, senza $g([a,b]) \subset [a,b]$, non basta per la convergenza globale.
• Contare le iterazioni della bisezione dall'errore: $c_n$ ha maggiorante $\frac{b-a}{2^{n+1}}$ (non $2^n$).
• Credere che $|g'(\alpha)| < 1$ garantisca la convergenza da ogni $x_0$: è una condizione locale, vale per $x_0$ abbastanza vicino ad $\alpha$.

Domande d'esame

1. Quante iterazioni di bisezione servono per ridurre l'errore sotto $10^{-8}$ partendo da un intervallo di lunghezza $2$? Risposta: $n \ge \log_2(2\cdot 10^8) - 1 = 26{,}6 \to 27$.
2. Enunciare il teorema delle contrazioni. Risposta: se $g([a,b]) \subset [a,b]$ e $|g'| \le L<1$, il punto fisso è unico e $x_{k+1} = g(x_k)$ converge per ogni $x_0$ con $|x_k - \alpha| \le \frac{L^k}{1-L}|x_1 - x_0|$.
3. Che ordine ha la convergenza di un'iterazione di punto fisso con $0<|g'(\alpha)|<1$? Risposta: lineare, con tasso $|g'(\alpha)|$.
4. Perché il residuo piccolo non basta come criterio di arresto? Risposta: l'errore è circa $|f(x_k)|/|f'(\alpha)|$, grande se $f'(\alpha)$ è piccola.