Calcolo della radice quadrata in modo iterativo
Questo algoritmo calcola una stima approssimativa della radice quadrata di un numero N tramite una funzione di trasformazione g(x), una stima iniziale x0 e un ciclo di iterazioni. $$ x_{n+1} = \frac{1}{2} \cdot ( x_n + \frac{N}{x_n} ) $$
L'algoritmo converge alla radice quadrata del numero N in modo quadratico.
Il computo termina quando il valore assoluto della differenza tra due stime successive è inferiore a una soglia che considero accettabile.
Ad esempio, 0,1 o altro.
$$ | x_{n+1} - x_n | < 0,1 $$
La soglia di accettazione determina l'accuratezza della stima e il numero di iterazioni.
Quanto più è bassa la soglia di accettazione, tanto più è accurata la stima e tante più iterazioni sono necessarie per ottenerla.
Nota. L'algoritmo che ho descritto è anche noto come "metodo babilonese" o "metodo di Erone". Erone di Alessandria era un matematico greco antico che ha descritto questo metodo nei suoi scritti, ma è noto che era utilizzato anche dai babilonesi molto tempo prima di lui.
Un esempio pratico
Considero il numero N=81
$$ 81 $$
Calcolo la stima della radice quadrata di 81 in modo iterativo
$$ \sqrt{81} $$
Fisso come soglia accettabile 0,1
$$ | x_{n+1} - x_n | < 0,1 $$
Scelgo come prima stima iniziale x0=1
$$ x_0 = 1 $$
Il valore iniziale x0 è a libera scelta. Qualsiasi valore va bene.
Nota. Tuttavia, dovendo scegliere un numero a caso è utile utilizzare la metà del radicando ossia x0=N/2 se N è pari oppure x0=(N+1)/2 se N è dispari. Spesso, questo fa risparmiare almeno un passaggio rispetto a x0=1.
A questo punto comincia il calcolo iterativo.
Iterazione 1
Calcolo la stima successiva x1 utilizzando la funzione di trasformazione g(x)
$$ x_1 = \frac{1}{2} \cdot ( x_0 + \frac{N}{x_0} ) $$
Dove N=81 e x0=1
$$ x_1 = \frac{1}{2} \cdot ( 1 + \frac{81}{1} ) $$
$$ x_1 = \frac{1}{2} \cdot 82 $$
$$ x_1 = 42 $$
La prima iterazione termina con una stima x1=42
La differenza rispetto alla stima iniziale x0=1 è ancora troppo grande
$$ | 42 - 1 | = 41 > 0,1 $$
Quindi, procedo con la seconda iterazione.
Iterazione 2
Nella seconda iterazione calcolo x2 a partire da x1=41
$$ x_2 = \frac{1}{2} \cdot ( x_1 + \frac{N}{x_1} ) $$
$$ x_2 = \frac{1}{2} \cdot ( 41 + \frac{81}{41} ) $$
$$ x_2 = 21,4878 $$
La differenza rispetto alla stima precedente x1=41 è ancora troppo grande.
$$ | 21,4878 - 41 | = 19.5122 > 0,1 $$
Quindi, procedo con la terza iterazione.
Iterazione 3
Questa volta calcolo x3 a partire da x2=21,4878
$$ x_3 = \frac{1}{2} \cdot ( x_2 + \frac{N}{x_2} ) $$
$$ x_3 = \frac{1}{2} \cdot ( 21,4878 + \frac{81}{21,4878} ) $$
$$ x_3 = 12,6287 $$
La differenza rispetto alla stima precedente x2=21,4878 si riduce ma è ancora grandeo.
$$ | 12,6287 - 21,4878 | = 8,8591 > 0,1 $$
Quindi, procedo con la quarta iterazione.
Iterazione 4
Nella quarta iterazione calcolo x4 a partire da x3=12,6287
$$ x_4 = \frac{1}{2} \cdot ( x_3 + \frac{N}{x_3} ) $$
$$ x_4 = \frac{1}{2} \cdot ( 12,6287 + \frac{81}{12,6287} ) $$
$$ x_4 = 9,5213 $$
La differenza rispetto alla stima precedente x3=12,6287 si riduce ulteriormente ma è ancora grande.
$$ | 9,5213 - 12,6287 | = 3.1074 > 0,1 $$
Quindi, procedo con la quinta iterazione.
Iterazione 5
Calcolo la stima x5 a partire da x4=9,5213
$$ x_5 = \frac{1}{2} \cdot ( x_4 + \frac{N}{x_4} ) $$
$$ x_5 = \frac{1}{2} \cdot ( 9,5213 + \frac{81}{9,5213} ) $$
$$ x_5 = 9,0143 $$
La differenza assoluta rispetto alla stima precedente x4=9,5213 è ancora troppo grande
$$ | 9,0143-9,513| = 0,4987 > 0,1 $$
Quindi procedo con la sesta iterazione.
Iterazione 6
Calcolo x6 a partire da x5=9,0143
$$ x_6 = \frac{1}{2} \cdot ( x_5 + \frac{N}{x_5} ) $$
$$ x_6 = \frac{1}{2} \cdot ( 9,0143 + \frac{81}{9,0143} ) $$
$$ x_6 = 9,000011 $$
Ora la differenza tra le stime x6 e x5 è sufficientemente piccola.
E' inferiore alla soglia di accettazione che ho fissato a 0,1
$$ | 9,000011 - 9,0143 | = 0,014 < 0,1 $$
L'algoritmo termina qui. L'ultimo valore trovato x6=9,000011 è il valore approssimato della radice quadrata di 81.
Nota. Riducendo ulteriormente la soglia di accettazione posso trovare un valore approssimato più accurato, al costo di compiere ulteriori iterazioni.
Osservazioni
Alcune osservazioni sull'algoritmo
- La stima iniziale x0 è a libera scelta. Qualsiasi valore va bene. Ad esempio, avrei potuto scegliere x0=1000 o x0=9999. In ogni caso l'algoritmo converge verso la stima della radice quadrata di N=81. Ovviamente, cambia il numero delle iterazioni necessarie per raggiungerla.
- Dovendo scegliere un numero a caso è utile utilizzare la metà del radicando ossia x0=N/2 se N è pari oppure x0=(N+1)/2 se è dispari.
- L'algoritmo calcola una stima della radice quadrata di qualsiasi numero. Ad esempio, se scelgo N=15129, l'algoritmo restituisce la sua radice quadrata, ossia 123, in dieci iterazioni.
E così via.
