Algoritmi di ricerca delle radici - Root-finding algorithms

In matematica e informatica , un algoritmo di ricerca delle radici è un algoritmo per trovare gli zeri , detti anche "radici", di funzioni continue . Uno zero di una funzione f , dai numeri reali ai numeri reali o dai numeri complessi ai numeri complessi, è un numero x tale che f ( x ) = 0 . Poiché, in genere, gli zeri di una funzione non possono essere calcolati esattamente né espressi in forma chiusa , gli algoritmi di ricerca delle radici forniscono approssimazioni agli zeri, espressi o come numeri in virgola mobile o come piccoli intervalli isolanti , o dischi per radici complesse (un intervallo o disco l'output è equivalente a un output approssimativo insieme a un limite di errore).

Risolvere un'equazione f ( x ) = g ( x ) equivale a trovare le radici della funzione h ( x ) = f ( x ) – g ( x ) . Pertanto, gli algoritmi di ricerca delle radici consentono di risolvere qualsiasi equazione definita da funzioni continue. Tuttavia, la maggior parte degli algoritmi di ricerca delle radici non garantisce che troveranno tutte le radici; in particolare, se tale algoritmo non trova alcuna radice, ciò non significa che non esista alcuna radice.

La maggior parte dei metodi numerici di ricerca della radice utilizza l' iterazione , producendo una sequenza di numeri che si spera convergano verso la radice come limite . Richiedono una o più ipotesi iniziali della radice come valori di partenza, quindi ogni iterazione dell'algoritmo produce un'approssimazione della radice successivamente più accurata. Poiché l'iterazione deve essere interrotta a un certo punto, questi metodi producono un'approssimazione alla radice, non una soluzione esatta. Molti metodi calcolano i valori successivi valutando una funzione ausiliaria sui valori precedenti. Il limite è quindi un punto fisso della funzione ausiliaria, che viene scelta per avere come punti fissi le radici dell'equazione originaria, e per convergere rapidamente a questi punti fissi.

Il comportamento degli algoritmi generali di ricerca delle radici è studiato nell'analisi numerica . Tuttavia, per i polinomi, lo studio della ricerca delle radici appartiene generalmente all'algebra del computer , poiché le proprietà algebriche dei polinomi sono fondamentali per gli algoritmi più efficienti. L'efficienza di un algoritmo può dipendere notevolmente dalle caratteristiche delle funzioni date. Ad esempio, molti algoritmi utilizzano la derivata della funzione di input, mentre altri lavorano su ogni funzione continua . In generale, non è garantito che gli algoritmi numerici trovino tutte le radici di una funzione, quindi non riuscire a trovare una radice non dimostra che non esiste una radice. Tuttavia, per i polinomi , esistono algoritmi specifici che utilizzano proprietà algebriche per certificare che non viene persa alcuna radice e per localizzare le radici in intervalli separati (o dischi per radici complesse) sufficientemente piccoli da garantire la convergenza dei metodi numerici (tipicamente il metodo di Newton ) alla radice univoca così situata.

Metodi di bracketing

I metodi di parentesi quadre determinano intervalli successivi più piccoli (parentesi) che contengono una radice. Quando l'intervallo è sufficientemente piccolo, è stata trovata una radice. In genere usano il teorema del valore intermedio , che afferma che se una funzione continua ha valori di segno opposto ai punti finali di un intervallo, allora la funzione ha almeno una radice nell'intervallo. Pertanto, richiedono di iniziare con un intervallo tale che la funzione assuma segni opposti ai punti finali dell'intervallo. Tuttavia, nel caso dei polinomi esistono altri metodi ( regola dei segni di Cartesio , teorema di Budan e teorema di Sturm ) per ottenere informazioni sul numero di radici in un intervallo. Conducono a algoritmi efficienti per l' isolamento della radice reale dei polinomi, che assicurano la ricerca di tutte le radici reali con una precisione garantita.

Metodo di bisezione

L'algoritmo di ricerca della radice più semplice è il metodo di bisezione . Lasciate f sia una funzione continua , per cui si sa un intervallo [ a , b ] tale che f ( a ) e f ( b ) hanno segni opposti (una staffa). Sia c = ( a + b )/2 il centro dell'intervallo (il punto medio o il punto che biseca l'intervallo). Allora o f ( a ) e f ( c ) , oppure f ( c ) e f ( b ) hanno segni opposti, e uno ha diviso per due la dimensione dell'intervallo. Sebbene il metodo di bisezione sia robusto, guadagna uno e solo un bit di precisione con ogni iterazione. Altri metodi, in condizioni appropriate, possono acquisire precisione più rapidamente.

Falsa posizione ( regula falsi )

Il metodo della falsa posizione , chiamato anche metodo regula falsi , è simile al metodo della bisezione, ma invece di utilizzare il centro dell'intervallo della ricerca di bisezione utilizza l' intercetta x della linea che collega i valori della funzione tracciata agli estremi dell'intervallo , questo è

La posizione falsa è simile al metodo secante , tranne per il fatto che, invece di mantenere gli ultimi due punti, si assicura di mantenere un punto su entrambi i lati della radice. Il metodo della falsa posizione può essere più veloce del metodo della bisezione e non diverge mai come il metodo secante; tuttavia, potrebbe non riuscire a convergere in alcune implementazioni ingenue a causa di errori di arrotondamento che possono portare a un segno errato per f ( c ) ; tipicamente, ciò può verificarsi se la velocità di variazione di f è grande nelle vicinanze della radice.

Metodo ITP

Il metodo ITP è l'unico metodo noto per racchiudere la radice con le stesse garanzie del caso peggiore del metodo di bisezione, garantendo al contempo una convergenza superlineare alla radice delle funzioni lisce come il metodo secante. È anche l'unico metodo noto garantito per superare il metodo di bisezione in media per qualsiasi distribuzione continua sulla posizione della radice (vedi ITP Method#Analysis ). Lo fa tenendo traccia sia dell'intervallo di bracketing che dell'intervallo minmax in cui ogni punto al suo interno converge alla stessa velocità del metodo di bisezione. La costruzione del punto c interrogato segue tre fasi: interpolazione (simile alla regula falsi), troncamento (aggiustando la regula falsi simile a Regula falsi § Miglioramenti in regula falsi ) e quindi proiezione sull'intervallo minmax. La combinazione di questi passaggi produce simultaneamente un metodo ottimale minmax con garanzie simili ai metodi basati sull'interpolazione per funzioni regolari e, in pratica, supererà sia il metodo di bisezione che i metodi basati sull'interpolazione sia per funzioni regolari che non regolari.

interpolazione

Molti processi di ricerca delle radici funzionano per interpolazione . Consiste nell'utilizzare gli ultimi valori approssimati calcolati della radice per approssimare la funzione con un polinomio di basso grado, che assume gli stessi valori in corrispondenza di queste radici approssimate. Quindi la radice del polinomio viene calcolata e utilizzata come nuovo valore approssimato della radice della funzione e il processo viene iterato.

Due valori consentono di interpolare una funzione con un polinomio di grado uno (ovvero approssimare il grafico della funzione con una retta). Questa è la base del metodo secante . Tre valori definiscono una funzione quadratica , che approssima il grafico della funzione con una parabola . Questo è il metodo di Muller .

Regula falsi è anche un metodo di interpolazione, che differisce dal metodo secante utilizzando, per l'interpolazione con una linea, due punti che non sono necessariamente gli ultimi due punti calcolati.

Metodi iterativi

Sebbene tutti gli algoritmi di ricerca della radice procedano per iterazione , un metodo di ricerca della radice iterativo generalmente utilizza un tipo specifico di iterazione, consistente nel definire una funzione ausiliaria, che viene applicata alle ultime approssimazioni calcolate di una radice per ottenere una nuova approssimazione. L'iterazione si ferma quando viene raggiunto un punto fisso ( fino alla precisione desiderata) della funzione ausiliaria, cioè quando il nuovo valore calcolato è sufficientemente vicino ai precedenti.

Metodo di Newton (e metodi simili basati su derivati)

Il metodo di Newton assume che la funzione f abbia una derivata continua . Il metodo di Newton potrebbe non convergere se avviato troppo lontano da una radice. Tuttavia, quando converge, è più veloce del metodo di bisezione ed è solitamente quadratico. Il metodo di Newton è importante anche perché si generalizza facilmente a problemi di dimensioni superiori. I metodi di tipo Newton con ordini di convergenza superiori sono i metodi del capofamiglia . Il primo dopo il metodo di Newton è il metodo di Halley con ordine di convergenza cubico.

metodo secante

Sostituendo la derivata nel metodo di Newton con una differenza finita , otteniamo il metodo secante . Questo metodo non richiede il calcolo (né l'esistenza) di un derivato, ma il prezzo è di convergenza più lenta (l'ordine è di circa 1,6 ( rapporto aureo )). Una generalizzazione del metodo secante in dimensioni superiori è il metodo di Broyden .

Il metodo di Steffensen

Se usiamo un fit polinomiale per rimuovere la parte quadratica della differenza finita usata nel metodo Secante, in modo che approssimi meglio la derivata, otteniamo il metodo di Steffensen , che ha convergenza quadratica, e il cui comportamento (sia buono che cattivo) è essenzialmente lo stesso del metodo di Newton ma non richiede una derivata.

Interpolazione inversa

La comparsa di valori complessi nei metodi di interpolazione può essere evitata interpolando l' inverso di f , ottenendo il metodo di interpolazione quadratica inversa . Di nuovo, la convergenza è asintoticamente più veloce del metodo secante, ma l'interpolazione quadratica inversa spesso si comporta male quando le iterazioni non sono vicine alla radice.

Combinazioni di metodi

Il metodo di Brent

Il metodo di Brent è una combinazione del metodo di bisezione, del metodo secante e dell'interpolazione quadratica inversa . Ad ogni iterazione, il metodo di Brent decide quale metodo tra questi tre è probabilmente quello migliore e procede facendo un passo secondo quel metodo. Questo dà un metodo robusto e veloce, che gode quindi di una notevole popolarità.

Metodo dei cavalieri

Il metodo di Ridders è un metodo ibrido che utilizza il valore della funzione nel punto medio dell'intervallo per eseguire un'interpolazione esponenziale alla radice. Ciò fornisce una rapida convergenza con una convergenza garantita di al massimo il doppio del numero di iterazioni del metodo di bisezione.

Radici di polinomi

Trovare le radici del polinomio è un problema di vecchia data che è stato oggetto di molte ricerche nel corso della storia. Una testimonianza di ciò è che fino al XIX secolo l' algebra significava essenzialmente teoria delle equazioni polinomiali .

Trovare la radice di un polinomio lineare (grado uno) è facile e richiede solo una divisione. Per i polinomi quadratici (grado due), la formula quadratica produce una soluzione, ma la sua valutazione numerica potrebbe richiedere una certa attenzione per garantire la stabilità numerica . Per i gradi tre e quattro esistono soluzioni in forma chiusa in termini di radicali , che generalmente non sono convenienti per la valutazione numerica, in quanto troppo complicate e che comportano il calcolo di più radici n- esime il cui calcolo non è più facile del calcolo diretto della radici del polinomio (ad esempio l'espressione delle radici reali di un polinomio cubico può coinvolgere radici cubiche non reali ). Per i polinomi di grado cinque o superiore il teorema di Abel-Ruffini afferma che non esiste, in generale, alcuna espressione radicale delle radici.

Quindi, ad eccezione di gradi molto bassi, la ricerca delle radici dei polinomi consiste nel trovare approssimazioni delle radici. Per il teorema fondamentale dell'algebra , si sa che un polinomio di grado n ha al massimo n radici reali o complesse, e questo numero si raggiunge per quasi tutti i polinomi.

Ne consegue che il problema della ricerca delle radici per i polinomi può essere suddiviso in tre diversi sottoproblemi;

  • Trovare una radice
  • Trovare tutte le radici
  • Trovare le radici in una specifica regione del piano complesso , tipicamente le radici reali o le radici reali in un dato intervallo (ad esempio, quando le radici rappresentano una grandezza fisica, interessano solo quelle reali positive).

Per trovare una radice, il metodo di Newton e altri metodi iterativi generali funzionano generalmente bene.

Per trovare tutte le radici, il metodo più antico è, una volta trovata una radice r , dividere il polinomio per xr , e riavviare iterativamente la ricerca di una radice del quoziente polinomio. Tuttavia, ad eccezione dei gradi bassi, questo non funziona bene a causa dell'instabilità numerica : il polinomio di Wilkinson mostra che una modifica molto piccola di un coefficiente può cambiare drasticamente non solo il valore delle radici, ma anche la loro natura (reale o complessa). Inoltre, anche con una buona approssimazione, quando si valuta un polinomio con una radice approssimata, si può ottenere un risultato che è lontano dallo zero. Ad esempio, se un polinomio di grado 20 (il grado del polinomio di Wilkinson) ha una radice prossima a 10, la derivata del polinomio alla radice può essere dell'ordine di questo implica che un errore di sul valore della radice può produrre un valore del polinomio alla radice approssimata che è dell'ordine di

Per evitare questi problemi, sono stati elaborati metodi che calcolano tutte le radici contemporaneamente, con la precisione desiderata. Attualmente il metodo più efficiente è il metodo Aberth . Una libera implementazione è disponibile sotto il nome di MPSolve . Questa è un'implementazione di riferimento, che può trovare abitualmente le radici di polinomi di grado maggiore di 1.000, con più di 1.000 cifre decimali significative.

I metodi per calcolare tutte le radici possono essere usati per calcolare le radici reali. Tuttavia, può essere difficile decidere se una radice con una piccola parte immaginaria è reale o meno. Inoltre, poiché il numero delle radici reali è, in media, il logaritmo del grado, è uno spreco di risorse del computer calcolare le radici non reali quando si è interessati alle radici reali.

Il metodo più antico per calcolare il numero di radici reali e il numero di radici in un intervallo risulta dal teorema di Sturm , ma i metodi basati sulla regola dei segni di Cartesio e le sue estensioni - i teoremi di Budan e Vincent - sono generalmente più efficienti. Per la ricerca delle radici, tutto procede riducendo la dimensione degli intervalli in cui vengono cercate le radici fino a ottenere intervalli contenenti zero o una radice. Quindi gli intervalli contenenti una radice possono essere ulteriormente ridotti per ottenere una convergenza quadratica del metodo di Newton alle radici isolate. I principali sistemi di computer algebra ( Maple , Mathematica , SageMath , PARI/GP ) hanno ciascuno una variante di questo metodo come algoritmo predefinito per le radici reali di un polinomio.

Un'altra classe di metodi si basa sulla conversione del problema della ricerca delle radici polinomiali nel problema della ricerca degli autovalori della matrice compagna del polinomio. In linea di principio, si può usare qualsiasi algoritmo di autovalori per trovare le radici del polinomio. Tuttavia, per ragioni di efficienza si preferiscono metodi che impieghino la struttura della matrice, cioè che possano essere implementati in forma matrix-free. Tra questi metodi vi è il metodo delle potenze , la cui applicazione alla trasposizione della matrice compagna è il classico metodo di Bernoulli per trovare la radice del modulo massimo. Il metodo della potenza inversa con spostamenti, che trova prima una radice più piccola, è ciò che guida la variante complessa ( cpoly ) dell'algoritmo di Jenkins-Traub e gli conferisce la sua stabilità numerica. Inoltre, è insensibile a più radici e ha una rapida convergenza con l'ordine (dove è il rapporto aureo ) anche in presenza di radici cluster. Questa rapida convergenza comporta un costo di tre valutazioni polinomiali per passo, risultando in un residuo di O (| f ( x )| 2+3 φ ) , ovvero una convergenza più lenta rispetto a tre passi del metodo di Newton.

Trovare una radice

Il metodo più utilizzato per il calcolo di una radice è il metodo di Newton , che consiste nelle iterazioni del calcolo di

partendo da un valore ben scelto

Se f è un polinomio, il calcolo è più veloce quando si utilizza il metodo di Horner o la valutazione con pre-elaborazione per calcolare il polinomio e la sua derivata in ogni iterazione.

La convergenza è generalmente quadratica , può convergere molto lentamente o addirittura non convergere affatto. In particolare, se il polinomio non ha radice reale, ed è reale, allora il metodo di Newton non può convergere. Tuttavia, se il polinomio ha una radice reale, che è più grande della radice reale più grande della sua derivata, allora il metodo di Newton converge quadraticamente a questa radice più grande se è più grande di questa radice più grande (ci sono modi semplici per calcolare un limite superiore della radici, vedere Proprietà delle radici polinomiali ). Questo è il punto di partenza del metodo Horner per il calcolo delle radici.

Quando è stata trovata una radice r , si può usare la divisione euclidea per rimuovere il fattore xr dal polinomio. Calcolare una radice del quoziente risultante e ripetere il processo fornisce, in linea di principio, un modo per calcolare tutte le radici. Tuttavia, questo schema iterativo è numericamente instabile; gli errori di approssimazione si accumulano durante le successive fattorizzazioni, cosicché le ultime radici sono determinate con un polinomio che si discosta ampiamente da un fattore del polinomio originario. Per ridurre questo errore, si può, per ogni radice trovata, riavviare il metodo di Newton con il polinomio originale, e questa radice approssimata come valore iniziale.

Tuttavia, non vi è alcuna garanzia che ciò consentirà di trovare tutte le radici. Infatti, il problema di trovare le radici di un polinomio dai suoi coefficienti è in generale altamente mal condizionato . Ciò è illustrato dal polinomio di Wilkinson : le radici di questo polinomio di grado 20 sono i primi 20 interi positivi; cambiando l'ultimo bit della rappresentazione a 32 bit di uno dei suoi coefficienti (uguale a –210) si produce un polinomio con solo 10 radici reali e 10 radici complesse con parti immaginarie maggiori di 0,6.

Strettamente correlati al metodo di Newton sono il metodo di Halley e il metodo di Laguerre . Entrambi usano il polinomio e le sue due prime derivazioni per un processo iterativo che ha una convergenza cubica . Combinando due passaggi consecutivi di questi metodi in un unico test, si ottiene un tasso di convergenza di 9, al costo di 6 valutazioni polinomiali (con regola di Horner). D'altra parte, la combinazione di tre fasi del metodo di Newton fornisce un tasso di convergenza di 8 al costo dello stesso numero di valutazioni polinomiali. Questo dà un leggero vantaggio a questi metodi (meno chiaro per il metodo di Laguerre, poiché ad ogni passo deve essere calcolata una radice quadrata).

Quando si applicano questi metodi a polinomi con coefficienti reali e punti di partenza reali, i metodi di Newton e Halley rimangono all'interno della linea dei numeri reali. Bisogna scegliere punti di partenza complessi per trovare radici complesse. Al contrario, il metodo Laguerre con una radice quadrata nella sua valutazione lascerà spontaneamente l'asse reale.

Trovare le radici in coppia

Se il dato polinomio ha solo coefficienti reali, si potrebbe voler evitare calcoli con numeri complessi. A tal fine, si devono trovare fattori quadratici per coppie di radici complesse coniugate. L'applicazione del metodo di Newton multidimensionale a questo compito risulta nel metodo di Bairstow .

La vera variante dell'algoritmo di Jenkins-Traub è un miglioramento di questo metodo.

Trovare tutte le radici in una volta

Il semplice Durand-Kerner e il metodo Aberth leggermente più complicato trovano simultaneamente tutte le radici usando solo l' aritmetica dei numeri complessi semplici . Algoritmi accelerati per la valutazione e l'interpolazione multipunto simili alla trasformata di Fourier veloce possono aiutare ad accelerarli per grandi gradi del polinomio. Si consiglia di scegliere un insieme di punti iniziali asimmetrici, ma uniformemente distribuiti. L'implementazione di questo metodo nel software gratuito MPSolve è un riferimento per la sua efficienza e la sua precisione.

Un altro metodo con questo stile è il metodo Dandelin-Gräffe (a volte anche attribuito a Lobachevsky ), che utilizza trasformazioni polinomiali per quadrare ripetutamente e implicitamente le radici. Ciò amplifica notevolmente le variazioni nelle radici. Applicando le formule di Viète si ottengono facili approssimazioni per il modulo delle radici e, con qualche sforzo in più, per le radici stesse.

Metodi di esclusione e inclusione

Esistono diversi test veloci che dicono se un segmento della linea reale o una regione del piano complesso non contiene radici. Limitando il modulo delle radici e suddividendo ricorsivamente la regione iniziale indicata da questi limiti, è possibile isolare piccole regioni che possono contenere radici e quindi applicare altri metodi per localizzarle esattamente.

Tutti questi metodi implicano la ricerca dei coefficienti delle versioni spostate e scalate del polinomio. Per grandi gradi, i metodi accelerati basati su FFT diventano praticabili.

Per le radici reali, vedere le sezioni successive.

L' algoritmo di Lehmer-Schur utilizza il test di Schur-Cohn per i cerchi; una variante, l'algoritmo di bisezione globale di Wilf utilizza un calcolo del numero di avvolgimento per regioni rettangolari nel piano complesso.

Il metodo del cerchio di divisione utilizza trasformazioni polinomiali basate su FFT per trovare fattori di grande grado corrispondenti a gruppi di radici. La precisione della fattorizzazione viene massimizzata utilizzando un'iterazione di tipo Newton. Questo metodo è utile per trovare le radici di polinomi di grado elevato a precisione arbitraria; ha una complessità quasi ottimale in questa impostazione.

Isolamento della radice reale

Trovare le vere radici di un polinomio con coefficienti reali è un problema che ha ricevuto molta attenzione dall'inizio del 19° secolo, ed è ancora un campo di ricerca attivo. La maggior parte degli algoritmi di ricerca delle radici può trovare alcune radici reali, ma non può certificare di aver trovato tutte le radici. I metodi per trovare tutte le radici complesse, come il metodo Aberth, possono fornire le radici reali. Tuttavia, a causa dell'instabilità numerica dei polinomi (vedi il polinomio di Wilkinson ), potrebbero aver bisogno di un'aritmetica di precisione arbitraria per decidere quali radici sono reali. Inoltre, calcolano tutte le radici complesse quando solo poche sono reali.

Ne segue che il modo standard di calcolare le radici reali è calcolare i primi intervalli disgiunti, detti intervalli isolanti , in modo che ciascuno contenga esattamente una radice reale, e insieme contengano tutte le radici. Questo calcolo è chiamato isolamento della radice reale . Avendo l'intervallo di isolamento, si possono usare metodi numerici veloci, come il metodo di Newton per migliorare la precisione del risultato.

Il più antico algoritmo completo per l'isolamento della radice reale risulta dal teorema di Sturm . Tuttavia, sembra essere molto meno efficiente dei metodi basati sulla regola dei segni di Cartesio e sul teorema di Vincent . Questi metodi si dividono in due classi principali, una che utilizza le frazioni continue e l'altra che utilizza la bisezione. Entrambi i metodi sono stati notevolmente migliorati dall'inizio del 21° secolo. Con questi miglioramenti raggiungono una complessità computazionale simile a quella dei migliori algoritmi per calcolare tutte le radici (anche quando tutte le radici sono reali).

Questi algoritmi sono stati implementati e sono disponibili in Mathematica (metodo della frazione continua) e Maple (metodo della bisezione). Entrambe le implementazioni possono trovare abitualmente le radici reali di polinomi di grado maggiore di 1.000.

Trovare più radici di polinomi

La maggior parte degli algoritmi di ricerca delle radici si comporta male quando ci sono più radici o radici molto vicine. Tuttavia, per i polinomi i cui coefficienti sono dati esattamente come interi o numeri razionali , esiste un metodo efficiente per fattorizzarli in fattori che hanno solo radici semplici e i cui coefficienti sono anch'essi dati esattamente. Questo metodo, chiamato fattorizzazione senza quadrati , si basa sul fatto che le radici multiple di un polinomio sono le radici del massimo comun divisore del polinomio e della sua derivata.

La fattorizzazione senza quadrati di un polinomio p è una fattorizzazione in cui ciascuno è 1 o un polinomio senza radici multiple e due diversi non hanno radice comune.

Un metodo efficiente per calcolare questa fattorizzazione è l'algoritmo di Yun .

Guarda anche

Riferimenti

  • premere, WH; Teukolsky, SA; Vetterling, WT; Flannery, BP (2007). "Capitolo 9. Determinazione della radice e insiemi non lineari di equazioni" . Ricette numeriche: L'arte del calcolo scientifico (3a ed.). New York: Cambridge University Press. ISBN 978-0-521-88068-8.
  1. ^ "Radici polinomiali - Radici MATLAB" . MathWorks . 2021-03-01 . Estratto 2021-09-20 .