15. Aritmetica con i Numeri in Virgola Mobile: Problemi e Limitazioni

I numeri in virgola mobile sono rappresentati nell’hardware del computer come frazioni in base 2 (binari). Ad esempio, la frazione decimale 0.625 ha valore 6/10 + 2/100 + 5/1000, e allo stesso modo la frazione binaria 0.101 ha valore 1/2 + 0/4 + 1/8. Queste due frazioni hanno valori identici, l’unica vera differenza è che la prima è scritta in notazione frazionaria in base 10, e la seconda in base 2.

Sfortunatamente, la maggior parte delle frazioni decimali non può essere rappresentata esattamente come frazioni binarie. Una conseguenza è che, in generale, i numeri in virgola mobile decimali che inserisci sono solo approssimati dai numeri in virgola mobile binari effettivamente memorizzati nella macchina.

Il problema è più facile da comprendere inizialmente in base 10. Considera la frazione 1/3. Puoi approssimarla come una frazione in base 10:

0.3

o, meglio,

0.33

o, meglio,

0.333

e così via. Non importa quanti cifre sei disposto a scrivere, il risultato non sarà mai esattamente 1/3, ma sarà una approssimazione sempre migliore di 1/3.

Allo stesso modo, non importa quante cifre in base 2 sei disposto a usare, il valore decimale 0.1 non può essere rappresentato esattamente come una frazione in base 2. In base 2, 1/10 è la frazione che si ripete all’infinito:

0.0001100110011001100110011001100110011001100110011...

Fermati a un qualsiasi numero finito di bit, e otterrai un’approssimazione. Nella maggior parte delle macchine odierne, i float sono approssimati usando una frazione binaria con il numeratore che utilizza i primi 53 bit a partire dal bit più significativo e con il denominatore come potenza di due. Nel caso di 1/10, la frazione binaria è 3602879701896397 / 2 ** 55 che è vicino ma non esattamente uguale al valore vero di 1/10.

Molti utenti non sono a conoscenza dell’approssimazione a causa del modo in cui i valori sono visualizzati. Python stampa solo un’approssimazione decimale al valore decimale vero dell’approssimazione binaria memorizzata nella macchina. Nella maggior parte delle macchine, se Python dovesse stampare il valore decimale vero dell’approssimazione binaria memorizzata per 0.1, dovrebbe visualizzare:

>>> 0.1
0.1000000000000000055511151231257827021181583404541015625

Questo è un numero di cifre maggiore di quelle che la maggior parte delle persone trova utile, quindi Python mantiene il numero di cifre gestibile visualizzando invece un valore arrotondato:

>>> 1 / 10
0.1

Ricorda solo, anche se il risultato stampato sembra il valore esatto di 1/10, il valore effettivamente memorizzato è la frazione binaria rappresentabile più vicina.

Interessante è che ci sono molti numeri decimali diversi che condividono la stessa frazione binaria approssimata più vicina. Ad esempio, i numeri 0.1 e 0.10000000000000001 e 0.1000000000000000055511151231257827021181583404541015625 sono tutti approssimati da 3602879701896397 / 2 ** 55. Poiché tutti questi valori decimali condividono la stessa approssimazione, qualsiasi di essi potrebbe essere visualizzato pur mantenendo l’invariante eval(repr(x)) == x.

Storicamente, il prompt di Python e la funzione integrata repr() avrebbero scelto quella con 17 cifre significative, 0.10000000000000001. A partire da Python 3.1, Python (sulla maggior parte dei sistemi) è ora in grado di scegliere il più corto di questi e visualizzare semplicemente 0.1.

Nota che questa è la natura stessa della virgola mobile binaria: non è un bug di Python e non è nemmeno un bug nel tuo codice. Vedrai lo stesso tipo di cosa in tutti i linguaggi che supportano l’aritmetica in virgola mobile del tuo hardware (anche se alcuni linguaggi potrebbero non visualizzare la differenza per impostazione predefinita o in tutte le modalità di output).

Per un output più gradevole, potresti voler usare la formattazione delle stringhe per produrre un numero limitato di cifre significative:

>>> format(math.pi, '.12g')  # give 12 significant digits
'3.14159265359'

>>> format(math.pi, '.2f')   # give 2 digits after the point
'3.14'

>>> repr(math.pi)
'3.141592653589793'

È importante rendersi conto che questo è, in un certo senso, un’illusione: stai semplicemente arrotondando la visualizzazione del vero valore della macchina.

Un’illusione può generare un’altra. Ad esempio, poiché 0.1 non è esattamente 1/10, sommare tre valori di 0.1 potrebbe non dare esattamente 0.3, nemmeno:

