Giugno 2025
Spiegazione dei principali metodi di ricerca di zeri di funzioni reali. Definizioni preliminari, e introduzione ai metodi di bisezione, Newton e quasi-Newton.
Consideriamo una funzione \(f\) reale. Un metodo iterativo è un algoritmo scritto con lo scopo della costruzione di una successione che converge allo zero cercato della funzione \(f\).
Per definizione di algoritmo, questo deve arrestarsi secondo un criterio di arresto semplice da verificare, non costoso ed efficace. Il criterio di arresto si basa su diversi fattori:
tolleranza, cioè una quantità che in modo assoluto o generale controlla l’errore commesso,
numero di iterazioni svolte,
costo della singola iterazione del metodo,
velocità di convergenza.
Cerchiamo allora di dare una definizione di velocità di convergenza. Sia \((x_k)_k\) la successione convergente a \(\hat x\) generata dal metodo iterativo in questione. Definiamo la successione \(e_k= |x_k - \hat x|\). Allora, poiché il metodo è convergente, la successione \(e_k\) tende a zero.
Diciamo che \((x_k)\) converge con ordine di convergenza \(p>1\) se
\[\lim_{k\rightarrow\infty} \frac{e_{k+1}}{(e_k)^p}=c\in (0, +\infty)\]
cioè se \[|x_{n+1}-\hat x| = O(|x_n-\hat x|^{p})\quad n\rightarrow \infty\].
Se invece \(p=1\), la convergenza si dice lineare, e si richiede che la costante asintotica sia \(c \in (0, 1)\), per il criterio del rapporto. Se così non fosse, infatti, la successione non sarebbe convergente. Si dice che la convergenza è quadratica quando \(p = 2\), e così via. La successione si dice a convergenza sublineare se \(c = 1 = p\): è il caso peggiore.
Notiamo che tale numero \(p\) non può essere minore di uno, per il semplice fatto che la successione del rapporto sarebbe divergente altrimenti.
Si parla di metodo di iterazione funzionale quando la successione costruita dal metodo è del tipo
\[x_{k+1} = \Phi(x_k)\]
Richiamiamo alcune definizioni dall’analisi matematica:
Si dice che \(x_0\) è punto fisso per la funzione \(\Phi\) se \(\Phi(x_0) = x_0\).
Diciamo che \(\Phi\) è una contrazione se \[\Phi:D\rightarrow D \subseteq\mathbb{R}\] e se \[\exists L\in[0, 1]: |\Phi(x) - \Phi(y)| \leq L|x-y| \quad\forall x,y\in D\]
Il lemma di contrazione, di cui ometto la dimostrazione (sebbene parte del programma di Analisi Numerica I, 2025, ma facilmente reperibile negli appunti del corso di Analisi I o nelle slides caricate dalla professoressa) assicura che
Se \(x_0\in D\) è un punto fisso per \(\Phi\), allora è il suo unico punto fisso.
La successione \((x_k)_k\) definita da \(x_{k+1} = \Phi(x_k)\) è tale che \(\{x_k\} \subseteq D\) e \((x_k)_k \rightarrow x_0\) punto fisso.
Questo risultato è utile per dimostrare il seguente teorema.
Teorema 1 (Ordine di convergenza per metodi a iterazione funzionale). Sia \(\Phi : D\rightarrow D\) contrazione su \(D\) connesso e aperto, sia \(\hat x \in D\) il suo punto fisso, e sia \(L\) la costante di contrazione di \(\Phi\). Se \(\Phi \in C^p(D)\) e \(p\in\mathbb{N}\) è tale che
\[\forall k = 1\dots (p-1) \quad \Phi^{(k)}(\hat x) = 0 \text{ e } \Phi^{(p)}\neq 0\]
(cioè \(\hat x\) è uno zero multiplo di molteplicità \(p\) per \(\Phi\)), allora la successione \((x_k)_k\) definita da \[x_{k+1} = \Phi(x_k)\]
ha ordine di convergenza \(p\).
Infatti, per la formula di Taylor con resto in forma di Lagrange per funzioni di una variabile, senza perdita di generalità, esiste un \(\eta_k \in (x_k, \hat x)\) tale che
\[x_{k+1} = \Phi(x_k) = \Phi(\hat x) + \frac{\Phi^{(p)}(\eta_k)(x_k - \hat x)^p}{p!}\]
Adottando la notazione dell’errore con \(e_k\) vista prima, otteniamo
\[\lim_k\frac{e_{k+1}}{e_k^p} = \frac{\Phi^{(p)}(\hat x)}{p!}\] ◻
Siamo interessati a cercare un modo di velocizzare un metodo iterativo a convergenza lineare. Supponiamo allora che \((x_k)_k\) sia una successione con convergenza lineare, cioè tale che \(x_{k+1} = O(x_k)\). Allora possiamo approssimare
\[\frac{x_{k+1}- \hat x}{x_k - \hat x} = \frac{x_k - \hat x}{x_{k-1}- \hat x}\]
e ottenere
\[\hat x= \frac{x_{k-1}x_{k+1} - x_k^2}{x_{k+1}-2x_k +x_{k-1}}\]
Supponiamo adesso che \(x_{k+1} = f(x_k)\). La tecnica di accelerazione di Aitken è un metodo iterativo che fornisce una successione data da
\[\begin{aligned} z_0 &= x_0 \\ \Phi_{A}(z)& = \frac{zf(f(z)) - f^2(z)}{f(f(z)) - 2f(z) + z} \\ \Phi_A(\hat x) &= \hat x \end{aligned}\]
E’ necessario definire \(\Phi_A(\hat x)\) poiché \(\hat x\) è l’unico punto in cui la successione di Aitken si annulla. Questa, insieme ad altre proprietà sono riassunte nel seguente teorema.
Teorema 2 (Aitken). La successione \(\Phi_A\) definita come sopra:
è ben definita
ha \(\hat x\) come unico punto fisso
è continua
La dimostrazione si fa per punti.
Supponiamo che \[f(f(z)) - f(z) - (f(z) - z) = 0\] Allora \[|f(f(z)) - f(z)| = |f(z) - z|\] e \[|f(f(z)) - f(z)| \leq L|f(z) - z|\] che implica che \(f(z) = z\) che è assurdo a meno che \(z= \hat x\).
Supponiamo che \(\hat x\) e \(z\) siano due punti fissi per \(\Phi_A\). Allora \[\begin{aligned} z &= \frac{zf(f(z))-f^2(z)}{f(f(z)) - 2f(z)+z} \\ z(f(f(z)) - 2f(z) +z) &= zf(f(z)) - f^2(z) \\ -2zf(z)+z^2 &=-f^2(z) \\ f(z) - z &= 0 \end{aligned}\] Ma se ciò accade, \(\Phi_A\) non è più ben definita.
Per il teorema di de l’Hopital,
\[\lim_{z\rightarrow \hat x} \frac{f(f(z))+zf'(f(z))\cdot f'(z)-2f(z)f'(z)}{1-2f'(z)+f'(f(z))f'(z)} = \lim\hat x\frac{1-2f'(z)+f'(f(z))f'(z)}{1-2f'(z)+f'(f(z))f'(z)} = \hat x\]
◻
Alcune proprietà fondamentali.
Se la successione generata dal metodo che si vuole ottimizzare ha convergenza lineare con costante asintotica \(c\in (0, 1)\), si ottiene \[|x_{k+2}-\hat x| \sim c|x_{k+1}-\hat x| \sim c^2|x_k -\hat x|\] dunque, la successione base degli indici pari ha ordine di convergenza comunque lineare, ma con costante asintotica più piccola. Si riesce a dimostrare, tuttavia, che l’ordine di convergenza della successione di Aitken a lei associata converge quadraticamente, seppur con il doppio dei flops per passo, dunque conviene.
Se la successione generata ha ordine di convergenza \(p>1\), non è più conveniente utilizzare Aitken, in quanto \[|x_{k+2}-\hat x| \sim c|x_{k+1}-\hat x|^p \sim c^{p+1}|x_k -\hat x|^{p^2}\] mentre si ha anche \[|z_{k+1} -\hat x| \sim c|z_k-\hat x|^{2p-1}\] con \(p^2>2p-1\).
Il metodo di bisezione è il primo metodo che vediamo di ricerca degli zeri di \(f\in C^1([a_0, b_0])\), e si basa fondamentalmente sul teorema di esistenza degli zeri per funzioni ad una variabile. La successione per ricorrenza è definita come
\[x_{k+1} = \frac{a_k +b_k}{2} \text{ dove } [a_{k+1}, b_{k+1}] = \begin{cases} [a_k, x_k] \text{ se } f(a_k)f(x_k)<0 \\ [x_k, b_k] \text{ altrimenti} \end{cases}\]
L’intervallo \([a_k, b_k]\) viene comunemente detto intervallo di confidenza, perché è quello nel quale ci aspettiamo che si trovi lo zero di funzione. L’errore che si commette è maggiorabile con
\[e_k \leq |x_k- \hat x| \leq \frac{b-a}{2^{k-1}}:=\eta_k\]
La convergenza è lineare, globale ma non monotona.
Per assicurare che il metodo si arresti, un primo tentativo è quello di imporre un bound assoluto (e indipendente dal caso particolare trattato) al termine \(\eta_k\). Questo controllo non tiene di conto della pendenza della funzione, e implica diversi problemi nel numero di iterazioni fatte. Se ad esempio la funzione risulta piatta in un intorno dello zero ricercato, l’errore commesso approssimando \(x_k\) con \(\hat x\) è molto alto.
Per ovviare a questo problema, si può controllare che \[|f(x_k)| \leq |f'(\hat x)| \cdot \text{tolx}\] dove si può approssimare \[f'(\hat x) = \frac{f(b_k)-f(a_k)}{b_k-a_k}\] ottenendo \[f(x) = f(\hat x) + f'(\hat x)(x-\hat x) \implies \kappa =\frac{1}{f'(\hat x)}\] dove \(\kappa\) è numero di condizionamento del problema.
E’ chiaro che il problema risulti malcondizionato se \(\hat x\) è uno zero multiplo. In tal caso possiamo scegliere di prendere il valore della prima derivata n-esima diversa da zero, oppure scegliere di calcolare tutte le iterazioni.
Il metodo di Newton si applica a funzioni \(f\in C^1(K)\) ed è un metodo a convergenza locale, definito dalla successione
\[x_{k+1} = x_k- \frac{f(x_k)}{f'(x_k)}\]
Questo metodo si generalizza facilmente al caso multidimensionale, e la formula si può interpretare sia in modo analitico, che in modo geometrico, come metodo delle tangenti. Abbiamo due risultati generali per zeri semplici o multipli che dimostriamo.
Teorema 3 (Convergenza locale di Newton per zeri semplici). Sia \(f\in C^2(U_{\hat x})\), intorno di \(\hat x\) zero semplice, e supponiamo che la successione di Newton associata sia ben definita e convergente. Allora se \(f''(\hat x) = 0\), \(p =2\), altrimenti \(p > 2\).
Si ha che \[\begin{aligned} 0 &= f(\hat x)=f(x_k)+f'(x_k)(\hat x-x_k) + \frac{f''(\xi_k)(\hat x-x_k)^2}{2}\\ &=f'(x_k)[\frac{f(x_k)}{f'(x_k)} + \hat x - x_k ] +\frac{f''(\xi_k)(\hat x-x_k)^2}{2}\\ &\implies f'(x_k)(x_{k+1}-\hat x) = \frac{f''(\xi_k)e_k^2}{2} \\ &\implies \frac{e_{k+1}}{e_k}=\frac{f''(\xi_k)}{f'(x_k)}\rightarrow\frac{f''(\hat x)}{f'(\hat x)} \end{aligned}\] ◻
Il prossimo teorema si applica anche con l’ipotesi più debole che \(f\in C^2(U_{\hat x})\), ma per semplicità vedremo il risultato quando \(f\in C^{M+1}\).
Teorema 4 (Convergenza locale di Newton per zeri multipli). Sia \(f\in C^{M+1}(U_{\hat x})\), intorno di \(\hat x\) zero multiplo con molteplicità \(M\), e supponiamo che la successione di Newton associata sia ben definita e convergente. Allora la convergenza è lineare con costante \(c = 1-\frac{1}{M}\).
Infatti, supponendno \(f\) sviluppabile in serie di Taylor, otteniamo
\[f(x) = (x-\hat x)^Mg(x)\]
Allora,
\[\begin{aligned} \frac{e_{k+1}}{e_k} &= \frac{x_{k+1} - \hat x}{e_k} = \frac{ e_k - \frac{f(x_k)}{f'(x_k)} }{e_k} \\ &= 1 - \frac{1}{e_k}\cdot \Bigg[ \frac{(x_k-\hat x)^Mg(x_k)}{(x_k -\hat x)^{M-1}[Mg(x_k) + (x-\hat xg'(x_k) ]} \Bigg]\\ &= 1 - \Bigg[ \frac{Mg(x_k)}{Mg(x_k) + (x-\hat xg'(x_k) } \Bigg] \rightarrow 1-\frac{1}{M} \end{aligned}\] ◻
I teoremi che abbiamo visto adesso assicurano, sotto opportune ipotesi, soltanto la convergenza locale. Esiste un risultato che dimostra la convergenza "globale" su intervalli abbastanza grandi.
Teorema 5. (Convergenza globale di Newton)
Supponiamo che \(f\in C^2(\mathcal U _{\hat x})\). Supponiamo che la funzione \(f\) sia positiva e convessa, oppure negativa e concava su \((a, \hat x)\), e sia \(f'(\hat x)\neq 0\) su \((a, x)\). Allora, la successione con punto d’innesco in \((a, \hat x)\) è crescente e limitata superiormente da \(\hat x\), dunque è convergente.
Supponiamo per semplicità che \(f'>0\), e che \(f, f'' <0\) sull’intervallo \((a, \hat x)\). Vogliamo mostrare che
\[x_k < x_{k+1} <\hat x \quad \forall k\in\mathbb{N}\]
Dimostriamolo per induzione, e supponiamo \(k=0\). \[\begin{aligned} x_1=x_0 - \frac{f(x_0)}{f'(x_0)} \implies x_1 > x_0 \end{aligned}\]
Supponiamo per semplicità che \(k = 1\). Allora \[\begin{aligned} 0 &= f(\hat x) \\ &= f(x_0)+f'(x_0) (x-x_0)+\frac{f''(\xi_0)(\hat x- x_0)^2}{2} \\ &= f'(x_0)\Bigg[ \hat x - x_0 + \frac{f(x_0)}{f'(x_0)} \Bigg]+\frac{f''(\xi_0)(\hat x- x_0)^2}{2} \\ &= f'(x_0) (\hat x-x_1)+\frac{f''(\xi_0)(\hat x- x_0)^2}{2}\\ \implies (x_1 - \hat x) &= \frac{f''(\xi_0)(\hat x- x_0)^2}{2} \\ (x_1 - \hat x) &= \frac{1}{2}\frac{f''(\xi_0)}{f'(x_0)} (\hat x- x_0)^2 < 0\\ \implies x_0 &< x_1 < \hat x \end{aligned}\] ◻
Il metodo di Newton può essere applicato al caso particolare in cui si sia interessati al calcolo della radice n-esima di un numero. La funzione di cui si vuole cercare lo zero è \(f(x) = x^n-\alpha\). Si ottiene:
\[\begin{aligned} x_{k+1} &= x_k-\frac{x_k^n - \alpha}{nx_k^{n-1}}=\frac{(n-1)x_k^n + \alpha}{nx_k^{n-1}} \\ &= \frac{1}{n}\bigg( (n-1)x_k + \frac{\alpha}{x_k^{n-1}} \bigg) \end{aligned}\]
I metodi che introdurremo adesso sono metodi di quasi-Newton, ottenuti apportando modifiche al metodo di Newton per renderlo più veloce o per poterlo applicare su funzioni non \(C^2\).
Poiché il calcolo della derivata può risultare particolamente impegnativo, cerchiamo di approssimarla in prima battuta con una costante \(\mu\). Il metodo che otteniamo è il metodo delle corde, un metodo a iterazione funzionale, dato dalla seguente espressione:
\[x_{k+1} = \Phi(x_k) = x_k - \frac{f(x_k)}{\mu}\]
e con
\[\Phi'(x_k) = 1-\frac{f'(x_k)}{\mu}\]
Il termine \(|\Phi'(x_k)|\) deve essere imposto \(<1\), perchè si vuole realizzare la condizione di convergenza. Infatti, se questo accade, otteniamo subito che \(\Phi\) è una contrazione, e quindi, dai risultati visti all’inizio, genera una successione di punti convergente. Così facendo, si ottiene:
\[\begin{aligned} -1 &< \frac{f'(\hat x)}{\mu} < 1 \\ -2 &< \frac{f'(\hat x)}{\mu} \\ 2 &> \frac{f'(\hat x)}{\mu} \end{aligned}\]
Il metodo delle corde è un metodo poco efficiente, a convergenza lineare, anche per zeri semplici, dunque non conviene.
Il metodo delle secanti è un metodo con memoria (a passo doppio) e approssima meglio la derivata prima, utilizzando informazioni già calcolate precedentemente. Il fattore \(\mu = \mu_k\) viene calcolato nuovamente ad ogni iterazione con la formula
\[\mu_k = \frac{f(x_k) - f(x_{k-1})}{x_k - x_{k-1}}\]
Il nome è giustificato, infatti \(\mu_k\) è proprio la pendenza della retta secante nei punti \((x_k, f(x_k))\) e \((x_{k+1}, f(x_{k+1}))\).
Se lo zero è semplice, l’ordine di convergenza è \[p = \frac{1+\sqrt{5}}{2}\]
Da un punto di vista computazionale, sono necessarie alcune osservazioni. Per arrestare il metodo, si utilizza un criterio di arresto misto, dove si controlla, per ogni iterazione, che
\[\frac{|x_{k+1} - x_k|}{\text{tol}_\text{abs} + \text{tol}_\text{rel}|x_{k+1}|} \leq 1\]
Inoltre, in aritmetica finita, è possibile riarrangiare i termini per rendere il calcolo più stabile, ottenendo
\[x_{k+1} = \frac{x_{k-1}f(x_k) - x_k f(x_{k+1})}{f(x_k) - f({x_{k+1}})}\]
Proponiamo infine un esempio per il calcolo della radice quadrata di un numero con il metodo delle secanti. La funzione di cui si vuole cercare lo zero è \(f(x) = x^2 - \alpha\). Si ha \[\begin{aligned} x_{k+1} &= x_k - \frac{x_k - x_{k-1}}{x_k^2-x_{k-1}^2}(x_k^2 - \alpha) \\ &= x_k - \frac{x_k^2 - \alpha}{x_k + x_{k+1}}\\ &= \frac{x_{k-1}x_k - \alpha}{x_{k-1}+x_k} \end{aligned}\]