Precisione e accuratezza nei calcoli a virgola mobile

Numero KB originale: 125056

Riepilogo

Esistono molte situazioni in cui precisione, arrotondamento e accuratezza nei calcoli a virgola mobile possono funzionare per generare risultati sorprendenti per il programmatore. Devono seguire le quattro regole generali:

  1. In un calcolo con precisione singola e doppia, il risultato non sarà in genere più accurato della precisione singola. Se è necessaria una precisione doppia, assicurarsi che tutti i termini del calcolo, incluse le costanti, siano specificati con precisione doppia.

  2. Non presupporre mai che un valore numerico semplice sia rappresentato in modo accurato nel computer. La maggior parte dei valori a virgola mobile non può essere rappresentata con precisione come valore binario finito. Ad esempio, .1 è .0001100110011... in binario (si ripete per sempre), quindi non può essere rappresentato con precisione completa in un computer usando l'aritmetica binaria, che include tutti i PC.

  3. Non presupporre mai che il risultato sia accurato fino all'ultima posizione decimale. Esistono sempre piccole differenze tra la risposta "vera" e ciò che può essere calcolato con la precisione finita di qualsiasi unità di elaborazione a virgola mobile.

  4. Non confrontare mai due valori a virgola mobile per verificare se sono uguali o meno. Si tratta di un corollario della regola 3. Ci sono quasi sempre piccole differenze tra i numeri che "dovrebbero" essere uguali. Verificare invece sempre se i numeri sono quasi uguali. In altre parole, verificare se la differenza tra loro è piccola o insignificante.

Ulteriori informazioni

In generale, le regole descritte in precedenza si applicano a tutti i linguaggi, inclusi C, C++ e assembler. Gli esempi seguenti illustrano alcune delle regole che usano FORTRAN PowerStation. Tutti i campioni sono stati compilati usando FORTRAN PowerStation 32 senza alcuna opzione, ad eccezione dell'ultimo, scritto in C.

Esempio 1

Il primo esempio illustra due aspetti:

  • Le costanti FORTRAN sono a precisione singola per impostazione predefinita (le costanti C sono a precisione doppia per impostazione predefinita).
  • I calcoli che contengono termini di precisione singoli non sono molto più accurati dei calcoli in cui tutti i termini sono precisione singola.

Dopo essere stato inizializzato con 1.1 (una singola costante di precisione), y è impreciso come una singola variabile di precisione.

x = 1.100000000000000  y = 1.100000023841858

Il risultato della moltiplicazione di un singolo valore di precisione per un valore di precisione doppia accurato è quasi altrettanto negativo della moltiplicazione di due singoli valori di precisione. Entrambi i calcoli hanno migliaia di volte l'errore che moltiplica due valori di precisione doppia.

true = 1.320000000000000 (multiplying 2 double precision values)
y    = 1.320000052452087 (multiplying a double and a single)
z    = 1.320000081062318 (multiplying 2 single precision values)

Codice di esempio

C Compile options: none

       real*8 x,y,z
       x = 1.1D0
       y = 1.1
       print *, 'x =',x, 'y =', y
       y = 1.2 * x
       z = 1.2 * 1.1
       print *, x, y, z
       end

Esempio 2

L'esempio 2 usa l'equazione quadratica. Dimostra che anche i calcoli a precisione doppia non sono perfetti e che il risultato di un calcolo deve essere testato prima che dipenda dal fatto che piccoli errori possano avere risultati drastici. L'input alla funzione radice quadrata nell'esempio 2 è solo leggermente negativo, ma non è ancora valido. Se i calcoli a precisione doppia non hanno avuto errori lievi, il risultato sarebbe:

Root =   -1.1500000000

Genera invece l'errore seguente:

Errore di run-time M6201: MATH

  • sqrt: errore DOMAIN

Codice di esempio

C Compile options: none

       real*8 a,b,c,x,y
       a=1.0D0
       b=2.3D0
       c=1.322D0
       x = b**2
       y = 4*a*c
       print *,x,y,x-y
       print "(' Root =',F16.10)",(-b+dsqrt(x-y))/(2*a)
       end

Esempio 3

L'esempio 3 dimostra che a causa di ottimizzazioni che si verificano anche se l'ottimizzazione non è attivata, i valori possono mantenere temporaneamente una precisione superiore al previsto e che non è consigliabile testare due valori a virgola mobile per l'uguaglianza.

In questo esempio due valori sono entrambi uguali e non uguali. Al primo if, il valore di Z è ancora nello stack del coprocessore e ha la stessa precisione di Y. X non è quindi uguale a Y e il primo messaggio viene stampato. Al momento del secondo IF, Z doveva essere caricato dalla memoria e quindi aveva la stessa precisione e valore di X, e viene stampato anche il secondo messaggio.

