Jump to content

Cuidado com a aritmética de ponto flutuante!


fredericopissarra

Recommended Posts

Posted

Diferente do que você pode achar, ponto-flutuante não funciona no domínio dos números reais (ℝ), mas sim no domínio dos racionais (ℚ). E isso faz uma grande diferença. Um valor em ponto flutuante é codificado para uma fração usando um número definido de bits (dependendo do tipo, se float, double, long double ou __float128). Para exemplificar, falarei do tipo float aqui (porque tem menos bits significativos), mas o texto aplica-se a todos os outros...

Um valor em ponto flutuante segue a formula:

png.latex.png.0222b95068a74d4d2d3e148db71fed5a.png

Onde s é um único bit (MSB), e tem 8 bits e f tem 23 bits. O valor f é inteiro e faz parte da estrutura do tipo float. Note a limitação de 23 bits. Isso significa que o valor máximo é 0x7fffff (ou 2²³-1)... O bit s nos dá o sinal do valor (1 = negativo) e o expoente e faz parte do fator de escala que faz o "ponto flutuar" para a esquerda ou direita (dai o nome).

Como o denominator é uma potência de 2, os únicos valores exatos que podem ser representados em float que tenham uma parte fracionária diferente de 0.0 (ou seja, o zero depois do "ponto" não conta), sempre terminará em 5 (deixo pro leitor tentar comprovar isso!)... Isso quer dizer que valores como 0.1, 0.2, 0.3, 0.4, 0.6, 0.7, 0.8 e 0.9 não podem ser representados exatamente - para citar exemplos com apenas um algarismo "decimal".

O valor 0.1, decimal, quando convertido para binário, resulta numa dizima periódica: 0.000110011001100110011... ad infinitum. Isso pode parecer estranho, mas considere a divisão de dois valores decimais inteiros exatos: 1/3. Algo similar acontece com valores binários.

Acontece que o valor codificado sempre terá 1.0 implícito (não está na estrutura do float), assim, o valor binário acima tem que sofrer o deslocamento do "ponto" em 4 bits para a direita e a "parte fracionária" tem que ter 23 bits... Ou seja, 1.10011001100110011001100*2¹²³⁻¹²⁷. No fim das contas teremos s=0, f=0x4ccccc e e=123. Acontece que f=0x4ccccd causa um erro menor (porque o próximo bit seria 1, mas ele não cabe no float).

Para exemplificar, eis um programinha que imprime a estrutura do float e o valor obtido à partir de 0.1f:

#include <stdio.h>

struct float_t {
  unsigned int f:23;
  unsigned int e:8;
  unsigned int s:1;
} __attribute__((packed));

void main( void )
{
  float v = 0.1f;
  struct float_t *vp = (struct float_t *)&v;

  printf("v=%.28f\n"
         "  s=%u, e=%u, f=%#x\n",
         v, vp->s, vp->e, vp->f);
}

Compilando e executando, obtemos isso:

$ cc -o float float.c
$ ./float
v=0.1000000014901161193847656250
  s=0, e=123, f=0x4ccccd

Repare nessa sequencia final de v... Repare no último algarismo (5, o zero à esquerda não conta, né?). E os valores de s, e e f são exatamente o que eu predisse, não é?

O tipo double tem f de 52 bits e e com 11 bits. O tipo long double tem f com 63 bits e e com 15 bits. O tipo __float128 (uma extensão do GCC?!) tem f com 112 bits e e também com 15 bits... Para saber o que isso significa em termos de precisão (a quantidade de "casas decimais" garantida) basta usar a formuletinha:

png.latex2.png.c4317f72f2f708a1fa129dd9ea172a8d.png

Isso dá uma ideia de grandeza do valor decimal, em termos de precisão (o log está ai porque queremos determinar a grandeza na base 10 com base no máximo possível na base 2, ou seja, 10^x = 2^y, onde sabemos y e queremos x). Aqui, Pdec é a precisão em algarismos decimais, Pbin é a quantidade de bits significativos (de f mais o bit implícito, no caso de float, double e __float128 - o tipo long double não tem um bit "implícito", mas ele também conta). Mas, para ficar mais fácil, a quantidade de algarismos decimais é aproximadamente igual a 30% da quantidade de bits significativos... Assim, um float possui 7 algarismos decimais de precisão, um double, 16 e um long double 19.

Note que falei "algarismos", não "casas depois da 'vírgula'"... Num valor como 3.1415926 o algarismo 6 final é incerto, já que temos 8 algarismos aqui... Aliás, este é o motivo porque printf imprime 6 "casas decimais depois da 'vírgula'" por default quanto usamos a coversão "%f", sem especificar a precisão (que tem um significado diferente da precisão do "ponto flutuante", como ficou claro, espero!)...

Mas, essa precisão decimal é enganosa. Como disse antes, a precisão é medida em quantidade de bits de f e do bit 1.0, implícito. E isso gera outro problema: O maior número inteiro que pode ser seguramente convertido para float só pode ter até 24 bits de tamanho! Considere o valor binário fracionário 1.111...111 ( com 23 bits fracionários). Se deslocarmos o ponto em 23 posições para a direita teremos 0xffffff, ou seja, 24 bits. O 25º bit, se houver um, simplesmente não cabe num float... Isso significa que algo como:

int x = 1234567891;
float f = x;

