Na aritmética de ponto flutuante, por que a imprecisão numérica resulta da adição de um pequeno termo a uma diferença de grandes termos?

13

Eu tenho lido o livro Computer Simulation of Liquids de Allen e Tildesley. A partir da página 71, os autores discutem os vários algoritmos usados ​​para integrar as equações de movimento de Newton em simulações de dinâmica molecular (MD). A partir da página 78, os autores discutem o algoritmo Verlet, que talvez seja o algoritmo de integração canônica no MD. Eles afirmam:

Talvez o método mais amplamente utilizado para integrar as equações de movimento seja o inicialmente adotado por Verlet (1967) e atribuído a Stormer (Gear 1971). Este método é uma solução direta da equação de segunda ordem mir¨i=fi . O método é baseado nas posições r(t) , acelerações a(t) e nas posições da etapa anterior. A equação para avançar as posições é a seguinte:r(tδt)

(3.14)r(t+δt)=2r(t)r(tδt)+δt2a(t).

Há vários pontos a serem observados sobre o eqn (3.14). Será visto que as velocidades não aparecem de todo. Eles foram eliminados pela adição das equações obtidas pela expansão de Taylor sobre :r(t)

r(t+δt)=r(t)+δtv(t)+(1/2)δt2a(t)+...

(3.15)r(tδt)=r(t)δtv(t)+(1/2)δt2a(t)....

Posteriormente (na página 80), os autores declaram:

Contra o algoritmo Verlet, ... a forma do algoritmo pode desnecessariamente introduzir alguma imprecisão numérica. Isso ocorre porque, na eqn (3.14), um pequeno termo ( ) é adicionado a uma diferença de termos grandes ( ), para gerar a trajetória. O ( δ t 0 )O(δt2)O(δt0)

Eu acho que o "termo pequeno" é , e a "diferença de termos grandes" é .δt2a(t)2r(t)r(tδt)

Minha pergunta é: por que a imprecisão numérica resulta da adição de um termo pequeno a uma diferença de termos grandes?

Estou interessado em uma razão conceitual bastante básica, já que não estou familiarizada com detalhes da aritmética de ponto flutuante. Além disso, você conhece alguma referência do tipo "visão geral" (livros, artigos ou sites) que me introduzisse idéias básicas de aritmética de ponto flutuante relacionadas a esta pergunta? Obrigado pelo seu tempo.

Andrew
fonte

Respostas:

9

Sua observação "a forma do algoritmo pode desnecessariamente introduzir alguma imprecisão numérica" ​​está correta. Mas sua explicação '' Isso ocorre porque, na eqn (3.14), um pequeno termo ( ) é adicionado a uma diferença de termos grandes ( ), a fim de gerar a trajetória. '' é falso.O(δt2)O(δt0)

A verdadeira razão para a leve instabilidade numérica do algoritmo de Verlet é que ele é apenas marginalmente estável, porque a equação da diferença (essencialmente o caso em que você negligencia em Verlet) tem uma solução parasitária proporcional a , que faz com que os erros introduzidos cresçam linearmente em enquanto que para um método multipasso totalmente estável aplicado a uma equação diferencial dissipativa, o crescimento do erro é limitado.xk+1=2xkxk1akk

Edit: Note que o livro é sobre simulação numérica de dinâmica molecular, e para obter uma precisão razoável das expectativas resultantes é preciso um grande número de passos, como a balança de precisão com única . (Geralmente, o intervalo de tempo está nos picossegundos para seguir a escala de oscilação intrínseca. Mas as escalas de tempo biologicamente relevantes estão em milissegundos ou maiores ( ), embora geralmente não se calcule até agora.)NO(N1/2)N109

Para mais detalhes, consulte http://en.wikipedia.org/wiki/Linear_multistep_method#Stability_and_convergence

Arnold Neumaier
fonte
10

Se você está procurando uma boa introdução, sugiro o que todo cientista da computação deve saber sobre aritmética de ponto flutuante de David Goldberg . Pode ser um pouco detalhado, mas está disponível online gratuitamente.

Se você possui uma boa biblioteca, sugiro Computação Numérica de Michael Overton com Aritmética de Ponto Flutuante IEEE ou os primeiros capítulos da Precisão e Estabilidade de Algoritmos Numéricos de Nick Higham .

O que Allen e Tildesley estão se referindo especificamente é o cancelamento numérico . O resumo é que, se você tem, digamos, apenas três dígitos e subtrai 100de 101, recebe 1.00(em três dígitos). O número parece ter precisão de três dígitos, mas, na verdade, apenas o primeiro dígito é verdadeiro e o final .00é lixo. Por quê? Bem, 100e 101são apenas representações inexatas de, digamos, 100.12345e 101.4321, mas você só pode armazená-las como números de três dígitos.

