Ponto flutuante: por que 0,1 + 0,2 não dá 0,3

Ponto flutuante: por que 0,1 + 0,2 não dá 0,3

Você escreve um teste unitário para uma calculadora de frete. A conta é simples: 0,1 mais 0,2, e o teste espera 0,3. O teste falha com uma mensagem que parece erro de digitação, mas não é. O computador somou direito, só que o resultado que ele guardou não é o número que você escreveu no papel. Esta aula explica de onde vem essa diferença, como o padrão IEEE 754 organiza um número real em bits e por que a diferença aparece em cálculos que nada têm de exótico, como dinheiro, medição e radar.

A ideia é mais velha que o computador. Quem já tentou escrever 1/3 em decimal sabe que a conta nunca termina: 0,3333… continua para sempre, e a régua para no que cabe no papel. O computador tem um problema parecido, só que na base 2, e ele tem menos papel. A aula parte da régua para chegar aos bits, mostra as contas à mão e confere cada uma com código.

Pré-requisito: saber converter binário e hexadecimal e entender o que é um bit, o que está em Bits e binário. A ideia de número com sinal e de overflow, de Inteiros, complemento de dois e overflow, ajuda nos exemplos de limite, e o padrão de bits aparece de novo em Texto, ASCII, Unicode e UTF-8.

Mapa da aula

A aula começa pelo problema de 0,1 e pelo que a régua de decimal não consegue medir. Depois vem a notação científica em binário, que é a ideia por trás de todo número de ponto flutuante. Em seguida, o formato IEEE 754 é desmontado campo por campo, com o caso normal, o subnormal e os valores especiais. Os formatos menores (16 e 8 bits) vêm na sequência, e então a parte que mais machuca na prática: a perda de precisão, o dinheiro, a FPU dentro da CPU e dois bugs históricos com dados verificados. Os erros comuns, os exercícios com resposta, o glossário, as fontes e o que não foi confirmado fecham a aula.

1. Por que 0,1 não cabe em binário

Pense numa régua de marcas de um décimo. Você quer medir 1/3 de uma unidade e a régua só tem marcas de 0,1. Você coloca a régua e vê que o objeto passa um pouco da marca de 0,3. Pode dividir a marca em 0,01, depois em 0,001, mas 1/3 nunca cai exatamente numa marca, porque a conta nunca termina em base dez. Essa limitação vem da régua, e o objeto continua inteiro.

Em binário a régua é ainda mais limitada em um caso específico: as marcas são potências de dois, 0,5, 0,25, 0,125 e assim por diante. O número 0,5 cabe exato, 0,25 também, mas 0,1 não. Para escrever 0,1 em binário é preciso uma dízima: 0,0001100110011001100110011… com o padrão 0011 se repetindo para sempre. Qualquer computador que guarde um número em um número fixo de bits vai cortar essa sequência em algum ponto, e o valor guardado é o mais próximo que cabe.

Essa é a origem do problema. Quando você escreve 0.1 numa linguagem comum, o computador não guarda o número que você digitou, e sim o vizinho binário mais próximo. O programa abaixo mostra o valor exato que fica guardado em cada número de 64 bits, com todas as casas decimais, e compara os dois vizinhos da soma 0,1 + 0,2.

from fractions import Fraction
from decimal import Decimal, getcontext
import struct

getcontext().prec = 80

def double_de_bits(b):
    return struct.unpack('>d', b.to_bytes(8, 'big'))[0]

def bits_de_double(x):
    return int.from_bytes(struct.pack('>d', x), 'big')

def decimal_exato(f):
    return Decimal(f.numerator) / Decimal(f.denominator)

soma = 0.1 + 0.2
print('0.1 + 0.2 =', repr(soma))
print('bits de 0,1 + 0,2 =', hex(bits_de_double(soma)))
print('bits de 0,3       =', hex(bits_de_double(0.3)))

exata = Fraction(0.1) + Fraction(0.2)           # soma exata dos dois doubles guardados
baixo = double_de_bits(0x3FD3333333333333)      # vizinho de baixo (0x...33)
cima = double_de_bits(0x3FD3333333333334)       # vizinho de cima (0x...34)
print('soma exata       =', decimal_exato(exata))
print('vizinho 0x...33  =', Decimal(baixo))
print('vizinho 0x...34  =', Decimal(cima))
print('distância até 0x...33 =', float(exata - Fraction(baixo)))
print('distância até 0x...34 =', float(Fraction(cima) - exata))
print('empate exato?', exata - Fraction(baixo) == Fraction(cima) - exata)
print('último bit da mantissa: 0x33 ->', 0x33 & 1, '| 0x34 ->', 0x34 & 1)
print('0x33 e 0x34 diferem em', bin(0x33 ^ 0x34).count('1'), 'bits')

A saída do programa é esta:

0.1 + 0.2 = 0.30000000000000004
bits de 0,1 + 0,2 = 0x3fd3333333333334
bits de 0,3       = 0x3fd3333333333333
soma exata       = 0.3000000000000000166533453693773481063544750213623046875
vizinho 0x...33  = 0.299999999999999988897769753748434595763683319091796875
vizinho 0x...34  = 0.3000000000000000444089209850062616169452667236328125
distância até 0x...33 = 2.7755575615628914e-17
distância até 0x...34 = 2.7755575615628914e-17
empate exato? True
último bit da mantissa: 0x33 -> 1 | 0x34 -> 0
0x33 e 0x34 diferem em 3 bits

Repare que 0,1 guardado é um pouco maior que 0,1, e 0,3 guardado é um pouco menor que 0,3. Somados, os dois valores guardados dão 0,3000000000000000166533453693773481063544750213623046875. Esse número não cai em nenhum dos dois vizinhos de 0,3 em double. Ele cai exatamente no meio entre eles, a 2,7756 vezes 10 elevado a -17 de cada um. Quando o valor exato fica no meio, o padrão IEEE 754 escolhe o vizinho cuja última casa da mantissa é par. Esse vizinho é o 0x3FD3333333333334, que aparece como 0,30000000000000004. A soma fez o arredondamento previsto pelo padrão. O que surpreende é a regra de empate, aplicada a uma conta cujo resultado exato fica entre duas marcas.

A figura abaixo mostra a diferença nos bits. O valor de 0,3 em double é 0x3FD3333333333333, e o de 0,1 + 0,2 é 0x3FD3333333333334. Os oito bytes são iguais, exceto o último.

Olhe a figura de cima para baixo. Os dois quadros têm oito bytes cada, mostrados em ordem de leitura, com o byte mais significativo primeiro. Os bytes 3F, D3 e cinco 33 aparecem nos dois, e só o último muda: 0,3 termina em 00110011 (0x33) e a soma termina em 00110100 (0x34). Os dois valores são vizinhos, e a distância entre eles é de uma unidade da última casa da mantissa, o que se chama 1 ulp (unidade da última casa). Em bits, porém, a troca toca três posições, porque 0x33 e 0x34 diferem em três bits. O traço tracejado marca que os dois números não são iguais. É por isso que a aula fala de “último byte”: a diferença é pequena, mas o teste de igualdade não perdoa nada.

Esse é o ponto que mais confunde quem chega aqui. O erro não cresce por acumulação quando a soma é de dois números: cada número já chega com um pequeno desvio, e a soma herda os dois. Quando você soma muitos números, os desvios se acumulam e a diferença pode ficar visível em casas decimais que você esperava certas. A seção 5 mostra isso com código.

2. Notação científica em binário

A forma de guardar números reais em computador é uma notação científica em base dois. Você já usa essa notação em decimal: 6.022 vezes 10 elevado a 23 é uma forma compacta de escrever um número grande, com uma mantissa, que são os dígitos significativos, e um expoente, que diz onde está a vírgula. Um computador faz o mesmo com base dois. O número é escrito como uma mantissa entre 1 e 2, multiplicada por dois elevado a um expoente inteiro. Por exemplo, 0,1 é aproximadamente 1,6 vezes 2 elevado a -4, porque 2 elevado a -4 vale 0,0625 e 1,6 vezes 0,0625 dá 0,1.

