Изменения

Перейти к: навигация, поиск

Adaptive precision arithmetic

1061 байт добавлено, 20:33, 30 декабря 2011
Нет описания правки
{{В разработке}}
==Мотивация==
Все вычисления, производимые компьютером во <tex>floating-point</tex>[http://en.wikipedia.org/wiki/Floating_point]в модели с ''плавающей точкой'' модели, имеют погрешность. При большом количестве арифметических действий она возрастает. Во многих случаях результирующая погрешность уже не устраивает, и требуется либо абсолютно точное вычисление, либо меньшая погрешность. Одним из решений данной проблемы является хранение чисел в виде рациональных дробей, в которых числитель и знаменатель представляется в виде длинного целого числа. Но работать с такими числами довольно "дорого" по времени и тяжело в реализации: необходимо писать факторизацию чисел, эффективно сокращать дроби. Для улучшения работы нужны определенные оптимизации. Одной из них и является использование ''адаптивной арифметики'' (англ. ''adaptive precision arithmetic'').
==BackgroundБазовые понятия=====Представление чисел===Большинство современных процессоров поддерживают числа с плавающей точкой в форме <tex> \pm \mathrm{significand } \times 2^{\mathrm{exponent}}</tex>. Значащая часть числа (мантисса) представляет собой <tex>p</tex>-битное двоичное число в форме <tex>b.bbb \dots</tex>, где каждое <tex>b</tex>
обозначает один бит. Также имеется один бит на знак.
Далее в этой статье фраза "что-либо представимо в <tex>p</tex> битах" будет означать представимость в <tex>p</tex> битах мантиссы, ''не'' учитывая знак и экспоненту.
==Базовые понятия=====РасширенияРазложения===
{{Определение
|definition=
Ноль не пересекается ни с одним другим числом.
'''Например''', числа <tex>1100</tex> и <tex>-10.1</tex> не пересекаются, а <tex>110</tex> и <tex>10</tex> {{--- }} пересекаются.
{{Определение
}}
'''Например''', числа <tex>1100</tex> и <tex>11</tex> {{--- }} смежные, а <tex>1100</tex> и <tex>1000</tex> {{--- }} нет.
Иногда для использовании точной арифметики может понадобиться больше, чем <tex>p</tex> бит для хранения величин. В связи с этим вводится одно из базовых форм хранения чисел для такой арифметики.
{{Определение
|definition=
'''РасширениемРазложением''' числа <tex>x</tex> называется такое его представление <tex>x = x_n + x_{n-1} + \dots + x_1</tex>, где каждое <tex>x_i</tex> выражено <tex>p</tex>-битным числом с плавающей точкой и называется '''компонентой''' этого расширенияразложения.
}}
{{Определение
|definition=
Расширение Разложение называется '''неперекрывающимся (несмежным)''', если все его компоненты взаимно не перекрываются (не являются смежными).
}}
Как правило, расширения разложения должны быть неперекрывающимися, а их компоненты должны быть упорядочены от большей к меньшей по величине (то есть <tex>x_n</tex> {{- --}} большая). Далее будут рассматриваться именно такая их форма.
<wikitex>Стоит отметить, что число может быть представлено несколькими возможными неперекрывающимися расширениямиразложениями: $1100 + -10.1 = 1001 + 0.1 = 1000 + 1 + 0.1$.
</wikitex>
Неперекрывающиеся расширения разложения нужны, например, для того, чтобы быстро вычислять знак выражения (смотрим знак большей по размеру компоненты), или для грубой оценки значения всего расширения разложения (берем большую по величине компоненту).
===Округление===
{{Определение
|definition=
'''Точное округление''' (англ. ''exact rounding'') {{- --}} такой вид округления, что:
* если точный результат ''может'' быть представлен в <tex>p</tex> битах, то результатом округления будет точное значение числа;
* если точный результат ''не может'' быть представлен в <tex>p</tex> битах, то результатом округления будет ближайшее <tex>p</tex>-битное значение.
}}
'''Например''', в 4-битной арифметике произведение <tex>111 \times 101 = 100011</tex> будет округлено в <tex>1.001 \times 2^5</tex>.
При вычислении результата может возникнуть ситуация, когда значение попадает в точности между двумя соседними <tex>p</tex>-битными значениями. Тогда требуется определить правило поведения в таком случае. Рассмотрим некоторые из них.
{{Определение
|definition=
'''Округление до ближайшего четного''' (англ. ''round-to-even'') {{- --}} правило, при котором округление в вышеописанном случае производится к ближайшему <tex>p</tex>-битному четному значению.
}}
{{Определение
|definition=
'''Округление к нулю''' (англ. ''round-toward-zero'') {{- --}} правило, при котором округление в вышеописанном случае производится <tex>p</tex>-битному значению, находящемуся между точным значением и нулем, а также ближайшему к точному значению.
}}
'''Например''', в 4-битной арифметике число <tex>10011</tex> будет округлено до <tex>1.010 \times 2^4</tex> по первому правилу, и до <tex>1.001 \times 2^4</tex> по второму.
Стоит отметить, что стандарт IEEE 754 использует округление до ближайшего четного по умолчанию.
Далее в этой статье символами <tex>\oplus, \ominus</tex> и <tex>\otimes</tex> будут обозначаться <tex>p</tex>-битные сложение, вычитание и умножение с точным округлением соответственно.
Из-за округления данные операции теряют некоторые важные свойства, например, ассициативностьассоциативность: <tex>(1000 \oplus 0.011) \oplus 0.011 = 1000</tex>, но <tex>1000 \oplus(0.011 \oplus 0.011) = 1001</tex>.
При анализе округления часто используют так называемый <tex>\operatorname{ulp}</tex>.
{{Определение
|definition=
'''ULP''' (англ. ''units in the last place'') {{- --}} эффективная величина самого младшего (<tex>p</tex>-ого) бита.
}}
'''Например''', <tex>\operatorname{ulp}(-1100) = 1, \operatorname{ulp}(1) = 0.001</tex> в <tex>p</tex>-битной арифметике.
Так же полезной нотацией является <tex> \operatorname{err}(a \circledast b) </tex>, которая обозначает ошибку округлении результата выполнения операции <tex>\circledast</tex>. Отметим, что если <tex>\operatorname{ulp}</tex> всегда беззнаковая величина, то <tex>\operatorname{err}</tex> может иметь знак. Для базовых операций (сложение, вычитание, умножение) <tex> a \circledast b = a \ast b + \operatorname{err}(a \circledast b)</tex>, и точное округление гарантирует, что <tex> |\operatorname{err}(a \circledast b)| \leqslant \frac{1}{2}\operatorname{ulp}(a \circledast b)</tex>.
==Свойства==
{{Лемма
|id=
lemma1
|statement=
Пусть <tex>a \oplus b = a + b + err(a \oplus b) </tex>.
Как следствие, верна следующая лемма:
{{Лемма
|id=
lemma2
|statement=
Ошибка округления <tex>err(a \oplus b)</tex> может быть представлена в <tex>p</tex> битах.
}}
 
{{Лемма
|id=
lemma3
|statement=
Пусть <tex>|a + b| \leqslant |b|</tex> и <tex>|a + b| \leqslant |a|</tex>. Тогда <tex>a \oplus b = a + b</tex>. (Аналогично для вычитания).
{{Лемма
|id=
lemma4
|statement=
Пусть <tex>b \in \left [ \frac{a}{2}, 2a \right ]</tex>. Тогда <tex>a \ominus b = a - b</tex>.
==Алгоритмы суммирования==
===Простое суммирование===
Важной базовой операцией во всех алгоритмах, основанных на представлении чисел в виде расширенийразложений, является сумма двух <tex>p</tex>-битных величин, результатом которой является расширение разложение длины два. Есть два алгоритма для выполнения этой задачи - алгоритмы Деккера и Кнута.
====Сумма по Деккеру====
{{Теорема
Dekker
|statement=
Пусть <tex>a</tex> и <tex>b</tex> есть <tex>p</tex>-битные числа, причем <tex>|a| \geqslant |b|</tex>. Тогда следующий алгоритм вернет неперекрывающееся расширение разложение <tex>x + y</tex> такое, что <tex>a + b = x + y</tex>, где <tex>x</tex> - приближение (аппроксимация) суммы <tex>a + b</tex>, а <tex>y</tex> представляет собой ошибку округления при вычислении <tex>x</tex>. (Аналогично для вычитания).
}}
<br clear="all" />
<wikitex>Проблема с использованием этой процедуры заключается в требовании <tex>|a| \geqslant |b|</tex>. Если это заранее не известно, то необходимо выполнить сравнение перед ее вызовом. В большинстве С компиляторов, возможно, самым быстрым переносимым способом реализовать эту проверку является выражение <tex>if ((a > b) == (a > -b))</tex>. На эту проверку уйдет некоторое время, но увеличение времени может быть на удивление большим из-за современных процессоров с суперскалярными и конвейерными архитектурами, в которых вызов условного оператора может сбросить ветку предсказаний.Это объяснение лишь гипотетическое и зависит от машины, но алгоритм $\mathrm {TwoSum}$, что будет описан ниже, избегает этого сравнения посредством трех дополнительных операций, что обычно на практике даже быстрее. Конечно же, $\mathrm {FastTwoSum }$ все же быстрее, если результат сравнения известен ''априори''.</wikitex>
====Сумма по Кнуту====
Knuth
|statement=
Пусть <tex>a</tex> и <tex>b</tex> есть <tex>p</tex>-битные числа, причем <tex>p > 3</tex>. Тогда следующий алгоритм вернет неперекрывающееся расширение разложение <tex>x + y</tex> такое, что <tex>a + b = x + y</tex>, где <tex>x</tex> {{- --}} приближение (аппроксимация) суммы <tex>a + b</tex>, а <tex>y</tex> представляет собой ошибку округления при вычислении <tex>x</tex>.
}}
<tex>7</tex> <tex>return (x, y)</tex>
<wikitex>Пример работы последнего алгоритма, когда $|a| < |b|$ и $|a| < |x|$. Сумма $11.11 + 1101$ будет представлена в виде расширения разложение $10000 + 0.11$.
[[Файл:Adaptive_10.jpg|слева]]
<br clear="all" />
{{Лемма
|statement=
Пусть $x$ и $y$ {{- --}} величины, возращенные алгоритмом $\mathrm {FastTwoSum}$ или $\mathrm {TwoSum}$. На машине с округлением до ближайшего четного $x$ и $y$ {{--- }} несмежные.
}}
</wikitex>
===Суммирование расширенийразложений===<wikitex>Базируясь на алгоритмах сложения двух <tex>p</tex>-битных чисел, описанных выше, можно предложить алгоритмы суммирования расширенийразложений. В этой статье их будет предложено два: $\mathrm {ExpansionSum }$ и $\mathrm {FastExpansionSum}$. Первый алгоритм прибавляет к $m$-элементному расширению разложению $n$-элементное расширение разложение за время <tex>O(mn)</tex>, в то время, как последний алгоритм делает это за <tex>O(n + m)</tex>.
Несмотря на такое различие в асимптотике, первый алгоритм на практике может оказаться быстрее на расширениях разложениях, чей размер мал и фиксирован, потому что программные циклы могут быть полностью развернуты, а косвенные расходы времени исчезают, так как можно отказаться от использования массива). Линейный же алгоритм имеет определенные условия, при которых подобные оптимизации невозможны.
ExpansionSum имеет свойство, что если входные данные были неперекрывающимися (несмежными), то и на выходе мы получим расширения разложения с соотв.соответствующим свойством.
FastExpansionSum имеет несколько серьезных недостатков.Первый из них - алгоритм не сохраняет свойств неперекрываемости/несмежности. Второй - алгоритм основывается на округлении до ближайшего четного, что делает его непереносимым.<wikitex>Как правило, общим недостатком алгоритмов работы с расширениями разложениями является то, что в расширении разложении на выходе могут быть нулевые элементы, даже если в исходном расширении разложении их не было. Например, если подать на вход расширениея разложениея $1111 + 0.0101$ и $1100 + 0.11$, то результатом будет $11100 + 0 + 0 + 0.0001$. К счастью, алгоритмы, описанные в этой статье хорошо справляются с этой проблемой.
</wikitex>
Перед описанием алгоритмов суммирования расширений разложений, приведем алгоритм прибавления к расширению разложению <tex>p</tex>-битного числа.
====Grow expansionСумма разложения и числа====<wikitex>{{Теорема
|statement=
Пусть $e = \sum^{m}_{i=1}e_i$ - неперекрывающееся $m$-компонентное расширение разложение; $b$ - $p$-битное число, где $p \geqslant 3$. Предполагается, что $e_1, e_2, \dots, e_m$ отсортированы в '''возрастающем''' порядке, причем все компоненты ненулевые. Тогда следующий алгоритм вернет такое расширение разложение $h$, что $h = \sum^{m + 1}_{i=1}h_i = e + b$, где компоненты $h$ также отсортированы в возрастающем порядке и не равны нулю. Алгоритм также сохраняет свойство несмежности/неперекрываемости.
}}
$3$ $(Q_i, h_i) \Leftarrow TwoSum(Q_{i-1}, e_i)$
$4$ $h_{m+1} \Leftarrow Q_m$
$5$ $return $ $h$
$Q_i$ - приближенное значение суммы $b$ и первых $i$ компонент $e$.
{{Лемма
|statement=
Каждая их первых $m$ компонент расширения разложения $h$ не больше чем соответствующая компонента из $e$, то есть $|h_i| \leqslant |e_i|$. Более того, $|h_1| \leqslant |b|$.
}}
</wikitex>
 ====Expansion sumПростое суммирование====<wikitex>{{Теорема
|statement=
Пусть $e = \sum^{m}_{i=1}e_i$ и $f = \sum^{n}_{i=1}f_i$ - неперекрывающиеся $m$- и $n$-компонентные расширения разложения, компоненты которых - $p$-битные числа, где $p \geqslant 3$. Предполагается, что компоненты обоих расширений разложений отсортированы в '''возрастающем''' порядке, причем все компоненты ненулевые. Тогда следующий алгоритм вернет такое расширение разложение $h$, что $h = \sum^{m + n}_{i=1}h_i = e + f$, где компоненты $h$ также отсортированы в возрастающем порядке и не равны нулю. Алгоритм также сохраняет свойство несмежности/неперекрываемости.
}}
<br clear="all" />
====Fast expansion sumБыстрое суммирование====<wikitex>В отличие от $\mathrm {ExpansionSum}$, $\mathrm {FastExpansionSum}$ не сохраняет свойства неперекрываемости и несмежности, но гарантирует, что если если входное расширение разложение было ''строго неперекрывающимся'', то на выходе получится расширение разложение с таким же свойством.
{{Определение
|definition=
Расширение Разложение называется '''строго неперекрывающимся''', если если его компоненты попарно неперекрываются, ни одна компонента не смежна никаким двум другим, а также любая пара смежных компонент состоит из степеней двойки.
}}
'''Например''', $11000 + 11$ и $10000 + 1000 + 10 + 1$ {{--- }} строго неперекрывающиеся расширения разложения, а $11100 + 11$ и $100 + 10 + 1$ {{--- }} нет.
Для расширения разложения эта характеристика означает, что ноль в записи расширения разложения должен появляться ''как минимум'' каждые $p + 1$ бит. Например, расширениеразложение, каждая компонента которого есть $p$-битное число, и максимальная величина которой равна $1111$, может быть не больше $1111.01111011110\dots$.
Каждое несмежное расширение разложение является строго неперекрывающимся, а также каждое строго неперекрывающееся расширение разложение есть неперекрывающееся. В общем случае, обратное не верно.
Утверждается, что после предположения о том, что все расширения разложения строго неперекрывающиеся, алгоритм, приведенный ниже, корректно работает на машинах с округлением до ближайшего четного.
{{Теорема
|statement=
Пусть $e = \sum^{m}_{i=1}e_i$ и $f = \sum^{n}_{i=1}f_i$ {{- --}} строго неперекрывающиеся $m$- и $n$-компонентные расширенияразложения, компоненты которых {{--- }} $p$-битные числа, где $p \geqslant 4$. Предполагается, что компоненты обоих расширений разложений отсортированы в '''возрастающем''' порядке, причем все компоненты ненулевые. Тогда следующий алгоритм,запущенный на машине с правилом округления до ближайшего четного, вернет такое расширение разложение $h$, что $h = \sum^{m + n}_{i=1}h_i = e + f$, где компоненты $h$ также отсортированы в возрастающем порядке и не равны нулю. Алгоритм также сохраняет свойство несмежности/неперекрываемости.
}}
<br clear="all" />
</wikitex>
 
==Ссылки==
* [http://en.wikipedia.org/wiki/Rounding Rounding]
* [http://www.google.ru/url?sa=t&rct=j&q=Adaptive%2Bprecision%2Barithmetic&source=web&cd=4&ved=0CDgQFjAD&url=http%3A%2F%2Fwww.cs.berkeley.edu%2F~jrs%2Fpapers%2Frobustr.pdf&ei=fZigTp6xHoWa-waGxKSOBQ&usg=AFQjCNGl9T4V-Y02Wi99TjgSDLFotO5xeg&sig2=pQ5SPC_lRGCtwoBZhMhGhA Jonathan Richard Shewchuk, "Adaptive Precision Floating-Point Arithmetic and Fast Robust Geometric Predicates"]
 
 
[[Категория: Вычислительная геометрия]]
419
правок

Навигация