Pedro
fonte
-1: onde está o cancelamento que você atribui à fórmula de Verlet? Tipicamente é pequeno, o que torna r ( \ t - δ t ) r ( t ) , sem o cancelamento resultante. Tente r ( t ) = 1 ! δtr(\tδt)r(t)r(t)=1
Arnold Neumaier
@ArnoldNeumaier: Sim, o exemplo de Allen e Tildesley não parece fazer muito sentido, eu só queria fornecer uma referência para os problemas que surgem quando "um pequeno termo [..] é adicionado a uma diferença de grandes termos", e é isso que o OP perguntou, não se é um problema no caso em questão.
Pedro Pedro
Mas adicionar um termo pequeno a um termo amplo é apenas um erro de arredondamento, nada perigoso. O cancelamento ocorre quando dois termos grandes quase iguais são subtraídos para obter um termo minúsculo. Isso se torna um problema apenas quando os intermediários subtraídos são muito maiores que o resultado final de uma computação ou quando o pequeno resultado intermediário afetado pelo cancelamento é dividido por outro pequeno elemento.
Arnold Neumaier
@ ArnoldNeumaier: Como, eu acho, é bastante óbvio pela minha resposta, eu estava me referindo ao problema de calcular a diferença, não a soma.
Pedro
1
@ ArnoldNeumaier: O ponto foi levantado, mas espero que você entenda que eu considero isso bastante mesquinho para um "-1".
Pedro Pedro
5

Para aplicar o exemplo de Pedro à equação , assuma que suas variáveis ​​estão armazenadas com os seguintes valores:(3,14)

r ( t - δ t ) = 100 δ t 2 a ( t ) = 1,49

r(t)=101
r(t-δt)=100
δt2uma(t)=1,49

De deve seguir-se que(3,14)

r(t+δt)=103,49

mas, como só podemos usar três dígitos, o resultado fica truncado para

r(t+δt)=103

Este erro será propagado, de modo que, após 20 etapas, assumindo que permaneça inalterado, você obtém r ( t + 20 δ t ) = 331 em vez de 433,90 ,uma(t)r(t+20δt)=331433,90

Igor F.
fonte
Mas o efeito é tão grande apenas na aritmética decimal de 3 dígitos.
Arnold Neumaier
3

Pedro já dá o fato importante, ou seja, o cancelamento. O ponto é que todo número com o qual você calcula tem uma precisão associada; por exemplo, um único número de ponto flutuante de precisão pode representar apenas itens com aproximadamente 8 dígitos de precisão. Se você tiver dois números quase exatamente iguais, mas diferentes no sétimo dígito, a diferença será novamente um número de ponto flutuante de precisão única de 8 dígitos e parece ter precisão de 8 dígitos, mas, na realidade, apenas o primeiro 1 ou 2 dígitos são precisos porque as quantidades que você calculou não são precisas além desses 1 ou 2 dígitos da diferença.

Agora, o livro que você cita é de 1989. Naquela época, os cálculos eram feitos com maior precisão e arredondamento e cancelamento eram problemas sérios. Hoje, a maioria dos cálculos é feita com precisão dupla com 16 dígitos de precisão, e isso é muito menos um problema hoje do que costumava ser. Eu acho que vale a pena ler os parágrafos que você cita com um grão de sal e levá-los no contexto de seu tempo.

Wolfgang Bangerth
fonte
o cancelamento na aritmética de precisão dupla pode ser um problema tão grande quanto na precisão única. Um caso em questão é a eliminação gaussiana sem pivotamento, que geralmente produz resultados muito ruins devido ao cancelamento, mesmo em dupla precisão.
Arnold Neumaier
-1: a fórmula Verlet normalmente mantém todos os dígitos da precisão, não apenas 1 ou 2 de 8 na precisão única.
Arnold Neumaier
@ ArnoldNeumaier: Claro, você pode obter o mesmo tipo de problemas em dupla precisão. Tudo o que eu disse é que não se encontra com tanta frequência.
Wolfgang Bangerth
Se você perder 6 dígitos três vezes em uma cadeia de cálculos, todos os dígitos serão perdidos, mesmo com precisão dupla. Algoritmos que sofrem de cancelamento geralmente são ruins, mesmo em dupla precisão. O algoritmo de Verlet é diferente, pois não cancelamento, mas um leve crescimento linear de erros. Portanto, a perda de precisão não pode se multiplicar, tornando-a adequada para tempos de integração muito mais longos. Isso certamente era conhecido por Allen & Tildesley.
Arnold Neumaier 04/04/12
Certo. Mas o que eu quero dizer é que, se você tiver um algoritmo sem cancelamento, ainda ocorrerá um erro da ordem de 1e-8 com precisão única e, se executar etapas de tempo 1e8, poderá ter um problema mesmo se tudo o mais for exato. 1e8 etapas de tempo é uma ordem de magnitude que você pode ter para a ODE. Por outro lado, com precisão dupla, sua imprecisão em cada etapa é 1e-16 e exigiria 1e16 etapas de tempo para obter uma perda completa de precisão. Essa é uma série de etapas que você não encontrará na prática.
Wolfgang Bangerth 5/09/12