Repare que a mantissa fica sempre entre 1 e 2 no caso normal. Isso é uma convenção, e ela tem uma vantagem: o primeiro bit da mantissa é sempre 1, então o padrão não precisa armazená-lo. Esse bit escondido é chamado de bit implícito, e você vai vê-lo de novo na seção 3.

A notação científica tem uma consequência que a régua de decimal não mostra. Num número grande, o espaçamento entre os valores que o computador consegue representar cresce junto com o número. Com 32 bits, o espaçamento perto de 1 é de cerca de um décimo de milionésimo, e perto de 1 milhão é de cerca de 0,06. A régua fica mais grossa conforme você anda para a direita, e a precisão relativa, que é o erro em relação ao tamanho do número, fica parecida em todas as faixas. É um pouco como o hábito de escrever medidas de distância: um erro de um milímetro é grave numa peça de relógio e irrelevante numa rodovia, mas o computador tenta manter o mesmo número de algarismos corretos em qualquer escala.

Para ver a ideia sem os 32 bits, imagine um formato de 6 bits com 1 de sinal, 3 de expoente e 2 de mantissa, com bias 3. O expoente é guardado como um número sem sinal: o valor real mais um deslocamento fixo, que é o bias. Com bias 3, o campo 011 representa o expoente real 0. O campo de expoente com todos os bits em 1 fica reservado para infinito e NaN, como a seção 3 detalha, e o campo zero marca os valores subnormais. Sobram sete campos de expoente, e cada um tem quatro mantissas possíveis (2 elevado a 2). São 7 vezes 4, ou 28 valores não negativos, e entre eles estão o zero e 27 valores positivos. A figura a seguir coloca os 27 positivos numa escala logarítmica. A escala linear, usada numa régua comum, concentraria quase todos os pontos no lado esquerdo, e a escala log mostra o padrão.

Olhe a figura da esquerda para a direita. O eixo é logarítmico, então a distância entre 0,25 e 0,5 é igual à distância entre 4 e 8, porque em ambos o valor dobra. Cada faixa de potência de dois, de 0,25 a 0,5, de 0,5 a 1, de 1 a 2 e assim por diante, tem quatro pontos, e os quatro são espaçados de forma parecida dentro da faixa. O 1 aparece marcado porque é representável sem erro em qualquer formato binário, e os demais pontos são os outros valores de cada faixa. A faixa abaixo de 0,25 tem três positivos, 0,0625, 0,125 e 0,1875. Eles são os valores subnormais, números com expoente zero, que a seção 3 explica. O zero não aparece porque a escala log não tem como mostrá-lo.

A régua é um bom modelo, mas não é o formato real. O formato real de 32 bits tem 2 elevado a 23 valores em cada faixa de potência de dois, e a figura a seguir mostra quanto espaço existe entre dois valores consecutivos ao longo dos expoentes.

O gráfico mostra uma reta em escala log, o que significa que o espaço entre valores consecutivos cresce em progressão geométrica. Cada vez que o expoente sobe 1, o espaço dobra. Os dois marcadores laranja dizem o essencial: perto de 1 (expoente 0) o espaço é 2 elevado a -23, que é 0,00000011920928955078125 em decimal, e perto de 1024 (expoente 10) o espaço é 2 elevado a -13, cerca de 0,000122. Isso quer dizer que um float32 perto de mil já não consegue distinguir valores separados por menos de um décimo de milésimo, e um número perto de 16 milhões não distingue valores separados por menos de dois.

3. O formato IEEE 754

A notação científica binária é uma ideia. O IEEE 754 é o padrão que diz como guardar essa ideia em bits, com tamanhos fixos, regras de arredondamento e valores especiais para casos que a matemática também trata à parte. O padrão existe porque, antes dele, cada fabricante tinha o seu formato, e um mesmo programa dava resultados diferentes em máquinas diferentes. Uma consequência prática é que o mesmo padrão de bits tem o mesmo valor em praticamente qualquer CPU moderna, e isso é o que permite comparar resultados entre máquinas.

Um número de ponto flutuante tem três campos. O primeiro é o sinal, um bit: 0 para positivo e 1 para negativo. O segundo é o expoente, que é o número de casas que a vírgula anda para a direita ou para a esquerda. O terceiro é a mantissa, que guarda os dígitos significativos. Um detalhe importante é que o expoente é guardado com um deslocamento chamado bias, em vez de ser guardado com sinal. Em float32 o bias é 127, então o expoente real é o campo menos 127. Isso permite comparar dois números positivos olhando apenas os bits de expoente e mantissa, na ordem em que aparecem, um truque que o hardware usa.

Os três formatos mais usados têm tamanhos diferentes e o mesmo desenho. O binary16 tem 1 bit de sinal, 5 de expoente com bias 15 e 10 de mantissa. O binary32, que o C chama de float e o Python não usa por padrão, tem 1 de sinal, 8 de expoente com bias 127 e 23 de mantissa. O binary64, que o C chama de double e o Python usa para todo float, tem 1 de sinal, 11 de expoente com bias 1023 e 52 de mantissa. Em todos eles o valor de um número normal é 1, ponto, mantissa, vezes 2 elevado a expoente menos bias.

Para ver o padrão na prática, decomponha o 0,1 em float32. O valor exato dos bits é 0x3DCCCCCD. O primeiro bit é 0, então o número é positivo. Os próximos oito bits formam o expoente, que é 01111011, ou 123 em decimal. Subtraindo o bias 127 dá menos 4, e 2 elevado a -4 vale 0,0625. A mantissa tem 23 bits, 10011001100110011001101, e com o bit implícito a mantissa fica 1,6000000238. Multiplicando, 1,6000000238 vezes 0,0625 dá 0,10000000149, que é o valor guardado. A figura mostra a mesma decomposição campo por campo.

Leia a figura da esquerda para a direita. O primeiro campo é o sinal, seguido do expoente, que tem oito bits. A figura quebra os bits em bytes, por isso o expoente aparece em dois campos: sete bits no primeiro byte e o último bit no segundo, junto com os sete primeiros bits da mantissa. Os 23 bits da mantissa ocupam o resto, e os dois últimos bytes aparecem em hexadecimal, CC e CD. O CD é 11001101, e o CC seria 11001100. A seção a seguir explica por que o último byte é CD e não CC.

O caminho inverso mostra de onde vêm esses bits. Para 0,1, escreva o número como 1,6 vezes 2 elevado a -4. A mantissa é a parte fracionária 0,6 em binário, obtida dobrando o número e tirando a parte inteira a cada passo: 0,6 vira 1,2 (bit 1), 0,2 vira 0,4 (bit 0), 0,4 vira 0,8 (bit 0), 0,8 vira 1,6 (bit 1), e o ciclo se repete. Os 23 primeiros bits são 1001 1001 1001 1001 1001 100. O bit 24 é 1, seguido de outros bits que não são todos zero, então o resto passa da metade do último bit guardado. Nesse caso o arredondamento para o mais próximo soma 1 ao último bit: 100 vira 101, e a mantissa termina em 1101, que é o byte CD. Se o resto fosse menor que metade, o byte seria CC. O valor guardado é um pouco maior que 0,1 por causa desse arredondamento para cima.

Os números normais cobrem quase todo o intervalo, mas o padrão reserva casos especiais. Quando o expoente é zero, o bit implícito deixa de existir e o número vira um subnormal: o valor é 0,mantissa vezes 2 elevado a -126, em float32. Os subnormais preenchem o espaço entre zero e o menor normal, e sem eles existiria um salto enorme entre zero e o menor número normal. Com eles, a régua continua com marcas cada vez mais próximas perto de zero, até chegar a zero. O menor subnormal positivo em float32 vale cerca de 1,4 vezes 10 elevado a -45, o menor normal vale cerca de 1,18 vezes 10 elevado a -38, e o maior finito vale cerca de 3,4028 vezes 10 elevado a 38.