>>> 0.1 + 0.1 + 0.1 == 0.3
False

Inoltre, poiché lo 0.1 non può avvicinarsi al valore esatto di 1/10 e lo 0.3 non può avvicinarsi al valore esatto di 3/10, il pre-arrotondamento con la funzione round() non può aiutare:

>>> round(0.1, 1) + round(0.1, 1) + round(0.1, 1) == round(0.3, 1)
False

Anche se i numeri non possono essere resi più vicini ai loro valori esatti previsti, la funzione math.isclose() può essere utile per confrontare valori inesatti:

>>> math.isclose(0.1 + 0.1 + 0.1, 0.3)
True

In alternativa, la funzione round() può essere utilizzata per confrontare approssimazioni grossolane:

>>> round(math.pi, ndigits=2) == round(22 / 7, ndigits=2)
True

L’aritmetica con numeri in virgola mobile binari riserva molte sorprese di questo tipo. Il problema con «0.1» è spiegato in dettaglio preciso di seguito, nella sezione «Errore di Rappresentazione». Vedi Examples of Floating Point Problems per un riepilogo piacevole di come funziona la virgola mobile binaria e dei tipi di problemi comunemente riscontrati nella pratica. Vedi anche The Perils of Floating Point per un resoconto più completo di altre sorprese comuni.

Come viene detto vicino alla fine, «non ci sono risposte facili». Tuttavia, non essere troppo diffidente nei confronti della virgola mobile! Gli errori nelle operazioni float di Python sono ereditati dall’hardware in virgola mobile e nella maggior parte delle macchine sono dell’ordine di non più di 1 parte su 2**53 per operazione. Questo è più che adeguato per la maggior parte dei compiti, ma è necessario tenere a mente che non è aritmetica decimale e che ogni operazione float può soffrire di un nuovo errore di arrotondamento.

Sebbene esistano casi patologici, per la maggior parte dell’uso occasionale dell’aritmetica in virgola mobile otterrai il risultato che ti aspetti alla fine se semplicemente arrotondi la visualizzazione dei tuoi risultati finali al numero di cifre decimali che ti aspetti. str() di solito è sufficiente, e per un controllo più fine vedi gli specificatori di formato del metodo str.format() in Format String Syntax.

Per utilizzi che richiedono una rappresentazione decimale esatta, prova ad usare il modulo decimal che implementa l’aritmetica decimale adatta per applicazioni contabili e applicazioni ad alta precisione.

Un’altra forma di aritmetica esatta è supportata dal modulo fractions che implementa l’aritmetica basata su numeri razionali (così i numeri come 1/3 possono essere rappresentati esattamente).

Se sei un utente assiduo delle operazioni in virgola mobile, dovresti dare un’occhiata al pacchetto NumPy e a molti altri pacchetti per operazioni matematiche e statistiche forniti dal progetto SciPy. Vedi <https://scipy.org>.

Python fornisce strumenti che possono aiutare in quelle rare occasioni in cui davvero vuoi conoscere il valore esatto di un float. Il metodo float.as_integer_ratio() esprime il valore di un float come una frazione:

>>> x = 3.14159
>>> x.as_integer_ratio()
(3537115888337719, 1125899906842624)

Poiché il rapporto è esatto, può essere utilizzato per ricreare senza perdita il valore originale:

>>> x == 3537115888337719 / 1125899906842624
True

Il metodo float.hex() esprime un float in esadecimale (base 16), dando ancora una volta il valore esatto memorizzato dal tuo computer:

>>> x.hex()
'0x1.921f9f01b866ep+1'

Questa rappresentazione esadecimale precisa può essere utilizzata per ricostruire il valore del float esattamente:

>>> x == float.fromhex('0x1.921f9f01b866ep+1')
True

Poiché la rappresentazione è esatta, è utile per trasferire valori in modo affidabile tra diverse versioni di Python (indipendenza dalla piattaforma) e scambiare dati con altri linguaggi che supportano lo stesso formato (come Java e C99).

Un altro strumento utile è la funzione sum() che aiuta a mitigare la perdita di precisione durante la somma. Usa la precisione estesa per i passaggi di arrotondamento intermedi man mano che i valori vengono aggiunti a un totale incrementale. Questo può fare la differenza in termini di precisione complessiva in modo che gli errori non si accumulino al punto da influenzare il totale finale:

>>> 0.1 + 0.1 + 0.1 + 0.1 + 0.1 + 0.1 + 0.1 + 0.1 + 0.1 + 0.1 == 1.0
False
>>> sum([0.1] * 10) == 1.0
True

