15. Aritmetică în virgulă mobilă: probleme și limitări

Numerele în virgulă mobilă (de la englezescul floating-point) sunt reprezentate în mașina fizică de calcul ca fracții scrise în baza de numerație 2 (fracții binare). Mai precis, fracția zecimală (adică, scrisă în baza de numerație 10) 0.625 are valoarea (exprimată cu puteri de exponent negativ ale lui 10) 6/10 + 2/100 + 5/1000 în timp ce, folosind același tip de exprimare cu puteri de exponent negativ, fracția binară 0.101 are valoarea 1/2 + 0/4 + 1/8. Acestor două fracții le corespunde, dpdv. numeric, aceeași valoare, singura deosebire dintre ele fiind aceea că prima a fost scrisă cu notații fracționare ale bazei de numerație 10 pe când cea de-a doua cu notații ale bazei de numerație 2.

Din păcate, cele mai multe fracții zecimale nu pot fi reprezentate exact sub formă de fracții binare. O consecință a acestui fapt este că, în general, numerele zecimale în virgulă mobilă pe care le tastați sunt doar aproximații ale numerelor binare în virgulă mobilă care vor fi stocate în mașina de calcul.

Dificultatea este mai ușor de înțeles dacă o abordăm pornind din baza 10. Să luăm în considerare fracția 1/3. O putem aproxima ca fracție în baza 10

0.3

sau, și mai bine,

0.33

sau, și mai bine,

0.333

șamd. Oricât de multe cifre ne-am strădui să tastăm, rezultatul nu va fi niciodată exact 1/3, chiar dacă vom ajunge la aproximații tot mai bune ale lui 1/3.

La fel, oricât de multe cifre în baza de numerație 2 ne-am strădui să întrebuințăm, valoarea zecimală 0.1 nu va putea fi reprezentată exact ca fracție în baza 2. Pentru că, în baza 2, 1/10 este fracția (periodică mixtă) cu număr infinit de cifre la dreapta punctului

0.0001100110011001100110011001100110011001100110011...

Dacă de oprim, din tastat, după indiferent câți biți (cifre binare), tot ce obținem este o aproximație. Pe majoritatea mașinilor de calcul din zilele noastre, numerele în virgulă mobilă (în englezește, ca jargon, ele sunt denumite float-uri) sunt aproximate folosind o fracție binară al cărei numărător utilizează primii 53 de biți, începând cu bitul cel mai semnificativ, iar al cărei numitor este o putere a lui 2. În cazul lui 1/10, fracția binară este 3602879701896397 / 2 ** 55, adică foarte aproape de valoarea exactă a lui 1/10 însă nu chiar egală cu ea.

Mulți utilizatori nu sesizează aproximația datorită modului în care valorile numerice sunt afișate. Astfel, Python-ul va printa doar o aproximație zecimală a valorii zecimale precise pe care o are aproximarea binară stocată în mașina de calcul. Pe cele mai multe din mașinile de calcul, dacă ar trebui să afișeze valoarea zecimală precisă a aproximării binare stocate a lui 0.1, atunci Python-ul ar avea de afișat:

>>> 0.1
0.1000000000000000055511151231257827021181583404541015625

Dat fiind că au apărut mai multe cifre decât le sunt de folos majorității utilizatorilor, Python-ul menține numărul cifrelor la valori ușor de manageriat afișând, în loc de cele de mai sus, o valoare rotunjită:

>>> 1 / 10
0.1

Să reținem, așadar, că, în pofida faptului că rezultatul afișat arată la fel ca valoarea exactă a lui 1/10, valoarea efectiv stocată este cea mai apropiată (de rezultat) fracție binară pe care mașina de calcul o poate reprezenta.

Ca o curiozitate, există mai multe numere zecimale (așadar, diferite) cu aceeași cea mai bună aproximație printr-o fracție binară. De exemplu, numerele 0.1 și 0.10000000000000001, precum și 0.1000000000000000055511151231257827021181583404541015625 sunt cu toatele aproximate de 3602879701896397 / 2 ** 55. Dat fiind că toate aceste valori zecimale au aceeași aproximație, oricare din ele ar putea fi afișată, atunci când se cere, păstrându-se invariantul eval(repr(x)) == x.

În trecut, funcția predefinită a Python-ului repr() ar fi ales numărul cu 17 cifre semnificative, și anume pe 0.10000000000000001. În zilele noastre, începând cu Python 3.1, Python-ul este capabil (pe majoritatea mașinilor de calcul) să aleagă numărul cel mai scurt și să afișeze, pur și simplu, 0.1.

Rețineți, deci, că situația descrisă mai sus ține de însăși natura numerelor binare în virgulă mobilă: cu alte cuvinte, ea nu este o eroare (de la englezescul bug) a Python-ului și nici o eroare a codului scris de dumneavoastră. Această complicație poate fi întâlnită în toate limbajele de programare care suportă aritmetica în virgulă mobilă a mașinii fizice de calcul (chiar dacă unele limbaje ar putea să nu publice diferențele în mod implicit, și nici în toate modurile de afișare).