Codice di esempio

C Compile options: none

       real*8 y
       y=27.1024D0
       x=27.1024
       z=y
       if (x.ne.z) then
         print *,'X does not equal Z'
       end if
       if (x.eq.z) then
         print *,'X equals Z'
       end if
       end

Esempio 4

La prima parte del codice di esempio 4 calcola la più piccola differenza possibile tra due numeri prossimi a 1,0. Questa operazione viene eseguita aggiungendo un singolo bit alla rappresentazione binaria di 1.0.

x   = 1.00000000000000000  (one bit more than 1.0)
y   = 1.00000000000000000  (exactly 1.0)
x-y =  .00000000000000022  (smallest possible difference)

Alcune versioni di FORTRAN arrotondano i numeri quando li visualizzano in modo che l'imprecisione numerica intrinseca non sia così ovvia. Questo è il motivo per cui x e y hanno lo stesso aspetto quando vengono visualizzati.

La seconda parte del codice di esempio 4 calcola la differenza più piccola possibile tra due numeri prossimi a 10,0. Anche in questo caso, questa operazione viene eseguita aggiungendo un singolo bit alla rappresentazione binaria di 10.0. Si noti che la differenza tra numeri prossimi a 10 è maggiore della differenza vicina a 1. Questo dimostra il principio generale che maggiore è il valore assoluto di un numero, minore è la precisione con cui può essere archiviato in un determinato numero di bit.

x   = 10.00000000000000000  (one bit more than 10.0)
y   = 10.00000000000000000  (exactly 10.0)
x-y =   .00000000000000178

Viene inoltre visualizzata la rappresentazione binaria di questi numeri per indicare che differiscono solo per 1 bit.

x = 4024000000000001 Hex
y = 4024000000000000 Hex

L'ultima parte del codice di esempio 4 mostra che i valori decimali semplici non ripetuti spesso possono essere rappresentati in binario solo da una frazione ripetuta. In questo caso x=1.05, che richiede un fattore ripetuto CCCCCCCC.... (Esadecimale) nella mantissa. In FORTRAN, l'ultima cifra "C" viene arrotondata per eccesso a "D" per mantenere la massima accuratezza possibile:

x = 3FF0CCCCCCCCCCCD (Hex representation of 1.05D0)

Anche dopo l'arrotondamento, il risultato non è perfettamente accurato. Si verifica un errore dopo la cifra meno significativa, che è possibile visualizzare rimuovendo la prima cifra.

x-1 = .05000000000000004

Codice di esempio

C Compile options: none

       IMPLICIT real*8 (A-Z)
       integer*4 i(2)
       real*8 x,y
       equivalence (i(1),x)

       x=1.
       y=x
       i(1)=i(1)+1
       print "(1x,'x  =',F20.17,'  y=',f20.17)", x,y
       print "(1x,'x-y=',F20.17)", x-y
       print *

       x=10.
       y=x
       i(1)=i(1)+1
       print "(1x,'x  =',F20.17,'  y=',f20.17)", x,y
       print "(1x,'x-y=',F20.17)", x-y
       print *
       print "(1x,'x  =',Z16,' Hex  y=',Z16,' Hex')", x,y
       print *

       x=1.05D0
       print "(1x,'x  =',F20.17)", x
       print "(1x,'x  =',Z16,' Hex')", x
       x=x-1
       print "(1x,'x-1=',F20.17)", x
       print *
       end

Esempio 5

In C le costanti mobili sono doppie per impostazione predefinita. Usare una "f" per indicare un valore float, come in "89.95f".

/* Compile options needed: none
*/ 

#include <stdio.h>

void main()
   {
      float floatvar;
      double doublevar;

/* Print double constant. */ 
      printf("89.95 = %f\n", 89.95);      // 89.95 = 89.950000

/* Printf float constant */ 
      printf("89.95 = %f\n", 89.95F);     // 89.95 = 89.949997

/*** Use double constant. ***/ 
      floatvar = 89.95;
      doublevar = 89.95;

printf("89.95 = %f\n", floatvar);   // 89.95 = 89.949997
      printf("89.95 = %lf\n", doublevar); // 89.95 = 89.950000

/*** Use float constant. ***/ 
      floatvar = 89.95f;
      doublevar = 89.95f;

printf("89.95 = %f\n", floatvar);   // 89.95 = 89.949997
      printf("89.95 = %lf\n", doublevar); // 89.95 = 89.949997
   }