La funzione math.fsum() va oltre e traccia tutti i «bit perduti» man mano che i valori vengono aggiunti a un totale incrementale in modo che il risultato abbia solo un singolo arrotondamento. Questo è più lento di sum() ma sarà più accurato nei casi non comuni in cui input di grande magnitudine si annullano a vicenda lasciando una somma finale vicina allo zero:

>>> arr = [-0.10430216751806065, -266310978.67179024, 143401161448607.16,
...        -143401161400469.7, 266262841.31058735, -0.003244936839808227]
>>> float(sum(map(Fraction, arr)))   # Exact summation with single rounding
8.042173697819788e-13
>>> math.fsum(arr)                   # Single rounding
8.042173697819788e-13
>>> sum(arr)                         # Multiple roundings in extended precision
8.042178034628478e-13
>>> total = 0.0
>>> for x in arr:
...     total += x                   # Multiple roundings in standard precision
...
>>> total                            # Straight addition has no correct digits!
-0.0051575902860057365

15.1. Errore di Rappresentazione

Questa sezione spiega l’esempio «0.1» in dettaglio, e mostra come puoi eseguire un’analisi esatta di casi come questo da solo. Si assume una familiarità di base con la rappresentazione in virgola mobile binaria.

Errore di rappresentazione si riferisce al fatto che alcune (anzi, la maggior parte) delle frazioni decimali non possono essere rappresentate esattamente come frazioni binarie (base 2). Questa è la principale ragione per cui Python (o Perl, C, C++, Java, Fortran e molti altri) spesso non visualizzeranno il numero decimale esatto che ti aspetti.

Perché questo? 1/10 non è esattamente rappresentabile come frazione binaria. Dal 2000 in poi, quasi tutte le macchine utilizzano l’aritmetica in virgola mobile binaria IEEE 754, e quasi tutte le piattaforme mappano i float di Python ai valori «doppia precisione» IEEE 754 binary64. I valori IEEE 754 binary64 contengono 53 bit di precisione, quindi in input il computer cerca di convertire 0.1 nella frazione più vicina della forma J/2**N dove J è un intero contenente esattamente 53 bit. Riscrivendo:

1 / 10 ~= J / (2**N)

come:

J ~= 2**N / 10

e ricordando che J ha esattamente 53 bit (è >= 2**52 ma < 2**53), il miglior valore per N è 56:

>>> 2**52 <=  2**56 // 10  < 2**53
True

Cioè, 56 è l’unico valore per N che lascia J con esattamente 53 bit. Il miglior valore possibile per J è quindi quel quoziente arrotondato:

>>> q, r = divmod(2**56, 10)
>>> r
6

Poiché il resto è più della metà di 10, la migliore approssimazione si ottiene arrotondando verso l’alto:

>>> q+1
7205759403792794

Pertanto la migliore approssimazione possibile di 1/10 in doppia precisione IEEE 754 è:

7205759403792794 / 2 ** 56

Dividendo sia il numeratore che il denominatore per due si riduce la frazione a:

3602879701896397 / 2 ** 55

Nota che poiché abbiamo arrotondato verso l’alto, questo è in realtà un po” più grande di 1/10; se non avessimo arrotondato verso l’alto, il quoziente sarebbe stato un po” più piccolo di 1/10. Ma in nessun caso può essere esattamente 1/10!

Quindi il computer non «vede» mai 1/10: quello che vede è la frazione esatta data sopra, la migliore approssimazione doppia precisione IEEE 754 che può ottenere:

>>> 0.1 * 2 ** 55
3602879701896397.0

Se moltiplichiamo quella frazione per 10**55, possiamo vedere il valore fino a 55 cifre decimali:

>>> 3602879701896397 * 10 ** 55 // 2 ** 55
1000000000000000055511151231257827021181583404541015625

il che significa che il numero esatto memorizzato nel computer è pari al valore decimale 0.1000000000000000055511151231257827021181583404541015625. Invece di visualizzare il valore decimale completo, molti linguaggi (compresi le vecchie versioni di Python), arrotondano il risultato a 17 cifre significative:

>>> format(0.1, '.17f')
'0.10000000000000001'

I moduli fractions e decimal rendono questi calcoli facili:

>>> from decimal import Decimal
>>> from fractions import Fraction

>>> Fraction.from_float(0.1)
Fraction(3602879701896397, 36028797018963968)

>>> (0.1).as_integer_ratio()
(3602879701896397, 36028797018963968)

>>> Decimal.from_float(0.1)
Decimal('0.1000000000000000055511151231257827021181583404541015625')

>>> format(Decimal.from_float(0.1), '.17')
'0.10000000000000001'