Calculando com o IEEE 754

As operações elementares, algoritmos e a propagação de erro. A falha do míssil Patriot e a explosão do Arianne 5. Os 10 mandamentos de uma computação segura.

Transcript

1. NumPy e o IEEE 754/2008 Prof. Paulo R. G. Bordoni UFRJ
2. Não bastam lápis e papel. Sem o poder de criação de um artista, nada como isto existiria! Qual?
3. http://www.leonardoda-vinci.org/ Ele foi cientista, matemático, inventor, ...
4. La Gioconda está no museu do Louvre, em Paris. Todas as pinturas e desenhos de Da Vinci aqui apresentados são cortesia da www.leonardoda-vince.org. Jamais existirá um sorriso mais enigmático que o dela.
5. L'Ultima Cena Afresco de Leonardo Da Vinci para a igreja de seu protetor, o Duque Lodovico Sforza. O trabalho se encontra no convento de Santa Maria delle Grazie, em Milão. Já vi muitas cópias desta. É A Última Ceia. Só não sabia que também era dele!
6. Uma besta gigante e uma metralhadora!
7. A anatomia dos homens e uma relação sexual!
8. Jovens, sugiro fortemente uma visita à Wikipedia para saber mais saber mais sobre a Renascença. Além de Da Vinci, conheçam mais sobre Michelangelo, Rafael Sanzio e outros fora da Itália. http://en.wikipedia.org/wiki/Renaissance_art
9. http://en.wikipedia.org/wiki/Renaissance_art Vou embarcar nessa viagem, Mestres!
10. Para criar essas “coisas” maravilhosas, não bastam pincéis e tintas! Sim Sherlock, além da criatividade, o próprio ato de pintar. Estamos falando da dinâmica do processo.
11. Nossa análise até agora mostrou que: Mais de 99,99% dos números no padrão IEEE, Single ou Double apresentam • erro absoluto da ordem da ULP • erro relativo da ordem de ε. Esse é o aspecto estático das coisas. E a dinâmica dos processos?
12. De fato Professora, a principal razão de existência dos números reside na possibilidade de computar com eles. “C’est la raison d’etre” dos computadores. Como se comportarão as operações algébricas perante tais erros?
13. Mestres. Há uma pergunta anterior a esta. Sim Filósofo, antes precisamos saber como são realizadas as operações elementares no computador.
14. As operações algébricas são funções definidas para números reais. A adição recebe um par de números reais e devolve um número real. Isto também vale para racionais. + Adição z = x + y ℝ (x,y) y x ℝ ℝ
15. ∗ Multiplicação O mesmo vale para a multiplicação. Quero enfatizar que elas são fechadas, no sentido que tanto a soma como o produto são números reais (idem para os racionais). (x,y) y x ℝ ℝ z = x ∗ y ℝ
16. É mesmo, Loirinha. Será que a soma de dois Single também será um Single? Será que elas também são operações fechadas para os números dos padrões IEEE 754?
17. A resposta é um sonoro NÃO. Mostre um contraexemplo no Toy 1 Professor . Vejam, a soma de x = 0.125 e y = 1.0, ambos números do Toy 1, é x+y = 1.125, que não está no Toy 1.
18. Bem Loirinha, o passo seguinte é obrigar “a ovelha desgarrada” retornar ao lar. Mas e aí, o que eu faço com o resultado?
19. Em outras palavras, precisamos transformar o resultado num número do tipo IEEE 754. Sim, e esse “modus operandi” se aplica a todos os tipos float do IEEE 754.
20. + Adição com IEEE 754 z = x + y ℚ (x,y) y x IEEE 754 IEEE 754 z* = fl(x + y) IEEE 754 fl A imagem a gravar é esta! Tanto para a adição como para a multiplicação.
21. Não Loirinha, quando somamos dois números, o computador já dá o resultado diretamente! Mas, depois de calcular x + y, eu vou ter que aproximar o resultado, isto é calcular float64(x+y) ?
22. Sim. Por exemplo nos casos em que a soma de dois números cai fora de ℝSingle ou de ℝDouble . Inclusive, o IEEE 754 já devolve informações do tipo +Inf, -Inf ou NaN nos casos excepcionais.
23. Mas, o que há de diferente entre as operações de adição por computador e a comum? A resposta é tão simples, que parece tola: para realizar a adição comum, “com lápis e papel”, (i.é, pela definição) precisamos “casar as vírgulas”.
24. 2345,6 + 1,023456 2346,623456 Por exemplo, para somar 2345,2 com 1,023456, basta “casar as vírgulas”, e somar dígito a dígito, levando em conta os “vai 1”, quando for o caso. Vejam:
25. Nesse caso, “casar as vírgulas”, envolve reajustar o menor dos expoentes e a fração correspondente. Entretanto, no computador, os números estão em ponto flutuante IEEE 754. 2345,6 = 2,3452 ∗ 10 3 1,023456 = 1,023456 ∗ 10 0
26. 2,3456 ∗ 10 3  2,3456 ∗ 10 3 1,023456 ∗ 10 0  0,001023456 ∗ 10 3 A ação: “casar as vírgulas”, envolve: 1. Igualar o menor dos expoentes ao maior. 2. Deslocar, para a direita, a fração do número de menor expoente, tantas “casas decimais” quantas for a diferença dos expoentes.
27. Percebam que: • Deslocar 3 casas para a direita é absolutamente igual a multiplicar por 10 -3. • Para compensar multiplicamos o número por 10 3. • O número continua o mesmo. 2,3456 ∗ 10 3  2,3456 ∗ 10 3 1,023456 ∗ 10 0  0,001023456 ∗ 10 3
28. 2,3456 ∗103 + 0,001023456 ∗103 2,345623456 ∗103 Bem, depois disso é só somar as frações:
29. Certo Professor, mas voltemos ao 1º exemplo. Tem um detalhe, Professora. No caso abaixo, precisamos efetuar um reajuste final, reposicionando a vírgula. 9,42135 ∗102 + 2,8102 ∗102 12,23155 ∗102 1,223155 ∗ 103
30. O aspecto crucial é que o campo da fração é limitado a 23 bits no Single. Algo entre 6 e 7 dígitos decimais. Então, parte da fração será perdida: 2,346623456  2,346623 2,3456 ∗103 + 0,001023456 *103 2,346623456 ∗103
31. Mestres, e o computador? Estamos fazendo um Cálculo Numérico teórico? Professor, pegue essa bola!
32. As operações de adição e multiplicação para números reais satisfazem as propriedades básicas: • x+(y+z) = (x+y)+z ....... Associativa • x+y = y+x .................... Comutativa • x + 0 = x ...................... Existência e unicidade do neutro • -x + x = 0 .................... Existência e unicidade do oposto Adição • x∗(y ∗ z) = (x ∗ y) ∗ z .... Associativa • x ∗ y = y ∗ x .................... Comutativa • x ∗ 1 = x ......................... Existência e unicidade do neutro • x ∗ x-1 = 1 ....................... Existência e unicidade do inverso, se x ≠ 0 Multiplicação
33. Tais questões são de importância capital porque tudo que os computadores sabem fazer são essas duas operações elementares. Sim Sherlock. Mais as operações lógicas, inclusive as de comparação. Será que as propriedade serão mantidas para os números IEEE 754?
34. No IEEE 754, tanto a adição como a multiplicação são comutativas. Além disso, pela própria construção, todo número x possui seu oposto -x. Basta trocar o bit do sinal. Uma análise rápida, mostrará que todas as outras propriedades falham. Vejam os contra-exemplos a seguir:
35. A associatividade da adição falha no computador. É o que mostra este programa.
36. No programa, eu sorteio valores para x e y no intervalo (0, 1) e testo se (x + y) – y = x Vejam alguns resultados!
37. Neste programinha, testamos se cada número x do IEE 754 possui seu inverso 1/x. É só testar a veracidade da igualdade x * 1/x = 1.
38. Nele, eu sorteio x no intervalo (0, 1), calculo 1/x e multiplico-os. Executei-o 3 vezes. Confira: em todas as execuções, para algum valor sorteado de x, o produto não é 1.
39. Eis aí uma consequência: Não vale a lei do cancelamento!
40. Este programa mostra que 1/x, não necessariamente, é o inverso de x. Nele faço x assumir, sucessivamente, os valores x = 0.1, 0.2, ..., 0.9. Em seguida calculo as potências xk e x- k e o produto xk ∗ x- k, para expoentes 1,2, 3 ... Quando xk ∗ x- k ≠ 1, interrompo.
41. Vejam, em muitas situações: u * 1/u ≠ 1. Sim, mas houveram potenciações também.
42. Cara, é impressionante! Só não falhou para 0.5 e 1. (que são potências de 2).
43. Problemas de precisão: Compare a 1ª soma com os três seguintes.
44. Novamente, problemas na precisão da soma. Agora em Single e em Half, para acentuar os erros.
45. Ah, Mestre, são errinhos pequenos, nem aparecem em Single! Mas, “de grão em grão, a galinha enche o papo!” Podemos ter graves consequências!
47. Outro problema grave: zero, não é o único número no IEEE 754 que satisfaz x + y = x. Vejam: Cruzes! x + y = x, com y≠0
48. É, problemas desse tipo aparecerão ao efetuarmos adições com números IEEE 754 que diferem muito em magnitude! Outro problema famoso com os números IEEE 754 é o “cancelamento catastrófico”.
49. No IEEE 754, quando dois números x e y são muito próximos, sua diferença x - y, fica muito poluída pelos erros de representação. Pode até não significar NADA! Quando os números são quase iguais suas frações, na representação binária IEEE 754, diferem apenas nos últimos dígitos (que normalmente já estão poluídos). Ao subtrairmos, eles se cancelam, e a diferença pode ser lixo.
50. Vejam, no float16 esses erros de cancelamento são grosseiros. Eles são evidentes até mesmo em float32!
51. Cruzes! x - y = 0, com x ≠ y Pois é Loirinha, imagine as consequências possíveis, em contas posteriores que dependam desse resultado!
52. Um “ERRINHO” na conversão float64 int16
54. Mestres, fiz este programa em Python, com o NumPy, para calcular impostos no meu estágio. As respostas dão erro!
55. Vejam! Meu Chefe riu, e mandou eu fazer com o Excel. Use float64 no lugar de float32!
56. Troquei em tudo. Quase funcionou! Pois é, outro dia, eu quase ganhei na Mega-sena. Quase!
57. Refiz o programa importando e usando o módulo decimal. Vejam:
58. Os resultados agora estão corretos. A desvantagem, face às planilhas eletrônicas é que a execução fica MUITO mais lenta.
59. Multiplicando, é muito mais rápido!
60. Com certeza Professora! Meu objetivo ao efetuar um milhão de adições foi mostrar o quanto a utilização do módulo decimal pode tornar tudo mais lento! Achei um truque sujo!
61. Agora precisamos analisar as operações lógicas de comparação. Pois é Filósofo, com todos essas falhas nas operações, aposto que também teremos problemas nas igualdades e desigualdades!
62. P/ex., no Toy 1, todos os números reais entre 1.0 e 1.125 são levados em 1. Por exemplo, fl(1.12) = 1.0 = fl(1.10). Pois é Professora, depreendemos daí que x > y ⇒ fl(x) ≥ fl(y)
63. Confira, Surfista, que para números reais, valem as igualdades z = w = x (para x > 0) No IEEE 754, isto nem sempre é verdade.
64. E o que eu faço, Mestra? Pulo da ponte Rio-Niterói ? Não querida, nem quando brigar com o namorado! Use comparações seguras, conforme mostraremos a seguir.
65. ULP(x)/2 x x* Loirinha, observe pela figura que: dado um número real x, sempre existirá um número x* do Toy 1 para o qual | x – x* | < ULP(x)/2 Sim Mestra, para x na faixa normal do Toy 1 e usando arredondamento.
66. ULP(x)/2 x x* = fl(x) Como, por definição, x* = fl(x), temos | x – fl(x) | < ULP(x)/2 para x na faixa normal do Toy 1 e usando arredondamento. E este é um resultado válido para qualquer sistema de ponto flutuante IEEE 754.
67. Portanto, se x e y são tais que fl(x) = fl(y), teremos: | x – y | ≤ | x – fl(x) | + | fl(x) – fl(y) | + | fl(y) – y | < ULP(x)/2 + 0 + ULP(y)/2 = ε∗[ 2**exp(x) + 2**exp(y) ]/2 < ε∗ 2**max { exp(x), exp(y) } = max { ULP(x), ULP(y) }
68. Portanto, a forma segura de efetuar testes de igualdade para floats do IEEE 754 é verificar a desigualdade | x – y | ≤ ULP(x,y), onde ULP(x,y) = max { ULP(x), ULP(y) }. E também que: • x – y > ULP(x,y) ⇒ x > y • y – x > ULP(x,y) ⇒ x < y
69. O “mesmo” programa, com comparação segura. Marquei também o erro na conclusão do outro.
70. Neste programa, também usamos a comparação segura.
71. Vejam alguns resultados. Usamos float16 para apenas para evidenciar mais facilmente os erros.
72. O “mesmo” programa sem a comparação segura. Assinalei a conclusão falsa
73. O programa calcula o Épsilon para IEEE Double”.
74. Eis seu resultado! Em Double, que é o padrão de Python. Eps
75. Os 10 mandamentos do IEE 754 Para utilizar eficientemente o computador como uma ferramenta para a Computação Científica, precisamos conhecer as limitações impostas pelo padrão IEEE 754. Algumas pecados podem ser evitados, outros devem apenas ser compreendidos e aceitos.
76. Utilizarás o Teorema Fundamental da Representação de Ponto Flutuante IEEE 754 para estimar a propagação do erro. 1º mandamento Sim Filósofo, faremos isto no futuro.
77. Pois é, já vimos que em Single, 99.99 % das conversões apresentam erro. Em Double muito mais. Lembrai-vos do erro de conversão de números decimais para binários e vice-versa (input/output). 2º mandamento
78. Não usarás mais que 7 dígitos decimais no padrão Single, nem 16 no Double. 3º mandamento Pois é, podemos pedir a apresentação mais dígitos na formatação de saída, mas a precisão é no máximo 7 (ou 16).
79. Depois do acidente do “Patriot”, jamais esquecerei Filósofo. Lembrarás sempre que as operações aritméticas básicas podem arruinar a precisão do resultado. 4º mandamento
80. Utilizarás float16 para evidenciar a propagação de erro. 5º mandamento Sim Filósofo, já utilizamos várias vezes nos exemplos.
81. Pois é, basta lembrar o exemplo do Ariane 5. Quero comentar ainda que: • Conversões de Single para Double não aumentam o número de bits significativos. • Conversões para tipos inteiros sempre truncam ( zero). Pode até acontecer de int(x) ≠ x quando x é um “inteiro em ponto flutuante”. Não efetuarás conversões de tipo impunemente. 6º mandamento
82. Garantir a segurança é fundamental. Efetuarás apenas comparações seguras 7º mandamento
83. Pois é, a Intel e a AMD podem ter implementado detalhes do padrão de forma diferente. Idem para acréscimos nas máquinas virtuais Java ou Python, que é herdeira de C padrão, que usa IEEE 754. Não acreditarás que implementações distintas do padrão IEEE 754 apresentarão os mesmos resultados. 8º mandamento
84. São mandamentos que envolvem novos conteúdos, a serem ainda apresentados! Aguardem ... 9º mandamento 10º mandamento
85. Antes de encerrar e, para que você não afirme, no futuro, que foi induzido a uma ideia errada ...
86. Tchau, até a próxima aula!