1234567891 é 0x499602D3. Usando apenas os 24 bits superiores à partir do MSB setado temos 1001001100101100000001011010011. Os bits em vermelho simplesmente não cabem num float e serão zerados. O valor original tem 31 bits significativos e o outro só pode ter 24, mas, como o primeirio 1 é implícito, teremos que deslocar o ponto em 30 posições para a direita (e=157) e obteremos 1.00100110010110000000101. Descartando o implícito temos f=0x132c05, mas como o bit seguinte seria 1 (o MSB vermelho), então f=0x132c06. No fim das contas, Isso nos dá um valor diferente para o float (1234567936, é só zerar os bits vermelhos acima e "arredondar" o último 0101 para 0110)! Use o programinha acima, se não acredita em mim...

A outra lição aqui é que, além dos valores em ponto flutuante serem uma representação no domínio dos racionais, ela quase sempre é um arredondamento para a fração mais próxima (o default é round to nearest). Mas, lembre-se: O numerador (f) só pode ter 23 bits (52 no double, 63 no long double e 112 no __float128).

Existem alguns outros pontos interessantes... quando e=255 a representação é a de um não-numero e existem dois: NaN (not-a-number) e infinitos. Uma tentativa de extrair a raiz quadrada de valores negativos, por exemplo, resulta num NaN... Uma divisão 0/0 resulta num NaN, mas uma divisão de um valor válido por zero resulta num +Inf ou -Inf, dependendo do sinal do numerador.

Quando e=0 temos uma classe de valores em ponto flutuante chamados de subnormais, onde o 1.0 implícito não existe. A especificação IEEE-754 chama isso de uma "maneira graciosa de underflow" e cria penalidades para o processador (operações aritméticas com subnormais são mais lentas!). Existe uma exceção: Se e=0 e f=0 então temos o valor 0.0... Mas, há algo estranho ai... O bit s nos dá o sinal e, portanto, existem +0 e -0.

Escrevi tudo isso (e é apenas um pedacinho) para pedir que evitem usar ponto flutuante sempre que possível... Eles não são uma panaceia para resolver problemas de precisão. Muitas vezes, buscar maneiras engenhosas de realizar operações com valores inteiros causam menos problemas...

Para mais informações sobre ponto flutuante, recomendo a leitura de What every computer scientist should know about floating point arithmetic.

[]s
Fred

Posted

Mais uma das consequências da estrutura do ponto-flutuante é que as propriedades comutativa e distributiva (leis essenciais da aritmética) são mais desobedecidas do que obedecidas. Um exemplo:

#include <stdio.h>

#define SHOW_RESULT(c) \
  printf("[%s]\n", ((c))?"yes":"no")

void testfp1 ( void )
{
  double x = 1.2;

  printf ( "x = 1.2 - 0.4 - 0.4 - 0.4; x = 0.0? " );

  x -= 0.4;
  x -= 0.4;
  x -= 0.4;

  /* x should be 0.0, right!? Wrong! */
  SHOW_RESULT ( x == 0.0 );
}

void testfp2 ( void )
{
  double x;
  double y;

  printf ( "x = (0.1 + 0.2) + 0.3; y = 0.1 + (0.2 + 0.3); x == y ? " );

  x = ( 0.1 + 0.2 ) + 0.3;
  y = 0.1 + ( 0.2 + 0.3 );

  /* x == y, right? Wrong! */
  SHOW_RESULT ( x == y );
}

void testfp3 ( void )
{
  double x;
  double y;

  printf ( "x = (0.1 * 0.2) * 0.3; y = 0.1 * (0.2 * 0.3); x == y? " );

  x = ( 0.1 * 0.2 ) * 0.3;
  y = 0.1 * ( 0.2 * 0.3 );

  /* x == y, right? Wrong! */
  SHOW_RESULT ( x == y );
}

void testfp4 ( void )
{
  double x;
  double y;

  printf ( "x = (0.1 + 0.2) * 0.3; y = (0.1 * 0.3) + (0.2 * 0.3); x == y? " );

  x = ( 0.1 + 0.2 ) * 0.3;
  y = ( 0.1 * 0.3 ) + ( 0.2 * 0.3 );

  /* x == y, right? Wrong! */
  SHOW_RESULT ( x == y );
}

int main ( void )
{
  testfp1();
  testfp2();
  testfp3();
  testfp4();

  return 0;
}

Compilando e rodando, obtemos:

$ cc -o fptest fptest.c
$ ./fptest
x = 1.2 - 0.4 - 0.4 - 0.4; x = 0.0? [no]
x = (0.1 + 0.2) + 0.3; y = 0.1 + (0.2 + 0.3); x == y ? [no]
x = (0.1 * 0.2) * 0.3; y = 0.1 * (0.2 * 0.3); x == y? [no]
x = (0.1 + 0.2) * 0.3; y = (0.1 * 0.3) + (0.2 * 0.3); x == y? [no]

 

Posted

Eis outro exemplo da quebra da lei comutativa:

#include <stdio.h>

void main( void )
{
  float a, b, c, r1, r2;

  a = 1e8; b = -1e8; c = 1;

  r1 = a + b + c;
  r2 = a + c + b;

  printf( "%f, %f\n", r1, r2 );
}

Compilando e rodando:

$ cc -o fptest2 fptest2.c
$ ./fptest2
1.000000 0.000000

Ué?! (a+b+c) é diferente de (a+c+b)? ?

PS: Mesmo habilitando otimizações dá na mesma por causa da "associatividade" dos operadores +.

Archived

This topic is now archived and is closed to further replies.

  • Recently Browsing   0 members

    • No registered users viewing this page.
×
×
  • Create New...