Quando o expoente é todo 1, o padrão reserva dois valores especiais. Se a mantissa é zero, o número é infinito, positivo ou negativo, o que aparece em divisões por zero e em estouros. Se a mantissa é diferente de zero, o valor é NaN, do inglês not a number, usado para resultados indefinidos como zero dividido por zero ou a raiz de um número negativo. O NaN é a única coisa que não é igual a si mesma: math.nan == math.nan dá falso, e por isso o teste correto é math.isnan.

Zero também tem dois valores. Com o sinal positivo, o zero é todo bit zero, e com o sinal negativo é apenas o bit de sinal aceso. Eles são iguais na comparação, mas o sinal sobrevive a algumas operações. Pelo padrão IEEE 754, 1 dividido por menos zero dá menos infinito. O Python puro não segue isso e levanta ZeroDivisionError, como a saída do programa abaixo mostra.

O programa abaixo reúne o que foi dito nesta seção. A função decodifica32 usa frações exatas, para que a própria conta não introduza o erro que queremos medir. Ele decodifica 0,1 campo por campo, confere os casos especiais e mostra os limites do float32.

import struct, math
from fractions import Fraction

def bits32(x):
    return struct.unpack('>I', struct.pack('>f', x))[0]

def decodifica32(b):
    sinal = b >> 31
    expo = (b >> 23) & 0xFF
    mant = b & 0x7FFFFF
    if expo == 0xFF:
        return sinal, expo, mant, ('NaN' if mant else ('-inf' if sinal else '+inf'))
    if expo == 0:   # subnormal: sem bit implícito, expoente efetivo -126
        valor = Fraction(mant, 2**23) * Fraction(1, 2**126)
    else:           # normal: bit implícito 1, expoente real = campo - 127
        valor = (1 + Fraction(mant, 2**23)) * Fraction(2) ** (expo - 127)
    if sinal:
        valor = -valor
    return sinal, expo, mant, float(valor)

b = 0x3DCCCCCD
s, e, m, v = decodifica32(b)
print(f'bits: 0x{b:08X} = {b:032b}')
print(f'sinal={s}  expoente={e:08b} ({e}) -> {e} - 127 = {e - 127}  mantissa={m:023b}')
print('valor =', repr(v))
print('o struct devolve:', repr(struct.unpack('>f', struct.pack('>I', b))[0]))
print('bits de 0,1 em float32 =', hex(bits32(0.1)))

print('padrões crescentes para valores crescentes (positivos):')
for x in [0.1, 0.5, 1.0, 2.0, 1024.0]:
    print(f'  {x:<6} -> 0x{bits32(x):08X}')

print('casos especiais:')
for nome, x in [('-0.0', -0.0), ('+inf', math.inf), ('NaN', math.nan)]:
    print(f'  {nome:6} bits=0x{bits32(x):08X}  copysign={math.copysign(1, x)}')
print('  -0.0 == 0.0 ?', -0.0 == 0.0, '| NaN == NaN ?', math.nan == math.nan)
print('  1e308*10 =', 1e308 * 10, '| inf - inf =', math.inf - math.inf)
try:
    1 / -0.0
except ZeroDivisionError as ex:
    print('  1 / -0.0 em Python:', ex)

for nome, b in [('menor subnormal', 0x00000001), ('menor normal', 0x00800000), ('maior finito', 0x7F7FFFFF)]:
    print(f'{nome:16} 0x{b:08X} -> {decodifica32(b)[3]!r}')

A saída do programa é esta:

bits: 0x3DCCCCCD = 00111101110011001100110011001101
sinal=0  expoente=01111011 (123) -> 123 - 127 = -4  mantissa=10011001100110011001101
valor = 0.10000000149011612
o struct devolve: 0.10000000149011612
bits de 0,1 em float32 = 0x3dcccccd
padrões crescentes para valores crescentes (positivos):
  0.1    -> 0x3DCCCCCD
  0.5    -> 0x3F000000
  1.0    -> 0x3F800000
  2.0    -> 0x40000000
  1024.0 -> 0x44800000
casos especiais:
  -0.0   bits=0x80000000  copysign=-1.0
  +inf   bits=0x7F800000  copysign=1.0
  NaN    bits=0x7FC00000  copysign=1.0
  -0.0 == 0.0 ? True | NaN == NaN ? False
  1e308*10 = inf | inf - inf = nan
  1 / -0.0 em Python: float division by zero
menor subnormal  0x00000001 -> 1.401298464324817e-45
menor normal     0x00800000 -> 1.1754943508222875e-38
maior finito     0x7F7FFFFF -> 3.4028234663852886e+38

A figura seguinte mostra a ordem dos testes que um decodificador faz para classificar um float32. Os dois testes do expoente não se sobrepõem: um expoente não pode ser zero e 255 ao mesmo tempo, então a ordem entre eles não muda o resultado. O que sobra depois deles é o caso normal.

Leia a figura da esquerda para a direita. A partir do losango da esquerda, o expoente igual a zero leva ao teste da mantissa, que separa zero com sinal de subnormal. O expoente diferente de zero segue para o losango de cima ou de baixo conforme ele seja 255. Com expoente 255 e mantissa zero, o número é infinito, e com mantissa diferente de zero, é NaN. Os dois têm o mesmo expoente e se separam só pela mantissa: 1 dividido por zero dá infinito, e zero dividido por zero dá NaN.

Um mesmo número existe em formatos de tamanhos diferentes com bits diferentes, e o valor 1,0 é o melhor exemplo para ver o bias mudando. A figura seguinte compara o 1,0 nos três formatos. Em todos eles a mantissa é zero, porque 1,0 é exatamente o bit implícito sozinho. O que muda é o expoente guardado, que é o bias mais o expoente real, e aqui o expoente real é zero.

Leia a figura de cima para baixo. No binary16 o expoente guardado é 01111, que é 15, e o bias também é 15, então o valor real é zero e o número é 1,0, codificado como 0x3C00. No binary32 o expoente guardado é 127, com oito bits, e o resultado é 0x3F800000. No binary64 o expoente é 1023, com 11 bits, e o padrão dá 0x3FF0000000000000. Os três têm mantissa zero. A diferença está no expoente guardado e no tamanho dos campos, e como a figura quebra os bits em bytes, o expoente do binary32 aparece em dois campos.

Para encerrar a seção, vale notar o que o padrão promete. Ele define exatamente como as quatro operações básicas e a raiz quadrada devem arredondar, e por isso cada uma dessas operações dá o mesmo resultado arredondado em máquinas que seguem o padrão. Um programa inteiro só repete o resultado se a linguagem ou o compilador não reordenarem as contas nem usarem precisão extra, como a seção 9 detalha. O padrão não promete que o resultado é o número real exato, e a seção 5 mostra que ele não pode prometer isso.

4. Formatos menores: 16 e 8 bits

Nem todo cálculo precisa de 32 ou 64 bits. Redes neurais, por exemplo, processam bilhões de números, e guardar cada um com menos bits libera memória e largura de banda. Por isso surgiram formatos menores, que trocam precisão por espaço. O custo é que a régua fica mais grossa, e o número de algarismos corretos cai.

O float16 é o binary16 da seção anterior, com 5 bits de expoente e 10 de mantissa. O bfloat16 é um formato de 16 bits usado em aprendizado de máquina: ele usa 8 bits de expoente, como o float32, e 7 bits de mantissa. A origem do formato no Google é a informação usual da área, e não foi conferida em fonte nesta aula. Como o expoente é do mesmo tamanho, o bfloat16 alcança os mesmos números grandes e pequenos que o float32, só que com menos casas significativas. Essa escolha é útil porque o bfloat16 é, na prática, o float32 com os 16 bits de baixo descartados, e o expoente já está no lugar certo. A conversão pega os 16 bits de cima e arredonda pelo que sobrou, porque cortar sem arredondar muda o resultado. Para 0,1, cortar dá 0x3DCC, que vale 0,099609375, e o arredondamento para o mais próximo dá 0x3DCD, que vale 0,10009765625.

A figura abaixo mostra o número 0,1 nos dois formatos de 16 bits. Os dois guardam o mesmo número, mas com valores aproximados diferentes, porque a divisão dos bits é diferente.

