В арифметике с плавающей запятой, почему числовая неточность возникает в результате добавления маленького члена к разнице больших членов?

13

Я читал книгу Аллена и Тилдесли « Компьютерное моделирование жидкостей ». Начиная со страницы 71, авторы обсуждают различные алгоритмы, которые используются для интеграции уравнений движения Ньютона в моделирование молекулярной динамики (МД). Начиная со страницы 78, авторы обсуждают алгоритм Верле, который, возможно, является каноническим алгоритмом интегрирования в MD. Они заявляют:

Возможно, наиболее широко используемый метод интегрирования уравнений движения - это метод, первоначально принятый Верле (1967) и приписанный Штормеру (Gear 1971). Этот метод является прямым решением уравнения второго порядка $m_i \ddot{\textbf{r}}_i = \textbf{f}_i$ . Метод основан на позициях $\textbf{r}(t)$ , ускорениях $\textbf{a}(t)$ и позициях $\textbf{r}(t - \delta t)$ из предыдущего шага. Уравнение для продвижения позиций выглядит следующим образом:

$\begin{matrix} (3.14) & r (t + δ t) = 2 r (t) - r (t - δ t) + δ t^{2} a (t) . \end{matrix}$ $\tag{3.14}\textbf{r}(t + \delta t) = 2\textbf{r}(t) - \textbf{r}(t - \delta t) + \delta t^2 \textbf{a}(t).$
Следует отметить несколько моментов, касающихся уравнения (3.14). Будет видно, что скорости не появляются вообще. Они были устранены путем сложения уравнений, полученных разложением Тейлора по $\textbf{r}(t)$ :

$r (t + δ t) = r (t) + δ t v (t) + (1 / 2) δ t^{2} a (t) + . . .$ $\textbf{r}(t + \delta t) = \textbf{r}(t) + \delta t \textbf{v}(t) + (1/2) \delta t^2 \textbf{a}(t) + ...$
$\begin{matrix} (3.15) & r (t - δ t) = r (t) - δ t v (t) + (1 / 2) δ t^{2} a (t) - . . . . \end{matrix}$ $\tag{3.15}\textbf{r}(t - \delta t) = \textbf{r}(t) - \delta t \textbf{v}(t) + (1/2) \delta t^2 \textbf{a}(t) - ... .$

Затем, позже (на странице 80), авторы утверждают:

Против алгоритма Верле, ... форма алгоритма может излишне вводить некоторую числовую неточность. Это возникает потому , что, в уравнении (3.14), маленький член ( ) добавляют к разности больших членов ( ), для того , чтобы генерировать траекторию. $\mathcal{O}(\delta t^2)$ $\mathcal{O}(\delta t^0)$

Я полагаю , что «малый термин» является , а «разность больших терминов» является . $\delta t^2 \textbf{a}(t)$ $2\textbf{r}(t) - \textbf{r}(t - \delta t)$

Мой вопрос: почему численная неточность является результатом добавления маленького термина к разнице больших?

Меня интересует довольно базовая концептуальная причина, поскольку я совсем не знаком с деталями арифметики с плавающей запятой. Кроме того, знаете ли вы какие-либо «обзорные» ссылки (книги, статьи или веб-сайты), которые познакомили бы меня с основными идеями арифметики с плавающей точкой, связанной с этим вопросом? Спасибо за ваше время.

floating-point computer-arithmetic Андрей
источник

9

Их наблюдение «форма алгоритма может без необходимости вносить некоторую неточность в числовом выражении» является правильным. Но их объяснение «» Это происходит потому, что в уравнении (3.14), небольшой срок ( ) добавляют к разности больших членов ( $O(δt^2)$ ), для тогочтобы генерировать траекторию. '' ложно. $O(δt^0)$

Истинная причина небольшой численной нестабильности алгоритма Верлета состоит в том, что он является лишь незначительно устойчивым, потому что разностное уравнение (по существу, случай, когда вы игнорируете в Verlet) имеет паразитное значение решение пропорционально , что приводит к линейному росту ошибок, вносимых в тогда как для полностью стабильного многошагового метода, применяемого к диссипативному дифференциальному уравнению, рост ошибок ограничен. $x_{k+1}=2x_k-x_{k-1}$ $a$ $k$ $k$

Изменить: Обратите внимание , что книга о численном моделировании молекулярной динамики, и для получения разумной точности полученных ожиданий нужно огромное число шагов, как точность весы с только. (Часто шаг по времени в пикосекундах следовать внутренней шкале колебаний. Но биологически соответствующие шкалы времени в миллисекундах или больше ( ), хотя , как правило , один не рассчитывает , что далек.) $N$ $O(N^{-1/2})$ $N\sim 10^9$

Для получения дополнительной информации см. Http://en.wikipedia.org/wiki/Linear_multistep_method#Stability_and_convergence

Арнольд Ноймайер
источник

10

Если вы ищете хорошее введение, я бы предложил статью Дэвида Голдберга « Что должен знать каждый компьютерщик об арифметике с плавающей точкой» . Это может быть слишком подробно, но это доступно онлайн бесплатно.

Если у вас есть хорошая библиотека, я бы предложил « Числовые вычисления Майкла Овертона» с арифметикой IEEE с плавающей точкой или первые несколько глав « Точность и стабильность числовых алгоритмов Ника Хайама» .