Pentru o afișare plăcută ochiului, ați putea utiliza șirurile de formatare ca să produceți un număr restrâns de cifre semnificative:

>>> format(math.pi, '.12g')  # printează 12 cifre semnificative
'3.14159265359'

>>> format(math.pi, '.2f')   # printează 2 cifre la dreapta punctului
'3.14'

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

Merită să înțelegem că aceasta este, de fapt, o iluzie: noi nu facem decât să rotunjim afișajul valorii stocate în mașina de calcul.

Nicio iluzie nu vine de una singură. Astfel, cum 0.1 nu este chiar 1/10, nici suma a trei valori ale lui 0.1 nu va da chiar 0.3:

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

De asemeni, cum acest 0.1 nu poate ajunge oricât de aproape de valoarea exactă a lui 1/10 și nici 0.3-ul nu se poate apropia oricât de mult de valoarea exactă a lui 3/10, pre-rotunjirea valorilor bazată pe funcția round() nu ne va fi de niciun ajutor:

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

Chiar dacă numerele nu se pot apropia oricât de mult de doritele lor valori exacte, totuși, funcția math.isclose() ne poate fi de folos la compararea de valori inexacte:

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

Ca alternativă, funcția round() este utilă la compararea unor aproximații grosiere:

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

Aritmetica binară în virgulă mobilă este plină de surprize de acest fel. Dificultatea privitoare la „0.1” va fi explicată detaliat mai jos, în secțiunea „Erori de reprezentare”. Vedeți Exemple de complicații în virgulă mobilă pentru un sumar captivant al modului de funcționare al operațiilor în virgulă mobilă și al tipurilor de dificultăți întâlnite frecvent în practică. Vedeți, de asemeni, Pericolele virgulei mobile pentru o dare de seamă și mai cuprinzătoare a altor suprize des întâlnite.

Și, după cum se spune la finalul acestui din urmă material, „nu există răspunsuri simple.” Însă nu vă lăsați speriați prea ușor de virgula mobilă! Erorile produse de operațiile în virgulă mobilă din Python sunt moștenite de la arhitectura de calcul în virgulă mobilă a mașinii de calcul iar pe cele mai multe din mașini erorile nu vor depăși ordinul de 1 pe 2**53 per operație. Aceasta este mai mult decât potrivit pentru majoritatea activităților, dar nu scăpați din vedere că nu lucrați cu aritmetică zecimală și nici că orice operație în virgulă mobilă poate suferi de propria sa eroare de rotunjire.

Deși cazurile patologice vor continua să existe, utilizările obișnuite ale aritmeticii în virgulă mobilă vor conduce la rezultatul scontat atunci când veți rotunji afișarea rezultatelor finale la numărul de zecimale dorit. Chiar dacă funcția str() este, de obicei, suficientă, pentru un control mai fin al printării consultați specificatorii de format ai metodei str.format() din Format String Syntax.

Pentru cazurile de întrebuințare care au nevoie de reprezentări zecimale exacte, vă recomandăm să folosiți modulul decimal în care este implementată o aritmetică zecimală potrivită atât pentru aplicațiile de contabilitate cât și pentru cele de înaltă precizie (numerică).

Altă formă de aritmetică exactă este oferită de modulul fractions care implementează o aritmetică bazată pe numere raționale (astfel încât numerele de forma 1/3 să fie reprezentate exact).

În cazul în care sunteți un utilizator asiduu al operațiilor în virgulă mobilă atunci n-ar strica să aruncați o privire asupra pachetului NumPy ori asupra oricăruia din numeroasele pachete de calcul matematic și statistică oferite de către proiectul SciPy. Vezi <https://scipy.org>.

Python-ul vă pune la dispoziție unelte software care vă vor putea fi de folos în rarele ocazii în care vă va interesa cu adevărat valoarea exactă a unui număr în virgulă mobilă. Metoda float.as_integer_ratio() exprimă valoarea unui număr în virgulă mobilă sub formă de fracție:

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

Deoarece fracția este exactă, o putem folosi pentru a reconstitui, fără pierderi, valoarea originală:

>>> x == 3537115888337719 / 1125899906842624
True

Metoda float.hex() exprimă float-urile în format hexazecimal (adică, în baza de numerație 16), returnând, și ea, exact valoarea stocată în calculatorul dumneavoastră:

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

Această reprezentare hexazecimală precisă poate fi întrebuințată la reconstrucția fără pierderi a valorii float-ului:

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

Fiind o reprezentare exactă, ea este utilă la portarea (transportul, de la englezescul porting) în siguranță a valorilor între diferitele versiuni de Python (oferind, așadar, independența de platforma de calcul) precum și la schimbul de valori dintre Python și alte limbaje de programare care suportă acest format (cum sunt Java și C99).

Altă unealtă software folositoare este funcția sum() care știe cum să amelioreze pierderile de precizie de pe parcursul sumărilor. Ea se bazează pe precizia extinsă la rotunjirile pe care le efectuează în pașii de calcul intermediari în care valori succesive (termenii sumei) sunt acumulate într-o valoare totală (totalul cumulativ). Această tehnică va conta pentru acuratețea generală, împiedicând acumularea erorilor până la o valoare care să afecteze totalul general:

>>> 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

Funcția math.fsum() merge chiar mai departe și detectează toate „cifrele pierdute” de pe parcursul adăugirilor la totalul cumulativ, ceea ce conduce la o singură rotunjire a rezultatului final. Deși este mai înceată decât funcția sum(), ea se dovedește mai precisă în situațiile-limită, acolo unde date de valori mari se anihilează unele pe altele în calcul făcând ca grand totalul să fie aproape zero:

>>> valori_mari = [-0.10430216751806065, -266310978.67179024, 143401161448607.16,
...                -143401161400469.7, 266262841.31058735, -0.003244936839808227]
>>> float(sum(map(Fraction, valori_mari)))   # Sumare exactă cu o singură rotunjire
8.042173697819788e-13
>>> math.fsum(valori_mari)                   # O singură rotunjire
8.042173697819788e-13
>>> sum(valori_mari)                         # Rotunjiri multiple în precizie extinsă
8.042178034628478e-13
>>> totalul = 0.0
>>> for x in valori_mari:
...     totalul += x                         # Rotunjiri multiple în precizie obișnuită
...
>>>                                          # Adunarea directă nu conduce la nici
>>> totalul                                  # măcar o cifră corectă!
-0.0051575902860057365

15.1. Erori de reprezentare

Secțiunea de față se ocupă de exemplul „0.1” în detaliu și vă prezintă ce aveți de făcut pentru a realiza o asemenea analiză de caz riguroasă chiar dumneavoastră. Presupunem că dispuneți de o minimă familiaritate cu reprezentările binare în virgulă mobilă ale numerelor.

Eroarea de reprezentare se referă la faptul că anumite (de fapt, majoritatea lor) fracții zecimale nu pot fi reprezentate exact ca fracții binare (în baza 2). Această situație este motivul fundamental pentru care Python-ul (ori Perl-ul, C-ul, C++-ul, Java-ul, Fortran-ul și diverse alte limbaje de programare) nu afișează, în multe cazuri, exact numărul zecimal la care v-ați fi așteptat.

De ce asta? Pentru că 1/10 nu este reprezentabil exact ca fracție binară. Încă din anul 2000, cel puțin, aproape toate mașinile de calcul folosesc aritmetica binară în virgulă mobilă a standardului IEEE 754 iar majoritatea platformelor de calcul mapează float-urile din Python pe valorile binare, date în dubla precizie pe 64 de biți (binary64) a standardului IEEE 754. Valorile binary64 ale lui IEEE 754 conțin 53 de biți de precizie, astfel că, atunci când îl primește pe 0.1 ca dată de intrare, calculatorul va încerca să-l convertească în cea mai apropiată ca valoare fracție de forma J/2**N, unde J este un număr întreg de exact 53 de biți. Rescriind relația

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

drept

J ~= 2**N / 10

și ținând seama de faptul că J are exact 53 de biți (adică, este >= 2**52 dar < 2**53), deducem că valoarea cea mai bună pentru N este 56:

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

Cu alte cuvinte, 56 este singura valoare a lui N pentru care J va avea exact 53 de biți. Astfel, cea mai bună alegere a lui J va fi cea dată de rotunjirea câtului:

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

Dat fiind că restul (la împărțirea cu 10) este mai mare decât jumătate din 10, cea mai bună aproximație este cea obținută cu rotunjire prin adaos:

>>> q+1
7205759403792794

Așadar, cea mai bună aproximație a lui 1/10 în dubla precizie oferită de standardul IEEE 754 este:

7205759403792794 / 2 ** 56

Împărțind atât numărătorul cât și numitorul la 2 (adică, simplificând fracția cu 2), fracția devine:

3602879701896397 / 2 ** 55

Remarcați faptul că, deoarece am realizat o rotunjire prin adaos, această fracție este puțin mai mare decât 1/10; dacă nu am fi folosit adaosul, atunci câtul ar fi fost puțin mai mic decât 1/10. Însă, în niciun caz, nu am fi ajuns exact la 1/10!

În concluzie, calculatorul nu-l „vede” pe 1/10 niciodată : tot ce poate vedea este exact fracția dată mai sus, ea fiind cea mai bună aproximație în dublă precizie pe care o poate construi urmând standardul IEEE 754:

>>> 0.1 * 2 ** 55
3602879701896397.0

Dacă înmulțim fracția in cauză cu 10**55, atunci îi vom vedea valoarea scrisă cu 55 de cifre zecimale:

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

ceea ce înseamnă că numărul stocat în calculator este egal cu valoarea zecimală 0.1000000000000000055511151231257827021181583404541015625. În loc să afișeze această valoare zecimală în întregime, multe limbaje de programare (și, printre ele, versiunile mai vechi ale Python-ului), o vor rotunji la 17 cifre semnificative:

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

Modulele fractions și decimal ușurează aceste calcule:

>>> 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'