Os dois quadros têm 16 bits, e o que muda é a fronteira entre expoente e mantissa. No float16 os cinco bits de expoente, 01011, dão 11 menos 15, que é menos 4. No bfloat16 são oito bits, 01111011, que dão 123 menos 127, também menos 4. O expoente é o mesmo valor real nos dois, mas o bfloat16 tem uma mantissa de 7 bits contra 10 do float16. Por isso o float16 guarda 0,1 com mais casas (0,0999755859375), e o bfloat16 com menos (0x3DCD, que vale 0,10009765625). Quem escolhe o formato escolhe onde quer perder precisão: o float16 perde no tamanho máximo, que é 65504, enquanto o bfloat16 mantém o alcance do float32 e perde na casa decimal.

O programa abaixo converte 0,1 e 3,14159 nos dois formatos de 16 bits e confere os limites do float16. Depois ele enumera os 256 bytes do FP8 para achar os maiores valores finitos e os códigos especiais de cada formato.

import struct

def bits32(x):
    return struct.unpack('>I', struct.pack('>f', x))[0]

def float32_de_bits(b):
    return struct.unpack('>f', struct.pack('>I', b))[0]

print('float16 de 0.1 =', struct.unpack('>e', struct.pack('>e', 0.1))[0])
print('maior float16 =', struct.unpack('>e', struct.pack('>H', 0x7BFF))[0])
try:
    struct.pack('>e', 65520.0)
except OverflowError as ex:
    print('float16 de 65520:', ex)

def bf16_trunca(x):
    return bits32(x) >> 16

def bf16_arredonda(x):
    b = bits32(x)
    return (b + 0x7FFF + ((b >> 16) & 1)) >> 16   # o bit 16 faz o empate ir para o par

for x in [0.1, 3.14159]:
    t, a = bf16_trunca(x), bf16_arredonda(x)
    print(f'bfloat16 de {x}: truncado 0x{t:04X} = {float32_de_bits(t << 16)!r} | arredondado 0x{a:04X} = {float32_de_bits(a << 16)!r}')
print('float32 de 3.14159 =', float32_de_bits(bits32(3.14159)))

def decodifica(b, bits_exp, bits_man, bias):
    s = b >> (bits_exp + bits_man)
    e = (b >> bits_man) & ((1 << bits_exp) - 1)
    m = b & ((1 << bits_man) - 1)
    tam = 1 << bits_man
    if e == 0:
        return (-1) ** s * (m / tam) * 2.0 ** (1 - bias)
    return (-1) ** s * (1 + m / tam) * 2.0 ** (e - bias)

for nome, eb, mb, bias in [('E4M3', 4, 3, 7), ('E5M2', 5, 2, 15)]:
    especiais = {}
    maximo = (0.0, None)
    for b in range(256):
        e = (b >> mb) & ((1 << eb) - 1)
        m = b & ((1 << mb) - 1)
        if e == (1 << eb) - 1:
            if nome == 'E4M3' and m == (1 << mb) - 1:
                especiais.setdefault('NaN', []).append(b)   # E4M3: só S.1111.111 é NaN
                continue
            if nome == 'E5M2':
                tipo = 'infinito' if m == 0 else 'NaN'
                especiais.setdefault(tipo, []).append(b)
                continue
        v = decodifica(b, eb, mb, bias)
        if v > maximo[0]:
            maximo = (v, b)
    print(f'{nome}: maior finito = {maximo[0]} no código 0x{maximo[1]:02X}')
    for tipo, codigos in especiais.items():
        print(f'  {tipo}: ' + ', '.join(f'0x{c:02X}' for c in codigos))
print('448 = 1,75 x 2^8 =', 1.75 * 2 ** 8, '| 57344 = 1,75 x 2^15 =', 1.75 * 2 ** 15)

A saída do programa é esta:

float16 de 0.1 = 0.0999755859375
maior float16 = 65504.0
float16 de 65520: float too large to pack with e format
bfloat16 de 0.1: truncado 0x3DCC = 0.099609375 | arredondado 0x3DCD = 0.10009765625
bfloat16 de 3.14159: truncado 0x4049 = 3.140625 | arredondado 0x4049 = 3.140625
float32 de 3.14159 = 3.141590118408203
E4M3: maior finito = 448.0 no código 0x7E
  NaN: 0x7F, 0xFF
E5M2: maior finito = 57344.0 no código 0x7B
  infinito: 0x7C, 0xFC
  NaN: 0x7D, 0x7E, 0x7F, 0xFD, 0xFE, 0xFF
448 = 1,75 x 2^8 = 448.0 | 57344 = 1,75 x 2^15 = 57344.0

Os formatos de 8 bits levam a mesma ideia um passo adiante. O padrão FP8, descrito no artigo arXiv 2209.05433 (2022), com autores ligados à NVIDIA, define dois formatos de um byte: o E4M3, com 4 bits de expoente e 3 de mantissa, e o E5M2, com 5 de expoente e 2 de mantissa. Segundo o artigo, o E4M3 usa bias 7 e tem como maior valor finito 448, enquanto o E5M2 usa bias 15 e tem como maior valor 57344. O E5M2 segue as convenções do IEEE 754 e guarda infinito e NaN. O E4M3 abre mão do infinito e reserva só a combinação S.1111.111 para NaN. Com isso, o expoente 1111, que num formato IEEE seria reservado, passa a valer para números. Sem essa mudança, o maior valor seria 240. Com ela, sobe para 448, e o artigo cita esse ganho como sete magnitudes a mais.

As contas dos dois máximos são curtas. No E4M3, o expoente 1111 vale 15, menos o bias 7 dá 8, e a mantissa 110 vale 1,75. Então 1,75 vezes 2 elevado a 8 é 448. No E5M2, o expoente 11110 vale 30, menos o bias 15 dá 15, e a mantissa 11 vale 1,75, o que dá 1,75 vezes 2 elevado a 15, ou 57344. A mantissa do E4M3 é 110 e não 111 porque S.1111.111 é o NaN, que não entra na faixa dos números.

A figura seguinte mostra os três limites no byte: o maior valor finito do E4M3, o maior do E5M2 e o código de infinito do E5M2.

Olhe a figura de cima para baixo. O primeiro quadro é o E4M3 com expoente 1111 e mantissa 110, que dão 448. O segundo é o código de NaN do E4M3, com mantissa 111. O terceiro é o E5M2 com expoente 11110 e mantissa 11, que dá 57344. O último é o código de infinito do E5M2, com expoente 11111 e mantissa zero, que só existe porque o E5M2 segue o padrão IEEE. Quem usa FP8 precisa lembrar que o E4M3 não tem esse código de infinito. Quando uma conta produz um valor finito maior que o máximo, o artigo descreve a saturação como o tratamento usual: o resultado vira o maior valor representável. Um modo não saturante fica como opção para quem precisa tratar o estouro com rigor. Já um infinito ou um NaN de um formato mais largo, convertido para o E4M3, vira NaN, porque o E4M3 não tem infinito.

5. Arredondamento, épsilon e perda de precisão

Até aqui o arredondamento apareceu como um detalhe. Ele é o centro da perda de precisão. Quando uma conta produz um valor que não cabe nos bits, o padrão escolhe o vizinho mais próximo, e em caso de empate escolhe o que tem a última casa par. Esse modo, chamado arredondamento para o par, evita que somas longas acumulem um viés para cima ou para baixo.

Para entender o efeito, use o número 1e16. Em double, a partir de 2 elevado a 53, que é cerca de 9 vezes 10 elevado a 15, o espaçamento entre valores consecutivos passa a ser 2. Então os vizinhos de 1e16 estão a 2 de distância, e o número 1e16 mais 1 cai exatamente no meio entre 1e16 e 1e16 mais 2. Pelo arredondamento para o par, o resultado é 1e16. O bit de mantissa de 1 não existe nesse intervalo, e a soma se perde. A figura mostra os bits.

Olhe a figura da esquerda para a direita. O quadro da esquerda é 1e16, com o último byte 00. O quadro de cima, depois da seta barrada, é 1e16 mais 1, e tem exatamente os mesmos bytes: a soma não mudou nada, e a seta barrada diz isso. O quadro de baixo é 1e16 mais 2, que é o vizinho real, com o último byte 01. A diferença entre 1e16 e o vizinho está no último bit, então um 1 que cai no meio não consegue aparecer. Quem escreve 1e16 + 1 - 1e16 em Python obtém 0,0, e 1e16 - 1e16 + 1 obtém 1,0, porque a ordem das contas decide quando o 1 é descartado.

O exemplo de 1e16 mostra um caso extremo. O mesmo efeito aparece em formatos menores de forma mais visível. Em float32, o número 2 elevado a 24, que é 16777216, tem espaçamento 2 na faixa seguinte, e 2 elevado a 24 mais 1 volta a 16777216. O conceito que resume isso é o épsilon da máquina. Ele é a distância entre 1 e o próximo valor representável acima de 1, e vale 2 elevado a -52 em double, cerca de 2,22 vezes 10 elevado a -16, e 2 elevado a -23 em float32, cerca de 1,19 vezes 10 elevado a -7. Em float32, o valor 1 mais 2 elevado a -24 volta a 1,0, porque esse acréscimo é a metade do épsilon e cai no arredondamento. Já 1 mais 0,75 vezes 2 elevado a -23 passa da metade do épsilon e arredonda para cima, dando 1,0000001192092896. Ou seja, somar a 1 um valor menor ou igual à metade do épsilon devolve 1, e um valor maior que a metade, mas menor que um épsilon, arredonda para 1 mais épsilon. O épsilon é útil para decidir se duas quantidades são iguais o bastante para o problema, e não para decidir se são iguais exatamente.

A perda mais traiçoeira é o cancelamento. Ela acontece quando você subtrai dois números quase iguais, e o resultado guarda só os bits que restaram da diferença, que são poucos. Um exemplo clássico é calcular sqrt(x + 1) menos sqrt(x) para um x grande. Os dois valores são quase iguais, e a subtração joga fora os algarismos iguais, deixando só o erro de arredondamento de cada raiz. A figura mostra o erro relativo dessa fórmula ingênua para vários valores de x.

Leia o gráfico da esquerda para a direita. Ele traz 13 valores calculados, de x igual a 100 até x igual a 10 elevado a 14, um para cada potência de dez. A linha liga esses pontos, então o trecho entre eles é interpolação, e não medição contínua. Para x igual a 100 o erro é cerca de 6,5 vezes 10 elevado a -15, e para x igual a 10 elevado a 14 ele chega a 5,8 vezes 10 elevado a -3. A curva não é uma reta, e há um recuo entre 10 elevado a 11 e 10 elevado a 12, cuja causa exata não foi investigada nesta aula. A forma certa de calcular a mesma diferença é transformá-la numa expressão sem subtração, 1 sobre a soma das raízes, que não cancela. Quando x é 10 elevado a 12, por exemplo, a forma ingênua dá um erro relativo perto de 7,6 vezes 10 elevado a -6.

A ordem das somas também muda o resultado. Somar 0,1 dez vezes num laço dá 0,9999999999999999, e a mesma conta com sum([0.1] * 10) dá 1,0 no Python 3.12, porque a partir dessa versão a função sum usa a soma compensada de Neumaier para floats, como anunciado nas notas de versão do Python 3.12 (gh-100425). Essa compensação guarda o erro de cada soma e o devolve no fim, e por isso ela é mais exata. A compensação reduz o erro, mas não o elimina: é uma técnica e não uma garantia de exatidão. O programa abaixo reúne as contas desta seção: épsilon, arredondamento, somas e a tabela de erro do cancelamento.

import struct, math, sys
from decimal import Decimal, getcontext

getcontext().prec = 60

def f32(x):
    return struct.unpack('>f', struct.pack('>f', x))[0]

eps64 = sys.float_info.epsilon
eps32 = 2.0 ** -23
print('épsilon double =', eps64, '= 2^-52:', eps64 == 2.0 ** -52)
print('épsilon float32 = 2^-23 =', eps32)
print('float32: 1 + 0,75 x eps =', repr(f32(1 + 0.75 * eps32)), '(passa da metade, arredonda para cima)')
print('float32: 1 + 0,5 x eps  =', repr(f32(1 + 0.5 * eps32)), '(empate, vai para o par)')
print('float32: 1 + 2^-24      =', repr(f32(1 + 2.0 ** -24)))
print('float32: 2^24 + 1       =', repr(f32(f32(2.0 ** 24) + 1)))
print('double: 1e16 + 1 =', 1e16 + 1, '| 1e16 + 2 =', 1e16 + 2)
print('double: 1e16 + 1 - 1e16 =', 1e16 + 1 - 1e16, '| 1e16 - 1e16 + 1 =', 1e16 - 1e16 + 1)

acumulado = 0.0
for _ in range(10):
    acumulado += 0.1
print('laço somando 0.1 dez vezes =', repr(acumulado))
print('sum([0.1] * 10) =', repr(sum([0.1] * 10)))
print('math.fsum([0.1] * 10) =', repr(math.fsum([0.1] * 10)))

def erro_relativo(x):
    ingenuo = math.sqrt(x + 1) - math.sqrt(x)
    X = Decimal(x)
    exato = 1 / (((X + 1).sqrt()) + X.sqrt())   # mesma diferença sem subtração, em 60 dígitos
    return abs(Decimal(ingenuo) - exato) / exato

for k in range(2, 15):
    x = 10.0 ** k
    print(f'x = 1e{k:<2}  erro relativo = {float(erro_relativo(x)):.6e}')

x = 1e12
print('sqrt(x+1)-sqrt(x) com x=1e12:', repr(math.sqrt(x + 1) - math.sqrt(x)))
print('1/(sqrt(x+1)+sqrt(x)) com x=1e12:', repr(1 / (math.sqrt(x + 1) + math.sqrt(x))))

A saída do programa é esta:

épsilon double = 2.220446049250313e-16 = 2^-52: True
épsilon float32 = 2^-23 = 1.1920928955078125e-07
float32: 1 + 0,75 x eps = 1.0000001192092896 (passa da metade, arredonda para cima)
float32: 1 + 0,5 x eps  = 1.0 (empate, vai para o par)
float32: 1 + 2^-24      = 1.0
float32: 2^24 + 1       = 16777216.0
double: 1e16 + 1 = 1e+16 | 1e16 + 2 = 1.0000000000000002e+16
double: 1e16 + 1 - 1e16 = 0.0 | 1e16 - 1e16 + 1 = 1.0
laço somando 0.1 dez vezes = 0.9999999999999999
sum([0.1] * 10) = 1.0
math.fsum([0.1] * 10) = 1.0
x = 1e2   erro relativo = 6.504280e-15
x = 1e3   erro relativo = 1.141119e-13
x = 1e4   erro relativo = 2.130466e-13
x = 1e5   erro relativo = 6.422379e-13
x = 1e6   erro relativo = 9.255885e-11
x = 1e7   erro relativo = 1.976994e-10
x = 1e8   erro relativo = 1.361766e-08
x = 1e9   erro relativo = 1.562450e-07
x = 1e10  erro relativo = 1.116631e-06
x = 1e11  erro relativo = 8.899265e-06
x = 1e12  erro relativo = 7.614494e-06
x = 1e13  erro relativo = 1.610906e-03
x = 1e14  erro relativo = 5.828381e-03
sqrt(x+1)-sqrt(x) com x=1e12: 5.00003807246685e-07
1/(sqrt(x+1)+sqrt(x)) com x=1e12: 4.999999999998749e-07

A regra prática que sai dessa seção é a seguinte. Some números de tamanhos parecidos, e prefira somar os menores antes dos maiores quando a faixa de valores é grande. Evite subtrair dois números quase iguais, e quando não houver escolha, reescreva a fórmula. E desconfie de qualquer igualdade entre floats que não seja feita com uma tolerância escolhida de propósito.

6. Por que não usar float para dinheiro

Dinheiro parece o caso mais simples para ponto flutuante, porque a conta é com centavos e a maioria das pessoas pensa em decimal. O problema é justamente esse: 0,1 e 0,01 não cabem exatos em binário, então qualquer valor monetário guardado como float carrega um pequeno erro desde o começo. Somar 0,01 um milhão de vezes em float dá 10000,000000171856, enquanto math.fsum dá 10000,0. Num extrato com muitas linhas, o erro se acumula e pode aparecer no total.

O caso típico é o preço de 19,99 reais multiplicado por 100 para virar centavos. Em float64, o resultado é 1998,9999999999998, e a conversão para inteiro, com int(), trunca para 1998. O cliente paga um centavo a menos, ou a conta fica um centavo fora. Arredondar com round() antes de converter resolve este caso, e quando o mesmo valor é tratado como texto ou como Decimal, a conta é exata e dá 1999 centavos. A figura abaixo mostra os dois caminhos lado a lado.

Os dois grupos têm a mesma sequência de operações, e a diferença está no tipo do número. No grupo de cima, o preço entra como float64, a multiplicação devolve 1998,9999999999998 e o int() joga fora a parte fracionária, que era quase um centavo inteiro. No grupo de baixo, o preço entra como texto ou Decimal, a multiplicação dá 1999,00 exato e o valor é guardado como 1999 centavos inteiros. O erro é o mesmo em qualquer linguagem que use float para decimal, e é por isso que sistemas financeiros guardam centavos como inteiros ou usam tipos decimais.

O programa abaixo mostra as contas no Python, com a mesma entrada, e a soma de 0,01 um milhão de vezes.

from decimal import Decimal
import math

print('0.1 * 3 =', 0.1 * 3)
print('19.99 * 100 =', 19.99 * 100)
print('int(19.99 * 100) =', int(19.99 * 100))
print('round(19.99 * 100) =', round(19.99 * 100))
print("Decimal('19.99') * 100 =", Decimal('19.99') * 100)
print('Decimal(19.99) * 100 =', Decimal(19.99) * 100)

soma = 0.0
for _ in range(1_000_000):
    soma += 0.01
print('1 milhão de vezes 0,01, somados em float:', soma)
print('math.fsum de 1 milhão de 0,01:', math.fsum([0.01] * 1_000_000))

A saída do programa é esta:

0.1 * 3 = 0.30000000000000004
19.99 * 100 = 1998.9999999999998
int(19.99 * 100) = 1998
round(19.99 * 100) = 1999
Decimal('19.99') * 100 = 1999.00
Decimal(19.99) * 100 = 1998.999999999999843680598133
1 milhão de vezes 0,01, somados em float: 10000.000000171856
math.fsum de 1 milhão de 0,01: 10000.0

A regra é clara, e a exceção confirma. Valores de medição, como temperatura, distância ou o resultado de um modelo estatístico, aceitam float porque já têm incerteza muito maior que o erro de representação. Valores exatos, como dinheiro, contagem de estoque ou qualquer coisa que precisa bater centavo a centavo, não aceitam. Quando a conta é em centavos inteiros, o erro some: 1999 mais 1 é 2000, sem discussão.

7. A FPU dentro da CPU

Toda CPU de uso geral sabe fazer conta de ponto flutuante, mas não com a mesma unidade que faz as somas de inteiros. A ULA (unidade lógica e aritmética) cuida dos inteiros, e a unidade de ponto flutuante, chamada FPU, é um bloco separado, com seus próprios registradores e seus próprios circuitos de soma, multiplicação e divisão. Essa separação costuma ser explicada pelo fato de que os circuitos de ponto flutuante são bem maiores: alinhar expoentes, normalizar a mantissa e arredondar exige mais portas que uma soma de inteiros. Esse tamanho relativo não foi medido nesta aula.

No x86-64 moderno, os números de ponto flutuante são guardados em registradores de 128 bits, chamados xmm, que aparecem como xmm0 até xmm15 no modo SSE. Uma instrução SSE opera sobre um ou mais floats desses registradores, e o resultado volta para outro registrador. A FPU antiga, herdada do coprocessador 8087, usava uma pilha própria e ainda existe por compatibilidade, mas o código novo usa SSE. A figura mostra como os blocos se conectam dentro de um núcleo.

Siga a figura pelo caminho de baixo. Um float chega da memória e é carregado num registrador de ponto flutuante, o bloco verde. Desse registrador os operandos vão para a unidade de ponto flutuante, que devolve o resultado para o mesmo registrador. O caminho de cima é o dos inteiros, que tem o mesmo desenho com outra unidade. Os dois conjuntos de registradores se encontram só na conversão, a seta rosa: a ULA não faz conta de float, e a FPU não faz conta de endereço. Converter um inteiro em float exige mover o dado de um conjunto de registradores para o outro.

Uma consequência prática é o custo de misturar tipos. Um laço que mistura inteiros e floats a cada passo tende a gastar tempo movendo dados entre os dois conjuntos de registradores. O custo exato desse movimento não foi medido nesta aula. Compiladores tentam evitar isso, e o programador evita quando pode, mas a regra vale para qualquer cálculo em que o tipo muda muitas vezes.

8. Bugs reais: Pentium FDIV e Patriot

Os dois casos a seguir são de máquinas reais, e não são do mesmo tipo. O do Pentium é um defeito de hardware numa tabela da unidade de divisão. O do Patriot é um erro de precisão numa conversão de tempo, que, segundo o relatório do GAO, os registradores de 24 bits não permitiam evitar. Nenhum dos dois é um defeito do formato IEEE 754 em si, e ambos mostram que a precisão do número é só parte da história.

O primeiro é o defeito de divisão do processador Pentium, da Intel. O divisor do chip usava uma tabela de consulta com 2048 células, das quais 1066 deveriam ter um valor. Segundo a Wikipedia, cinco células que deveriam conter o valor +2 continham zero. No exemplo clássico, a divisão que passa por elas sai errada já na quarta casa decimal, mas esse é um caso particular. A mesma fonte diz que o pior caso aparece já no quinto algarismo significativo, e que o erro típico costuma aparecer no nono ou décimo. O professor de matemática Thomas Nicely, do Lynchburg College, descreveu o problema, e a Intel anunciou em 17 de janeiro de 1995 uma baixa contábil antes do imposto de 475 milhões de dólares para trocar os processadores afetados.

A conta que mostra o defeito é 4195835 dividido por 3145727. O resultado correto, com 14 algarismos, é 1,33382044913624. A Wikipedia dá o valor do chip defeituoso como 1,333739068902037589, que aparece aqui como 1,33373906890. O programa em C abaixo faz a divisão. Ele foi rodado num Mac com processador Apple, que tem arquitetura arm64 e não passa pelo divisor do Pentium, então o valor errado não aparece nessa máquina. O que o programa mostra é o resultado correto que o hardware atual dá.

/* Aula 1.3: o teste de divisão do Pentium. Num processador correto, o quociente é conhecido com precisão. */
#include <stdio.h>
int main(void) {
    volatile double x = 4195835.0, y = 3145727.0;   /* volatile impede o compilador de dobrar a conta */
    double q = x / y;
    printf("4195835 / 3145727 = %.18f\n", q);
    printf("resto do teste (x/y)*y - x = %.6f\n", (q * y) - x);
    return 0;
}

Compilado com cc -O0 fdiv.c, a saída é esta:

4195835 / 3145727 = 1.333820449136241093
resto do teste (x/y)*y - x = 0.000000

A figura abaixo mostra o caminho da divisão com a tabela completa e com a tabela defeituosa.

Leia a figura da esquerda para a direita. A divisão entra no bloco da esquerda, e dali há dois caminhos. O de cima consulta a tabela completa e chega ao quociente correto, em verde. O de baixo consulta a tabela com as cinco células zeradas e chega ao quociente errado, em rosa. O defeito está numa consulta que devia devolver um valor e devolveu zero, e não na aritmética de ponto flutuante em si. Por isso o defeito era difícil de achar: uma divisão só erra se o divisor cair numa das células zeradas, e os testes que não passam por elas dão a resposta certa.