То, на что конкретно ссылаются Аллен и Тилдесли, это численное аннулирование . Суть в том, что если у вас есть, скажем, только три цифры, и вы вычитаете 100из 101, вы получите 1.00(в три цифры). Число выглядит так, как будто оно состоит из трех цифр, но на самом деле только первая цифра является истинной, а конечные .00- мусором. Почему? Ну, 100и 101являются только неточными представлениями, скажем, 100.12345и 101.4321, но вы можете хранить их только в виде трехзначных чисел.

Pedro
источник

δ t

$\delta t$

r (\t - δ t) \approx r (t)

$r(\t-\delta t)\approx r(t)$

r (t) = 1

$r(t)=1$ !

Arnold Neumaier

@ArnoldNeumaier: Да, пример Аллена и Тилдесли, кажется, не имеет особого смысла, я только хотел дать некоторую справку по проблемам, возникающим, когда «маленький термин [..] добавляется к разнице больших терминов», вот что ОП спросила, не является ли это проблемой в данном случае.

Педро

Но добавление маленького термина к большому - просто ошибка округления, ничего опасного. Отмена - это когда два почти равных больших члена вычитаются, чтобы получить крошечный термин. Это становится проблемой только тогда, когда вычитаемые промежуточные продукты намного больше, чем конечный результат вычисления, или когда небольшой промежуточный результат, на который влияет отмена, делится на другой маленький элемент.

Арнольд Ноймайер

@ArnoldNeumaier: Как я думаю, из моего ответа вполне очевидно, что я имел в виду проблему вычисления разности, а не суммы.

Педро

1

@ArnoldNeumaier: Точка занята, но я надеюсь, вы понимаете, что я считаю это довольно мелким для «-1».

Педро

5

Чтобы применить пример Педро к уравнению $(3.14)$ предположим, что ваши переменные хранятся со следующими значениями:

р (T) знак равно 101

$\mathbf{r}(t) = 101$

р (T - δ T) знак равно 100

$\mathbf{r}(t - \delta t) = 100$

δ T^{2} a (T) знак равно 1,49

$\delta t^2 \mathbf{a}(t) = 1.49$

Из $(3.14)$ это должно следовать за этим

р (T + δ T) знак равно 103,49

$\mathbf{r}(t + \delta t) = 103.49$

но, поскольку мы можем использовать только три цифры, результат становится усеченным до

р (t + δ t) = 103

$\mathbf{r}(t + \delta t) = 103$

This error will propagate, so that after 20 steps, assuming $\mathbf{a}(t)$ remains unchanged, you get $\mathbf{r}(t + 20\delta t) = 331$ instead of $433.90$ ,

Igor F.
источник

But the effect is that large only in 3-digit decimal arithmetic.

Arnold Neumaier

3

Pedro already gives the important fact, namely cancellation. The point is that every number you compute with has an associated accuracy; for example, a single precision floating point number can only represent things up to approximately 8 digits of accuracy. If you have two numbers that are almost exactly the same but differ in the 7th digit, then the difference will again be an 8-digit single precision floating point number and it looks like it is accurate to 8 digits, but in reality only the first 1 or 2 digits are accurate because the quantities you computed it from are not accurate beyond this first 1 or 2 digits of the difference.

Now, the book you cite is from 1989. Back then, computations were most often done in single precision and round-off and cancellation were serious problems. Today, most computations are done using double precision with 16 digits of accuracy, and this is far less a problem today than it used to be. I think it is worthwhile reading the paragraphs you cite with a grain of salt and take them in the context of their time.

Wolfgang Bangerth
источник

cancellation in double precision arithmetic can be as big a problem as in single precision. A case in point is Gaussian elimination without pivoting, which often gives very poor results due to cancellation, even in double precision.

Arnold Neumaier

-1: The Verlet formula typically retains all digits of accuracy, not just 1 or 2 of 8 in single precision.

Arnold Neumaier

@ArnoldNeumaier: Sure, you can get the same kind of problems in double precision. All I said is that one doesn't encounter them quite as frequently.

Wolfgang Bangerth

If you lose 6 digits three times in a chain of computations you have lost all digits even in double precision. Algorithms suffering from cancellation will usually be poor even in double precision. Verlet's algorithm is different since there is no cancellation but a mild linear growth of errors. Thus the loss of accuracy cannot multiply, making it suitable for much longer integration times. This was surely known to Allen & Tildesley.

Arnold Neumaier

Правильно. Но я имею в виду, что если у вас есть алгоритм без отмены, вы по-прежнему сталкиваетесь с ошибкой порядка 1e-8 с одинарной точностью, и если вы делаете 1e8 временных шагов, у вас могут возникнуть проблемы, даже если все остальное точно. 1e8 временных шагов - это порядок, который вы можете иметь для ODE. С другой стороны, с двойной точностью ваша погрешность на каждом шаге составляет 1e-16, и для полной потери точности потребуется 1e16 временных шагов. Это ряд шагов, которые вы не встретите на практике.

Вольфганг Бангерт

В арифметике с плавающей запятой, почему числовая неточность возникает в результате добавления маленького члена к разнице больших членов?

Ответы: