41 Aritmética em Ponto Flutuante IEEE

Esse capítulo descreve funções para examinar a representação de números em ponto flutuante e controle do ambiente onde ocorrem os eventos envolvendo ponto flutuante de seu programa. As funções descritas nesse capítulo são declaradas no arquivo de cabeçalho ‘gsl_ieee_utils.h’.

41.1 Representação de números em ponto flutuante

O padrão IEEE para Aritmética em Ponto Flutuante no Formato Binário define formatos binários para números em precisão simples e em precisão dupla. Cada número é composto de três partes: um bit de sinal (s), um expoente (E) e uma fração (f). Os valores numéricos da combinação (s,E,f) é dado pela seguinte fórmula,

O formato para números em precisão simples usa 32 bits divididos da seguinte forma,

seeeeeeeefffffffffffffffffffffff
    
s = sign bit, 1 bit
e = exponent, 8 bits  (E_min=-126, E_max=127, bias=127)
f = fraction, 23 bits

O formato para precisão dupla usa 64 bits divididos da seguinte forma,

seeeeeeeeeeeffffffffffffffffffffffffffffffffffffffffffffffffffff

s = sign bit, 1 bit
e = exponent, 11 bits  (E_min=-1022, E_max=1023, bias=1023)
f = fraction, 52 bits

É muitas vezes útil estar apto a investigar o comportamento de um cálculo a nível de bit e a biblioteca fornece funções para mostrar as representações IEEE de uma forma legível a humanos.

Function: void gsl_ieee_fprintf_float (FILE * stream, const float * x)

Function: void gsl_ieee_fprintf_double (FILE * stream, const double * x)

Essas funções mostram uma versão formatada do número em ponto flutuante IEEE apontado por x para o fluxo stream. Um apontador é usado para informar o número indiretamente para evitar qualquer promoção indesejada de float para double. A saída toma uma das seguintes formas,

NaN

o símbolo de Não-um-Número (83)

Inf, -Inf

infinito positivo ou negativo

1.fffff...*2^E, -1.fffff...*2^E

um número em ponto flutuante normalizado

0.fffff...*2^E, -0.fffff...*2^E

um número em ponto flutuante denormalizado

0, -0

zero positivo ou negativo

A saída pode ser usada diretamente no modo GNU Emacs Calc precedendo essa saída com 2# para indicar binário.

Function: void gsl_ieee_printf_float (const float * x)

Function: void gsl_ieee_printf_double (const double * x)

Essas funções mostram uma versão formatada do número IEEE em ponto flutuante apontado por x para o fluxo stdout.

O seguinte programa demonstra o uso das funções mostrando as representações da precisão simples e da precisão dupla da fração 1/3. Para fins de comparação a representação do valor promovido de precisão simples para precisão dupla é também mostrado.

#include <stdio.h>
#include <gsl/gsl_ieee_utils.h>

int
main (void) 
{
  float f = 1.0/3.0;
  double d = 1.0/3.0;

  double fd = f; /* promote from float to double */
  
  printf (" f="); gsl_ieee_printf_float(&f); 
  printf ("\n");

  printf ("fd="); gsl_ieee_printf_double(&fd); 
  printf ("\n");

  printf (" d="); gsl_ieee_printf_double(&d); 
  printf ("\n");

  return 0;
}

A representação binária de 1/3 é 0.01010101... . a saída abaixo mostra que o formato IEEE normaliza essa fração para fornecer um dígito lider de 1,

 f= 1.01010101010101010101011*2^-2
fd= 1.0101010101010101010101100000000000000000000000000000*2^-2
 d= 1.0101010101010101010101010101010101010101010101010101*2^-2

A saída também mostra que um número em precisão simples é promovido a precisão dupla adicionando zeros na representação binária.

41.2 Ajustando seu Ambiente IEEE

O padrão IEEE define muitos modos para controle do comportamento de operações em ponto flutuante. esses modos especificam as importantes propriedades da aritmética computacional: a direção usada no arredondamento (e.g. se números devem ser arredondados para cima, para baixo ou para o número vizinho mais próximo), a precisão de arredondamento e como o programa deve tratar excessões aritméticas, tais como divisão por zero.

Muitos desses recursos podem ser agora controlados via funções padrão tais como fpsetround, que devem ser usadas sempre que estiverem disponíveis. Desafortunadamente no passado não houve uma API universal para controlar seu comportamento—cada sistema teve seu próprio caminho de baixo nível de tratar essas excessões. Para ajudar você a escrever programas portáveis a GSL permite especificar modos em uma forma que independe da plantaforma usando a variável de ambiente GSL_IEEE_MODE. A biblioteca então assume todas as inicializações específicas de máqina necessárias quando você chamar a função gsl_ieee_env_setup.

Function: void gsl_ieee_env_setup ()

Essa função lê a variável de ambiente GSL_IEEE_MODE e tenta ajustar o correspondente especificado modo IEEE. A variável de ambiente pode ser uma lista de palavras chave, separadas por vírgulas, como segue,

GSL_IEEE_MODE = "palavrachave,palavrachave,..."

onde palavrachave é um dos seguintes nomes de modo,

• single-precision
• double-precision
• extended-precision
• round-to-nearest
• round-down
• round-up
• round-to-zero
• mask-all
• mask-invalid
• mask-denormalized
• mask-division-by-zero
• mask-overflow
• mask-underflow
• trap-inexact
• trap-common

Se GSL_IEEE_MODE for vazia ou indefinida então a função retorna imediatamente e não tenta fazer mudanças no modo IEEE do sistema. Quando os modos da GSL_IEEE_MODE forem ajustados a função mostra uma mensagem curta mostrando os novos ajustes para relembrá-lo que os resultados do programa irão ser afetados.

Se o modo requisitado não for suportado pela plataforma sendo usada então a função chama o controlador de erro e retorna um código de erro de GSL_EUNSUP.

Quando opções forem especificadas usando esse método, o modo resultante é baseado em ajustes padronizados de maior precisão disponível (precisão dupla ou precisão extendida, dependendo da plataforma) no modo arredondar -para-o-mais-próximo, com todas as excessões habilitadas separadamente a partir da excessão INEXACT. A excessão INEXACT é gerada sempre que arredondamentos ocorrerem, de forma que a excessão INEXACT geralmente está desabilitada em cálculos científicos típicos. Todas as outras excessões em ponto flutuante são habilitadas por padrão, incluindo underflows e o uso de números denormalizados, por segurança. Essas excessões habilitadas podem ser desabilitadas com o ajuste individual mask- ou conjuntamente usando mask-all.

A seguinte ajustada combinação de modos é conveniente para muitos propósitos,

GSL_IEEE_MODE="double-precision,"\
                "mask-underflow,"\
                  "mask-denormalized"

Essa escolha ignora qualquer erros relativamente a pequenos números (ou denormalizados, ou underflowing para zero) mas captura overflows, divisão por zero e operações inválidas.

Note que na série x86 de processadores essa função ajusta ambos o modo original x87 e o mais recente modo MXCSR, que controla operações SSE em ponto flutuante. As unidades SSE em ponto flutuante não possuem um bit de controle de precisão, e sempre trabalha em precisão dupla. as palavras chave single-precision e extended-precision não causam efeito nesse caso.

Para demonstrar os efeitos de diferentes modos de arredondamento considere o seguinte programa que calcula e, a base dos logaritmos naturais, usando o somatório de uma série rapidamente decrescente,

#include <stdio.h>
#include <gsl/gsl_math.h>
#include <gsl/gsl_ieee_utils.h>

int
main (void)
{
  double x = 1, oldsum = 0, sum = 0; 
  int i = 0;

  gsl_ieee_env_setup (); /* read GSL_IEEE_MODE */

  do 
    {
      i++;
      
      oldsum = sum;
      sum += x;
      x = x / i;
      
      printf ("i=%2d sum=%.18f error=%g\n",
              i, sum, sum - M_E);

      if (i > 30)
         break;
    }  
  while (sum != oldsum);

  return 0;
}

Aqui está os resultados de executar o programa no modo round-to-nearest. Esse é o padrão IEEE de forma que não é realmente necessário especificar isso aqui,

$ GSL_IEEE_MODE="round-to-nearest" ./a.out 
i= 1 sum=1.000000000000000000 error=-1.71828
i= 2 sum=2.000000000000000000 error=-0.718282
....
i=18 sum=2.718281828459045535 error=4.44089e-16
i=19 sum=2.718281828459045535 error=4.44089e-16

Após dezenove termos o somatório converge para dentro de 4 \times 10^-16 do valor correto. Se agora mudarmos o modo de arredondamento para round-down o resultado final é menos preciso,

$ GSL_IEEE_MODE="round-down" ./a.out 
i= 1 sum=1.000000000000000000 error=-1.71828
....
i=19 sum=2.718281828459041094 error=-3.9968e-15

O resultado é sobre 4 \times 10^-15 abaixo o valor correto, uma ordem de magnitude pior que o resultado obtido no modo round-to-nearest.

Se mudarmos o modo de arredondamento para round-up então o resultado final é maior que o valor correto (quando adicionamos cada termo para o somatório o resultado final é sempre arredondado para cima,o que aumenta o somatório de pelo menos um tick até que o termo adicionado sofra um underflows para zero). Para evitar isso o problema irá precisar usar um critério de convergência segura, tai como while (fabs(sum - oldsum) > epsilon), com uma escolha adequada do valor de epsilon.

Finalmente podemos ver o efeito de calcular a soma usando arredondamento single-precision, no modo padrão round-to-nearest. Nesse caso o programa pensa que ainda está usando números em precisão dupla mas a CPU arredonda o resultado de cada operação em ponto flutuante para exatidão de single-precision. Isso simula o efeito de escrever o programa usando single-precision nas variáveis float ao invés de em variáveis double. A iteração para após cerca da metade do número de iterações e o resultado final é muito menos preciso,

$ GSL_IEEE_MODE="single-precision" ./a.out 
....
i=12 sum=2.718281984329223633 error=1.5587e-07

com um erro de O(10^-7), que corresponde a exatidão de precisão simples (cerca de 1 parte em 10^7). Continuar as iterações adicionais não diminui o erro pelo fato de todos os resultados subsequêntes arredondarem para o mesmo valor.

41.3 Referências e Leituras Adicionais

A referência para o padrão IEEE encontra-se em,

• ANSI/IEEE Std 754-1985, IEEE Standard for Binary Floating-Point Arithmetic.

Uma introdução mais pedagógica sobre o padrão pode ser encontrada no seguinte artigo,

• David Goldberg: What Every Computer Scientist Should Know About Floating-Point Arithmetic. ACM Computing Surveys, Vol. 23, No. 1 (March 1991), pages 5–48.

  Corrigendum: ACM Computing Surveys, Vol. 23, No. 3 (September 1991), page 413. and see also the sections by B. A. Wichmann and Charles B. Dunham in Surveyor’s Forum: “What Every Computer Scientist Should Know About Floating-Point Arithmetic”. ACM Computing Surveys, Vol. 24, No. 3 (September 1992), page 319.

Um livro texto detalhado sobre aritmética IEEE e seu uso prático está disponível da editora SIAM,

• Michael L. Overton, Numerical Computing with IEEE Floating Point Arithmetic, SIAM Press, ISBN 0898715717.

Esse documento foi gerado em 23 de Julho de 2013 usando texi2html 5.0.