O segundo caso é o do sistema de defesa antimísseis Patriot, usado pelo exército dos Estados Unidos na Guerra do Golfo. Um relatório do GAO, o órgão de auditoria do Congresso americano, de 4 de fevereiro de 1992, com o identificador B-247094 e o código IMTEC-92-26, descreve o que aconteceu em Dhahran, na Arábia Saudita. Em 25 de fevereiro de 1991, uma bateria Patriot que já operava havia mais de 100 horas seguidas não interceptou um míssil Scud, que atingiu um alojamento e matou 28 soldados americanos. O software modificado, que corrigia o erro, chegou a Dhahran em 26 de fevereiro, um dia depois.

O GAO atribui o erro ao relógio interno. O relatório diz que o tempo era convertido de inteiro para real, e que os registradores de 24 bits limitam a precisão dessa conversão. O relatório não traz a conta da deriva, nem o número de tiques, nem o valor de 0,1 truncado. Esses detalhes são reconstrução desta aula, feita a partir da tabela do Apêndice II. O efeito, pelo que o relatório mostra, é que após horas ligado o relógio atrasa o suficiente para a porta de alcance deixar de cobrir o alvo. A figura resume a tabela do Apêndice II.

A linha da figura sobe quase de forma reta, o que mostra que o erro cresce em proporção ao tempo ligado. A tabela do Apêndice II, lida no PDF do relatório, dá 7 metros após 1 hora, 55 após 8 horas, 137 após 20 horas, 330 após 48 horas, 494 após 72 horas e 687 após 100 horas. Para 100 horas, o relatório dá 0,3433 segundo de erro de tempo. O GAO também diz que a porta de alcance passa a errar quando o deslocamento chega a 50% da porta, o que acontece após cerca de 20 horas de operação contínua. O marcador de 20 horas na figura mostra esse limiar, com 137 metros. O relatório descreve o Scud como voando a cerca de Mach 5, ou 3750 mph, o que dá 1676 metros por segundo. Com essa velocidade, 0,3433 segundo de erro desloca o alvo em cerca de 575 metros, e não em 687. Para o relatório, 687 metros em 0,3433 segundo implicam cerca de 2001 metros por segundo. Essa diferença não foi reconciliada nesta aula.

A conta do erro de relógio é reconstrução desta aula, e não do GAO. Se 0,1 for truncado em 23 casas binárias depois da vírgula, o valor cortado é 0,00011001100110011001100, com erro de 9,5367 vezes 10 elevado a -8 segundo por tique. O relógio gera 36000 tiques por hora, e o erro dá 0,0034 segundo em 1 hora e 0,3433 segundo em 100 horas, como a tabela. Com 24 casas, o erro por tique cai para 3,5763 vezes 10 elevado a -8, e daria 0,1287 segundo em 100 horas, que não bate com a tabela. O registrador tem 24 bits, mas a reconstrução só bate com 23 casas fracionárias. A razão desse bit a menos não está no relatório e não foi confirmada nesta aula. O programa abaixo faz a conta com as duas quantidades de casas.

from fractions import Fraction
from math import floor

def truncado(x, casas):
    """Corta a fração x nas primeiras `casas` casas binárias depois da vírgula."""
    return Fraction(floor(x * 2 ** casas), 2 ** casas)

decimo = Fraction(1, 10)
print('0,1 truncado em 23 casas: 0.' + format(floor(decimo * 2 ** 23), '023b'))

for casas in (23, 24):
    erro_tique = decimo - truncado(decimo, casas)
    print(f'{casas} casas: erro por tique = {float(erro_tique):.4e} s')
    for horas in (1, 8, 20, 48, 72, 100):
        erro_s = erro_tique * 36000 * horas     # 36000 tiques de 0,1 s por hora
        print(f'   {horas:>3} h -> {float(erro_s):.4f} s')

velocidade_mach5 = 3750 * 0.44704                    # 3750 mph em metros por segundo
print('3750 mph =', round(velocidade_mach5, 1), 'm/s')
print('0,3433 s a essa velocidade =', round(velocidade_mach5 * 0.3433, 1), 'm')
print('687 m em 0,3433 s implicam', round(687 / 0.3433, 1), 'm/s')

A saída do programa é esta:

0,1 truncado em 23 casas: 0.00011001100110011001100
23 casas: erro por tique = 9.5367e-08 s
     1 h -> 0.0034 s
     8 h -> 0.0275 s
    20 h -> 0.0687 s
    48 h -> 0.1648 s
    72 h -> 0.2472 s
   100 h -> 0.3433 s
24 casas: erro por tique = 3.5763e-08 s
     1 h -> 0.0013 s
     8 h -> 0.0103 s
    20 h -> 0.0257 s
    48 h -> 0.0618 s
    72 h -> 0.0927 s
   100 h -> 0.1287 s
3750 mph = 1676.4 m/s
0,3433 s a essa velocidade = 575.5 m
687 m em 0,3433 s implicam 2001.2 m/s

9. Erros comuns e boas práticas

O primeiro erro é comparar floats com igualdade. Se o resultado de uma conta passa por uma soma, uma multiplicação ou uma divisão, a chance de ser exatamente igual a outro valor calculado de outra forma é pequena. Compare com uma tolerância escolhida para o problema, como abs(a - b) < 1e-9 * max(1, abs(a), abs(b)), ou use funções como math.isclose. A tolerância deve refletir a precisão real dos dados, e não um número arbitrário.

O segundo erro é tratar dinheiro como float. A solução é guardar centavos como inteiros, ou usar Decimal com contexto explícito, e a seção 6 mostra o motivo. O terceiro é subtrair dois números quase iguais sem perceber o cancelamento, como na seção 5. Quando a fórmula parece simples mas o resultado parece errado em valores grandes, procure a subtração perigosa antes de procurar um bug no compilador.

O quarto erro é somar muitos números de tamanhos muito diferentes na ordem em que aparecem. Quando um número grande entra no meio de uma soma de pequenos, os pequenos perdem bits antes de serem somados. Ordene por magnitude, use math.fsum ou sum do Python 3.12 quando a precisão importa, e teste com valores que estressam a conta.

O quinto erro é confiar em NaN para sinalizar erro sem testar. NaN não é igual a si mesmo e se propaga por qualquer operação. Um laço que produz NaN e continua a somar devolve NaN no fim, e o valor parece válido se ninguém verificar. Use math.isnan e trate o caso na entrada, e não na saída.

O sexto erro é assumir que o mesmo código dá o mesmo resultado em qualquer máquina. O padrão IEEE 754 garante o arredondamento das operações básicas, mas linguagens e compiladores podem reordenar contas, usar instruções com mais precisão intermediária ou trocar a função de biblioteca. Em teste de regressão numérica, fixe a precisão e o compilador, e compare com tolerância.

10. Exercícios

Os exercícios usam as contas da aula. Tente fazer à mão antes de conferir a resposta.

1. Decodifique o padrão 0x44800000 como float32. Qual é o sinal, o expoente real e o valor? Confira com struct.unpack('>f', bytes.fromhex('44800000')).

Resposta: o primeiro bit é 0, então o número é positivo. O campo de expoente é 137, e 137 menos 127 dá 10. A mantissa é zero, então o valor é 1 vezes 2 elevado a 10, ou 1024.

2. Por que 1e16 mais 1 dá 1e16 em double, mas 1e16 mais 2 dá um número diferente de 1e16?

Resposta: em double, 1e16 está entre 2 elevado a 53 e 2 elevado a 54, então o espaçamento entre valores vizinhos é 2. O número 1e16 mais 1 cai exatamente no meio entre 1e16 e 1e16 mais 2. O arredondamento para o par escolhe 1e16. O número 1e16 mais 2 é o vizinho exato seguinte, então ele é representável e diferente de 1e16.

3. O mini-float de 6 bits da seção 2 tem 1 bit de sinal, 3 de expoente com bias 3 e 2 de mantissa. Quantos valores não negativos finitos ele tem, e quantos são positivos?

Resposta: o expoente com todos os bits em 1 fica reservado para infinito e NaN, então sobram 7 valores de expoente, cada um com 4 mantissas. São 28 valores não negativos, contando o zero. Retirando o zero, sobram 27 positivos.

4. Calcule 19,99 vezes 100 de duas formas, uma com float e outra com Decimal('19.99'), e diga qual valor serve para guardar centavos.

Resposta: com float, o resultado é 1998,9999999999998, e int() dá 1998. Com Decimal, o resultado é 1999,00 exato. Para guardar centavos, o valor certo é o de Decimal, convertido para inteiro: 1999.

5. Em float32, o valor 1 mais 2 elevado a -24 é igual a 1? E 1 mais 2 elevado a -23?

Resposta: 2 elevado a -24 é metade do épsilon do float32, então o acréscimo cai no arredondamento e 1 mais 2 elevado a -24 dá 1,0. Já 1 mais 2 elevado a -23 é exatamente o épsilon, que é representável, e o resultado é 1,0000001192.

6. Codifique 0,75 em float32. Qual é o padrão de bits em hexadecimal? Confira com struct.pack('>f', 0.75).hex().

Resposta: 0,75 é 1,1 em binário vezes 2 elevado a -1, porque 1,1 em binário vale 1,5, e 1,5 vezes 0,5 dá 0,75. O expoente guardado é -1 mais 127, que dá 126, ou 01111110. A mantissa é 1 seguido de zeros, porque a parte fracionária 0,5 é o bit 1 logo depois da vírgula. O padrão é sinal 0, expoente 01111110 e mantissa 10000000000000000000000, ou seja, 0x3F400000.

Glossário

Ponto flutuante: forma de guardar um número real como sinal, expoente e mantissa, com a vírgula deslocada pelo expoente.

Mantissa: parte que guarda os dígitos significativos do número, com o bit implícito no caso normal. O IEEE 754 chama essa parte de significando.

Expoente: campo que diz a potência de dois do número, guardado com bias para ficar sem sinal.

Bias: valor somado ao expoente real para que o campo guardado seja sempre positivo. Vale 127 em float32, 1023 em float64 e 15 em float16.

Bit implícito: o 1 à esquerda da mantissa no caso normal, que o padrão não armazena.

Subnormal: número com expoente zero, sem bit implícito, que preenche o espaço entre zero e o menor normal.

Épsilon da máquina: distância entre 1 e o próximo valor representável acima de 1. Vale 2 elevado a -52 em double e 2 elevado a -23 em float32.

Ulp: unidade da última casa, a distância entre um valor representável e o vizinho seguinte, medida na última posição da mantissa.

Arredondamento para o par: regra do IEEE 754 que, em caso de empate, escolhe o valor cuja última casa é par.

Cancelamento: perda de precisão ao subtrair dois números quase iguais, porque os algarismos comuns somem.

NaN: valor de não número, usado para resultados indefinidos. Não é igual a si mesmo.

Soma compensada: técnica que guarda o erro de cada soma para devolvê-lo ao final, usada por sum em floats no Python 3.12.

bfloat16: formato de 16 bits com 8 bits de expoente e 7 de mantissa, o float32 com a parte fracionária encurtada. Usado em aprendizado de máquina.

E4M3 e E5M2: os dois formatos de 8 bits do FP8, com 4 e 5 bits de expoente, respectivamente.

ULA: unidade lógica e aritmética, o bloco do processador que faz as contas com inteiros.

FPU: unidade de ponto flutuante do processador, separada da ULA de inteiros, com seus próprios registradores.

Fontes

Consultadas em 1 de novembro de 2026.

• Notas de versão do Python 3.12, em docs.python.org/3.12/whatsnew/3.12.html, para a função sum com soma de Neumaier e o número gh-100425. O texto diz que sum() passou a usar soma de Neumaier para floats, com melhoria de exatidão.
• Artigo “FP8 Formats for Deep Learning”, arXiv 2209.05433, para os formatos E4M3 (bias 7, máximo 448) e E5M2 (bias 15, máximo 57344), para a regra de que o E5M2 segue o IEEE 754 com infinito e NaN, para o máximo de 240 sem a mudança do E4M3 e para a saturação e o modo não saturante descritos na seção de escala de tensores e na de tipos especiais.
• Relatório do GAO IMTEC-92-26, “Patriot Missile Defense: Software Problem Led to System Failure at Dhahran, Saudi Arabia”, com identificador B-247094 e data de 4 de fevereiro de 1992. O texto foi lido na camada de texto do PDF (gao.gov/assets/imtec-92-26.pdf) e conferido na página de apresentação (gao.gov/products/imtec-92-26). O Apêndice II traz a tabela de deslocamento, o texto traz a morte de 28 soldados, a velocidade de Mach 5 (3750 mph), o limiar de 50% após cerca de 20 horas, e a chegada do software em 26 de fevereiro de 1991.
• Verbete “Pentium FDIV bug” na Wikipedia em inglês, para as 2048 células, as 1066 preenchidas, as cinco células com zero, o valor de 1,333739068902037589, a baixa de 475 milhões de dólares anunciada em 17 de janeiro de 1995, e a informação sobre o algarismo em que o erro aparece. Thomas Nicely é citado como quem descreveu o problema.
• Os valores de bits e os limites de formato foram conferidos pelos scripts desta aula, rodados com Python 3.12.8 no macOS em 1 de novembro de 2026. O programa em C foi compilado com cc num Mac arm64. Essa conferência é a fonte principal desta aula para números.

O padrão IEEE 754 não foi lido diretamente nesta aula. As definições de campo, bias e casos especiais seguem a convenção conhecida do padrão e foram conferidas por decodificação, não pelo texto oficial.

O que não foi confirmado ou medido

O valor do Pentium defeituoso, 1,33373906890, vem da Wikipedia e não foi reproduzido nesta aula, porque o Mac usado tem processador arm64 sem o defeito. A reconstrução só cobre o caminho correto da divisão. A frase sobre o algarismo em que o erro costuma aparecer também vem da Wikipedia e não foi conferida em outra fonte.

No Patriot, a contagem de bits não foi reconciliada. O GAO fala em registradores de 24 bits, e a reconstrução do erro de relógio só reproduz a tabela do Apêndice II com 23 casas fracionárias. O mecanismo do 0,1 truncado e dos 36000 tiques por hora é reconstrução desta aula, e não está no relatório do GAO. A velocidade implícita de 2001 metros por segundo, para os 687 metros, também não bate com Mach 5, e essa diferença não foi explicada.

A data em que o defeito do Pentium foi descoberto varia entre as fontes, e a aula não a cita. O anúncio de 17 de janeiro de 1995 e a baixa de 475 milhões de dólares vêm da Wikipedia, e não de documento primário da Intel.

Os formatos FP8 foram descritos segundo o artigo, e não foram testados em hardware nem em biblioteca de aprendizado de máquina. A saturação e o modo não saturante também vêm do artigo, sem teste. A descrição do bfloat16 como corte dos 16 bits de baixo do float32, com arredondamento, foi conferida nos valores de 0,1 e de 3,14159, mas não em um acelerador real. A origem do bfloat16 no Google é a informação usual da área e não foi conferida em fonte nesta aula.

A afirmação de que, antes do IEEE 754, cada fabricante tinha o seu formato é a informação usual da área, sem fonte primária consultada nesta aula.

A descrição da FPU segue o x86-64 com SSE. Não foi aprofundada a pilha x87 do coprocessador antigo, nem a latência real das instruções em um processador específico. O custo de conversão entre inteiros e floats e o tamanho relativo da FPU não foram medidos; a aula os descreve só de forma qualitativa.

Os números de tamanho de faixa do float32, como 1,18 vezes 10 elevado a -38 para o menor normal e 3,4028 vezes 10 elevado a 38 para o maior, foram conferidos pelo struct do Python. A tabela de erro do cancelamento na figura 10 foi calculada com referência exata em decimal, e os valores são deste cálculo, não de um benchmark. A causa do recuo do erro entre 10 elevado a 11 e 10 elevado a 12 não foi investigada.

A figura 13 do Pentium e a figura 14 do Patriot são representações visuais de dados citados. Elas não reproduzem o hardware.

Relacionados

• Trilha de Hardware
• CPU
• Bits e binário
• Inteiros, complemento de dois e overflow
• Texto, ASCII, Unicode e UTF-8
• Álgebra booleana e portas lógicas