Comments 165
всё ради «быстродействия», «совместимости» или продаж нового оборудования!
А если честно то такое реализовать на FPGA, ASIC и новых чипах будет ещё сложнее чем флоаты (которые уже реализованы в аппаратных блоках в старших версиях). И самое худшее — приведёт к сильному удлинению конвеера CPU/GPU, ещё к более сложному управляющему устройству, к ещё более сложным алгоритмам предсказаний переходов инструкций и будущих данных и к ещё большим продажам новых IP ядер несовместимых с старыми. (может быть в этом и цель? как на авто: сделать «эко» режим и объявить остальное нормально ездящее вне закона? тем самым создав иллюзию повышения продаж)
IEEE (тут даже одинаковые выичлсения на одной и той же системе могут дать разные результаты)
Можно примерчик?
habr.com/ru/company/intel/blog/205970
И вот вам пример.
habr.com/ru/company/intel/blog/205970/#comment_7117974
Сорри за занудство, но где все-таки пример одинаковых вычислений?
Все ваши ссылки поясняют тривиальную вещь: разный порядок вычисления выражения дает разные результаты.
Сорри за занудство, но где все-таки пример одинаковых вычислений?Вот, например: "Math.Round returns different values for the same input.". ¯\_(ツ)_/¯
Why, you ask, can that happen? Good question; thanks for asking. Here’s the answer (with emphasis on the word “often”; the behavior depends on your hardware, compiler, etc.): floating point calculations and comparisons are often performed by special hardware that often contain special registers, and those registers often have more bits than a double. That means that intermediate floating point computations often have more bits than sizeof(double), and when a floating point value is written to RAM, it often gets truncated, often losing some bits of precision.
Это не про IEEE754. Это про то, как инженеры, делавшие конкретную платформу, хотели как лучше, а устроили всем геморрой. С 8087, например, такая история.
Остальные будем разбирать, или сами осилите?
Вот тут на sql.ru обсуждалось, что Oracle pl/sql даже при перестановке мест множителя дает различные результаты
sql.ru Oracle форум
А как воспроизвести ?
Я попробовал вот так в консоле браузера — ноль консоль логов произошло
for (let i = 0 ; i < 10e6; i++) {
let x = Math.random() * 10e8;
let y = Math.random() * 10e8;
let z = Math.random() * 10e8;
if (x*(y + z) !== (z + y) *x) {
console.log({x,y})
}
}
Ассоциативность:
(1e-200*1e200)*1e200 => 1e+200
1e-200*(1e200*1e200) => Infinity
Дистрибутивность:
1e200*(1e200-1e200) => 0
(1e200*1e200)-(1e200*1e200) => NaN
for (let i = 0 ; i < 1e6; i++) {
let x = Math.random() * 10e8;
let y = Math.random() * 10e8;
let z = Math.random() * 10e8;
if (x*(y + z) !== x*z + y*x) {
console.log({x,y})
}
}и ассоциативность (с
if (x*(y*z) !== (x*y)*z) {
Результат:
node check_math.js |wc -l
353838и
node check_math2.js |wc -l
350279(Примерно треть считает не совсем точно.)
Самое простое: перенос вычислений на NVDIA CUDA
Это означает, что величины с небольшой экспонентой, которые чаще всего используются, могут иметь большую точность, а реже используемые очень крупные и очень мелкие числа будут иметь меньшую точность.
И вот на этом возникает смутное ощущение, что равенство 32-битного Posit и 64-битного IEEE во всех-всех задачах есть только на бумаге...
Например, лёгкая манипуляция из того же отчёта LLNL:
Глядите, как всё хорошо — в трёх случаях из четырёх Posit даёт лучшую точность, чем IEEE single float. Но: давайте подумаем, а в какой реальной задаче будет нужно расстояние до звезды и какие масштабы будут с ним соседствовать? И получаем, что в реальной задаче рядом, вероятно, будет, либо астрономическая единица, либо порядок диаметра галактики. Поэтому все четыре примера, по сути, не отражают реальную сферу применения.
При этом: с привычным IEEE float у меня будет моя относительная погрешность 10-7 вне зависимости от того, решу я расстояние выражать в световых годах или километрах. А адаптивное число бит в экспоненте выглядит как много-много радости по поводу выбора масштабов в реальных вычислительных задачах.
На этом слайде 32-битные представления сравниваются.
На 64 битах Posit в обоих последних примерах, видимо, будет хуже, т.к. в IEEE экспонента всего 11 бит.
Ну вот они и показывают, что в single float будет уже бред, а posit хоть как-то, но держится. Другой вопрос — кому реально понадобится в программе расстояние между звёздами в планковских длинах выражать?
Лучше справляется – это не разрушает цепочку вычислений, деградируя до NaN или ±Infinity.
Возьмите любую игру: на зуме вылезают глюки Z-buffer'а из-за проблем с точностью вычислений на больших дистанциях (из-за того что вся геометрия вгоняется в один юнит в итоге после матричных перемножений. В смысле область размером 1x1x1 юнит)
То же самое у софтовых рендереров, когда, скажем, у горе-ахренетектора, (который не знает и не должен знать ничего о флоатах), здание находится в сцене в 2 млн сантиметров от нуля(по плану!) и рендерится «зубьями» (ошибка нахождения пересечения луча с треугольником, находящимся куда Макар гусей не гонял).
OpenEXR с его half-float спасибо тоже скажет.
PS: Кому-то выше отвечал, так получилось.
1 парсек = 3,2615637772… светового года
100,5*3,2615637772… = 331,6 (100,5*3.3) светового года
Хотя «точное» произведение будет равно 327,7871596…
А вот 100,500 парсек, где 500 — значащие цифры из измерений, будет уже 327,831 светового года.
Аналогично с очень большими экспонентами — при умножении мантиссы складываются, а точность дробной части равна точности наименьшей из перемножаемых.
То есть, мы по исходным данным точно знаем, что расстояние никак не равно 331.6 светового года — оно меньше.
Не знаю, что он хотел доказать, но показал только лишь то, что правила нужно не только зубрить, но и понимать.
Общее количество значащих цифр к ошибке, о которой он говорит, не имеет отношения. Если точность измерений одна десятая парсека, то она и будет одна десятая парсека, неважно 100 там парсеков или 1000000. И погрешность после умножения соответственно будет одна и та же.
Другое дело, что умножение флоатов в фиксированной разрядной сетке с этим тоже не связано. Там да, все значащие цифры важны. Пока влазим в разрядную сетку, относительные величины множителей вообще не имеют значения.
И погрешность после умножения соответственно будет одна и та же.Нет, не одна и та же. Точность, с которой известно второе слагаемое тоже влияет. Причём если мы урожаем, то влияет именно относительная точность (назад в школу, изучать как связаны умножение, сложение и логарифмы).
1.5 +- 0.1
δx = 0.1/1.5 ~= 0.0667
100.5 +- 0.1
δx = 0.1/100.5 ~= 0.001
1000000.5 +- 0.1
δx = 0.1/1000000.5 ~= 0.0000001
3.2615637772
δx ~= 0%
1.5 (0.0667) * 3.2615637772 (0) = 4.8923456658 (0.0667)
Δx = x * δx = 4.8923456658 * 0.0667 = 0.32631945590886 ~= 0.33
100.5 (0.001) * 3.2615637772 (0) = 327.7871596086 (0.001)
Δx = x * δx = 327.7871596086 * 0.001 = 0.3277871596086 ~= 0.33
1000000.5 (0.0000001) * 3.2615637772 (0) = 3261565.4079818886 (0.0000001)
Δx = x * δx = 3261565.4079818886 * 0.0000001 = 0.32615654079818886 ~= 0.33
В целом да, но вроде бы знаки отбрасываются после точного умножения, а не до.
> И получаем, что в реальной задаче рядом, вероятно, будет, либо астрономическая единица, либо порядок диаметра галактики.
Что-то мне подсказывает, что диаметрами галактики вряд ли будут считать, скорее уж метрами. Хотя, конечно, примеры «меньше единицы» надо было привести.
Я про то, что они привели сферический парсек в вакууме как пример.
posit, судя по описанию, работает через размен порядка на мантиссу. Поэтому числа вблизи 1 имеют больше значащих разрядов в мантиссе, а чем дальше, тем меньше.
Вопрос — если в реальной задаче за базовую единицу измерения взят парсек, то какие там будут встречаться характерные цифры? Насколько мне известно, это будет либо порядка 10-5пк (т.е. астрономической единицы), либо порядка 103-6пк (диаметр галактики).
Вывод: парсек как характерный масштаб вообще-то не очень полезен для применения с posit-числами. Только вот обычный float прожуёт и симуляцию планетной системы с расстояниями в парсеках, и симуляцию галактики, и не подавится. При этом программисту не нужно прилагать никаких усилий, чтобы специально обезразмерить задачу.
Так что представление это очень и очень нишевое, как мне кажется.
> Только вот обычный float прожуёт и симуляцию планетной системы с расстояниями в парсеках, и симуляцию галактики, и не подавится. При этом программисту не нужно прилагать никаких усилий, чтобы специально обезразмерить задачу.
Это до тех пор пока у вас в одном из вычислений не вылезет ±Infinity или NaN. Или накопление ошибки до неузнаваемости исказит результат. В этом смысле как раз posit с большей вероятностью прожует симуляцию галактики: у него другая стратегия underflow/overflow так что необычайно маленькое расстояние между двумя телами не приведет к делению на 0 и в конце вычислений на вас не будет смотреть NaN.
так что необычайно маленькое расстояние между двумя телами не приведет к делению на 0 и в конце вычислений на вас не будет смотреть NaNВот только ответ, который дадут поситы может на пару порядков отличаться от точного математического. И никакой возможности сказать — это более-менее точный ответ или практически случайное число, нет. Обычные флоаты тоже страдают от данной проблемы, но в гораздо меньшей степени.
Еще поситы обсуждали на gamedev.ru. Лично я для себя вынес, что поситы непригодны для сколь-нибудь серьезных применений.
Обычные флоаты честно выдерживают требуемую точность (исключение — underflow), либо выдают «не шмогла».Выдерживают требуемую точность они ровно в рамках одной операции, а дальше накопление ошибки может завести вас буквально куда угодно.
Поситы же дают какой-то ответ всегда.Не всегда. Деление на 0 – это ошибка (±inf), операция над ±inf – это тоже ±inf.
деградируя при этом точность до нуляЭто голословное утверждение. Нужно рассматривать конкретные примеры. На практике выходит так что там где float из-за одной из тысячи последовательных операций (например, a1*b1 + a2*b2 + ...) свалился в NaN/inf/0, posit даст приемлемый с точки зрения точности результат.
Вы же скорее описали ситуацию когда тупо неверно выбрана размерность для вычислений. Тут вам никакая магия не поможет.
Вот только ответ, который дадут поситы может на пару порядков отличаться от точного математического.Может, как и в случае с float.
И никакой возможности сказать — это более-менее точный ответ или практически случайное число, нет.Если вам нужно не только произвести вычисление, но и гарантировать его точность – был придуман тип Valids (одна из вариаций Posit), которая на любую операцию возвращает не просто значение, а диапазон, в котором реальное значение гарантированно лежит.
Правда, эта идея не получила широкой поддержки, поскольку оказалось что там где требуется гарантировать точность, практичнее использовать Quire (другой формат операций над Posit), которые просто не теряют точность во время операций суммирования/умножения.
Так что с одной стороны не рекомендую делать выводы только по gamedev.ru, с другой – позиты действительно на сегодняшний день бесполезны, если только вы не планируете производить свое железо. Аппаратной поддержки нигде нет, а софтварная эмуляция не сравнима с апаратной реализацией float. Так что пока можно применять разве что для эфективного сжатия чисел с потерей точности.
Выдерживают требуемую точность они ровно в рамках одной операции, а дальше накопление ошибки может завести вас буквально куда угодно.Накопление погрешности — это отдельная проблема. И как раз анализ на численную неустойчивость в случае поситов на порядок сложнее из-за их плавающей точности.
Вы же скорее описали ситуацию когда тупо неверно выбрана размерность для вычислений. Тут вам никакая магия не поможет.Вот как раз про это я и пишу. Обычные флоаты нечувствительны к небольшой неидеальности выбора погрешности, часто можно хоть на 10 порядков вверх или вниз сдвигать без какой-либо деградации результата. У поситов же даже один десятичный порядок сдвига от единицы может привести к потере точности.
Аппаратной поддержки нигде нетИ неудивительно. Я на gamedev.ru как раз разбирал этот вопрос.
Всё верно, но...
Поситы — это чтобы меньше потерять при сохранении некого float-результата в меньшем кол-ве бит.
В сравнении с 754 поситы теряют они действительно меньше, но иначе.
При неправильном использовании, ни 754, ни поситы не спасут.
Однако, при правильном использовании поситы позволяют потерять меньше.
А статья дурная, начиная с заголовка.
Интерферометр тоже за уши притянут — при их объемах и бюджете можно ASIC-и использовать для delta-кодирования и т.п.
Обычные флоаты нечувствительны к небольшой неидеальности выбора погрешности, часто можно хоть на 10 порядков вверх или вниз сдвигать без какой-либо деградации результата. У поситов же даже один десятичный порядок сдвига от единицы может привести к потере точности.Не понимаю откуда вы это взяли. У позитов есть стандартные варианты по аналогии с fp16, fp32, fp64 – p8e0, p16e1, p32e2. Если вы тупо замените fp32 на p32e2 вряд ли у вас будет деградация результата из-за того что вы что-то там не откалибровали. Можем проверить.
Не понимаю откуда вы это взяли.Да прямо из исходной статьи про поситы. Там есть красивый график точности для разных вариантов представления чисел. У стандартных флоатов ровный забор, а у поситов треугольник — в районе 1 они, действительно, точнее, но шаг влево или вправо от «вершины» и они начинают проигрывать.
Как я понял, DrSmile всё-таки прав.
8 бит на показатель в es2 требуются начиная с 11111000, и это (5-1)*22 + 0 = 16, а 216 — это 65536. Начиная с 220 точность ниже, чем у float. А ещё на один бит она станет хуже после 224, это далеко не 4 миллиарда.
Как я уже говорил, это плата за большую точность внутри диапазона. В центре диапазона у вас будет 4 дополнительных бита точности, так что пока большинство ваших значений лежит в интервале 1/2^20..2^20, а это диапазон в 12 десятичных порядков, точность ваших вычислений будет выше всегда.
На практике можно работать в гораздо более широком диапазоне, потому что точность деградирует достаточно медленно, и эта деградация компенсируется большей точностью на остальных числах. Здесь уже не имея конкретной задачи невозможно говорить какой подход проявит себя лучше.
Меня волнует то, в первую очередь, что распространение ошибки при фиксированной длине мантиссы достаточно просто проанализировать.
У меня как-то Inf-ы и NaN-ы вылезали, в основном, из-за catastrophic cancellation при вычитании (Inf от x/0 и NaN от 0/0). А от него не видно, как posit систематическим образом спасает.
Отсутствие инвариантности по масштабу напрягает. С Float32/64 я могу знать, что (2.0 * x) / (3.0 * x) отличается от 2.0 / 3.0 разве что последним битом (при разумных значениях x, естественно). У posit с этим как-то непонятно.
И главный вопрос, который у меня возник — я верно понял, что удвоение порядка — это плюс два бита в порядок и, соответственно, минус два из мантиссы? Если да, то более широкая мантисса по сравнению с IEEE там всё-таки в довольно узком диапазоне значений будет, и как-то ей хвастаться некорректно...
У меня как-то Inf-ы и NaN-ы вылезали, в основном, из-за catastrophic cancellation при вычитании (Inf от x/0 и NaN от 0/0). А от него не видно, как posit систематическим образом спасает.На самом деле как раз спасает систематически отличной от float стратегией underflow/overflow. Вы никогда не получите ±Infinity складывая, вычитая или умножая числа, отличные от ±Infinity. Вы никогда не получите ноль, поделив любое число отличное от нуля на любое другое число. На первый взгляд это звучит дико, но нужно понимать что понятия чисел в системе posit отличаются от чисел в системе float.
С Float32/64 я могу знать, что (2.0 * x) / (3.0 * x) отличается от 2.0 / 3.0 разве что последним битом (при разумных значениях x, естественно). У posit с этим как-то непонятно.Точность будет выше, поскольку числа 2 и 3 имеют порядок близкий к единице. Причем, как количество идентичных результатов будет выше, так и погрешность в случае ошибки будет меньше.
И главный вопрос, который у меня возник — я верно понял, что удвоение порядка — это плюс два бита в порядок и, соответственно, минус два из мантиссы?Скорее всего не верно. В рамках одной и той же системы в posit мантисса может занимать разное количество бит, и даже 0 бит. Это определятся значением regime – количеством первых последовательных бит, которые равны друг другу. Это немного сложная концпеция, которую сходу не объяснить, но основная суть в том что для разных значений мантисса будет иметь разную размерность.
Cуммирование:
Деление:
В вашем случае ситуация еще лучше, потому что 2 и 3 — удобные числа для posit.
И прежде чем писать про 1%-1.5% процента разберитесь что такое decimal error. Ошибка в 1% это decimal error = 0.004. И разберитесь что здесь находится по оси X.
Ошибка в 1% это decimal error = 0.004Сорри, пропустил нолик. Правильное утверждение: «кого волнует, что ошибка 15%, а не целых 10%?».
Если вы хотите увидеть там какой-то сюрприз, его там не будет. Точность float8 это 3 бинарных порядка мантисы (1/16), а значит в наихудшем случае decimal error = log10 (1+1/16)/1 = 0.02632. Ровно в это значение флоаты и упираются на графике.
Это графики для восьмибитных чисел.Ну тогда это еще одна подтасовка, ибо 8-битных IEEE-float, в принципе, не существуют. И нигде нету сравнений с реальными 32/64-битными флоатами. Почему они взяли 3-битную мантиссу и 4-битную экспоненту, а не наоборот? Наверное, потому что иначе графики не такие красивые получались.
Почему они взяли 3-битную мантиссу и 4-битную экспонентуЧтобы порядок экспонент был хоть как-то сравним.
Держите сравнение f32 и p32e2. Построено методом Монте-Карло на одном миллионе сэмплов.
Тот же график в другом масштабе, для полноты картины:
Я здесь вижу обратную сторону медали: у позитов наростание ошибок очень плавное и напоминает ужатую гиперболу. Это означает что аномально большие ошибки будут встречаться гораздо реже. У флоатов же 20% суммирований заканчивается с рандомной точностью.
А вы не знаете, случаем, почему Густафсон предлагает настолько жмотиться и делать стандартными p32es2?
Я тут прикинул, es3 лучше почти по всем параметрам, кроме потери одного бита в диапазоне от 1/16 до 16. Но блин, это уже нужно очень сильно заморочиться, чтобы в алгоритме важные переменные за этот диапазон гарантированно не вышли.
А дальше там более-менее прилично — точность не ниже float32 на диапазоне 18 десятичных порядков, и потом не так радикально падает.
"Не так радикально", правда, это число незначащих цифр O(|log2x|), как и в денормализованных float, только константа в es раз меньше (т.к. порядок в унарной форме записывается по сути).
Для es=3 у нас шаг не 16 бит, а 256 (2^2^es).
Умножение:
Деление:
Кажется тут позиты уделывают флоаты однозначно.
А что вы тут конкретно сравниваете?
Если вы берёте два позита, перемножаете их и сравниваете с "точным" результатом перемножения исходных позитов, — это не имеет смысла, т.к. даже в "точном" произведении значащих бит не больше, чем в исходных сомножителях.
На самом деле сведение к f64 дает неотлимый график:
(p1 * p2) as f64 / (p1 as f64 * p2 as f64)
Так это бессмысленно.
Вы берите сначала дроби, перемножайте как рациональные числа и сравнивайте результат с произведениями этих же дробей, округлённых в float и в posit.
Эти графики в первую очередь демонстрируют как позиты или флоаты сохраняют точность гуляя по полю собстенных значений при выполнении операций. Тут важно не столько то как выполняется операция, как то какие значения нам доступны в поле чисел.
Вас же скорее интересует точность преобразования рациональных чисел во флоаты/позиты. Такой график тоже можно построить, но это про другое.
Ну посмотрите, что ли, разницу (ab)c — a(bc). Это же входит в категорию "гуляют по полю собственных значений"?
Ожидаемо, график флоатов посередине улетает в бесконечность, потому что для случайных флоатов a*b*c в половине случаев приводит к underflow/overflow.
Вообще не потеряли точность флоаты в 35% случаев, позиты – в 67%.
В защиту флоатов я тут скажу, что без переполнения они таки гарантируют точность 2-23.
А дальше… В общем, Густафсону указывают на проблему: перемножение очень большого числа на очень маленькое вызывает потерю значащих цифр в произведении. Он этот аргумент просто игнорирует. А в физике такой юзкейс сплошь и рядом — тераваттный лазер с импульсом в наносекунду, тяга Saturn V поделить на массу Saturn V, произведение постоянной Планка на терагерцовую частоту, гравитационная постоянная на массу Земли, число Авогадро на постоянную Больцмана и т.д. Наиболее патологический для позитов случай будет (x * y) / y, где y — порядка 1025-30.
В лазерах сейчас, например, начали играться с такой штукой — как два импульса по мишени взаимодействуют между собой. И вполне могут бить сначала микросекунду мегаваттом, потом пикосекунду тераваттом. Произвольно менять единицы времени и мощности — ата-та, у нас же ещё есть плотность мишени, размер пятна, их тогда тоже придётся перемасштабировать. А если сложная симуляция, в которой подключается взаимодействие ЭМ поля с веществом — нужна ещё, как минимум, скорость света.
Когда мы два числа умножаем, то число значащих бит в мантиссе вырастает вдвое (в случае переменной битности — будет суммой числа бит обоих чисел). Вот, например, ежели два 32-битных целых числа умножить — будет 64-битное число. Итого, чтобы при умножении позитов потерялось точности больше, нежели у флоатов, надобно, чтобы значащих бит у них было в два раза меньше, нежели у флоатов.
Другими словами, даже позит 2,4е51 из картинки с примерами при умножении на 2е-52:
— будет иметь точность выше, чем у 32-битного флоата
— будет иметь точность всего на один бит меньше, чем максимальная точность соответствующих позитов
— при том, что во флоатах ни одно исходное число вовсе невозможно выразить!
Быть может, потому Густафсон такие аргументы и игнорирует, ибо для него это наверняка самоочевидно.
На картинке, оказывается, не совсем позиты, а другое представление плавающей точки тоже с динамической мантиссой, но дельта-кодированим показателя. У него точность выше, чем у IEEE, в ещё более узком интервале. Зато динамический диапазон шире, конечно.
У позитов чуть-чуть шире диапазон, где они сопоставимы с IEEE, но гораздо уже динамический диапазон и быстрее деградация точности.
IEEE — нормальная система. Jack of all trades, king of none. Что-то там с выжиманием лишних пары бит получается соорудить только для половинной точности без серьёных недостатков.
Да, ещё:
Вот, например, ежели два 32-битных целых числа умножить — будет 64-битное число.
Это так, но если мы работаем с плавающей точкой — за пределами размера мантиссы исходных сомножителей биты в произведении не имеют смысла.
Пример: чему равно 1,1x101x1,2x101 с точностью до 3 значащих цифр? 1,32x102? А вот фигвам, 1,25x102. Почему? Потому что 1,1x101 — это у меня было округлённое 10,6, а 1,2x101 — округлённое 11,8.
Пример, естественно, просто чтобы показать, что вычисление произведения с точностью выше точности сомножителей только маскирует то, что мы его точно не знаем, но не делает ответ более точным.
Представление IEEE нацелено на то, чтобы ситуация "биты в мантиссе есть, а толку от них уже нет" появлялась только при сложении чисел разного знака. Во всех представлениях с переменной длиной мантиссы это возникает как при сложении с разным знаком, так и при перемножении с порядками разных знаков.
Что-то я не пойму, к чему этот пример? Ежели у нас есть 11 и 12, то мы никак не можем догадаться, что на самом деле имеются в виду 10,6 и 11,8, хоть с позитом, хоть с флоатом. Для этого уже ванга нужна.
> так и при перемножении с порядками разных знаков
Насколько я всегда себе представлял, умножение производится так — мантиссы умножаем, а экспоненты складываем. Итого у мантиссы получается в два раза больше бит, а у экспоненты — либо столько же, либо на один больше. Итого точность нигде не теряется, кроме как в той операции, что мы должны поместить результат в число той же битности, что и операнд. Конечно, ежели мы возьмём позиты, у которых в мантиссе остались считанные биты, умножим и получим приблизительно единицу, то число значащих бит в ней будет меньше, нежели хотелось бы, но все эти биты мы потеряли ещё тогда, когда получали исходные числа, а умножение здесь ни при чём. И бит в ней будет столько же, сколько в двух исходных числах в сумме, то бишь в некотором смысле в точности мы приобретём.
Ежели у нас есть 11 и 12, то мы никак не можем догадаться, что на самом деле имеются в виду 10,6 и 11,8, хоть с позитом, хоть с флоатом. Для этого уже ванга нужна.
Что это сразу Ванга? Из-за масштабной неинвариантности мантиссы в позитах 1,06 при умножении на 10 может потерять точность и превратиться в 11 (утрирую, но не сильно — с es2 уже на масштабе 16 потеряется один бит). А 1,18, скажем, при делении на 10 превратится в 0,12. При обратном перемножении 11 на 0,12 масштаб будет снова около 1, и мы "возвращаем" третью значащую цифру, т.е. результат будет 1,32. Но он не равен округлённому до 3 знаков произведению исходных 1,06 и 1,18! А округление 1,3 до двух цифр — будет честным. Т.е позиты при умножении и делении "сжирают" информацию. Потом биты в мантиссу могут вернуться, а толку-то, если там уже мусор?
Эти графики в первую очередь демонстрируют как позиты или флоаты сохраняют точность гуляя по полю собстенных значений при выполнении операций.Я тут подумал и понял, что это неверно. Эти графики показывают вносимую погрешность операций, никак не соотнося ее с погрешностью исходных операндов. Т. е. если исходная погрешность 10-7, то дополнительной погрешностью операции в районе 10-9 можно просто пренебречь.
На практике, при анализе накопления ошибок, вносимой погрешностью операции можно пренебречь почти всегда, по сравнению с погрешностью операндов, кроме случая денормалов для обычных флоатов и умножений/делений поситов. Наверное, надо строить графики относительной вносимой погрешности: отношения погрешности, вносимой операцией, к погрешности из-за дискретизации исходных операндов.
Эти графики показывают вносимую погрешность операций, никак не соотнося ее с погрешностью исходных операндов.Согласен, я это и имел ввиду – как позиты сохраняют точность во время выполнения операций.
Я понял вашу мысль. Погрешность a + b не может быть ниже чем сумма погрешностей a и b. Но вот это мне не понятно:
надо строить графики относительной вносимой погрешности: отношения погрешности, вносимой операцией, к погрешности из-за дискретизации исходных операндов.Что нам даст отношение этих погрешностей? Обе погрешности могут быть очень высоким и очень низким, а их отношение – одинаковым, хотя второй вариант очевидно предпочтительней. Пока мне приходит в голову брать максимальную погрешность из двух для каждого случая. Вроде бы это правильных подход, потому что каждая погрешность вносит свой вклад независимо относительно истинного значения, а не последовательно.
Что нам даст отношение этих погрешностей?Ухудшение погрешности в результате операции. Погрешности операндов присутствуют до операции и про них есть отдельный график, интересно же узнать, насколько неточность операции увеличивает уже существующую погрешность.
В общем, такая относительная погрешность — это характеристика баланса точности до и после операции. Если такая погрешность много больше единицы, то это значит, что точность результата недостаточна, происходит потеря информации. Если много меньше — наоборот, точность результата излишняя, биты тратятся впустую.
Линия позитов идет ниже линии флоатов для любого X.
таки на графике деления (второй в вашем сообщении), линия позитов выглядит выше в диапазоне 12~14 000
То что линия позитов идет выше линии флоатов в том месте при делении говорит о том, что при делении флаты дают чуть больше точных результатов чем позиты. Но если результат неточен (а это большинство случаев), то позиты дают меньшую ошибку.
Возможно, наибольшими возможностями posit будут обладать в области машинного обучения, где 16-битные числа можно использовать для обучения, а 8-битные – для проверки.
Мне кажется, в области машинного обучения лучше вообще использовать fixed-point значения. Там числа находятся недалеко от нуля и имеют значащие разряды с десятых долей, зачем там астрономические или микроскопические величины, где экспонента экономит биты.
А теперь немного критики:
https://marc-b-reynolds.github.io/math/2019/02/06/Posit1.html
https://hal.inria.fr/hal-01959581v4
Из последней ссылки:
When posits are better than floats of the same size, they provide one or two extra digits of accuracy. When they are worse than floats, the degradation of accuracy can be arbitrarily large. This simple observation should prevent us from rushing to replace all the floats with posits.
Так же следует учитывать что posit справляется со значительным количеством случаев, где float не выдает адекватного результата вообще. Это стоит учитывать при сравнении. Но всегда можно найти такой пример, в котором любая разумная система исчисления будет лучше другой. В этом смысле действительно не корректно говорить что posit лучше чем float.
Главная проблема: точность posit хуже на сложении маленьких чисел.
Но ведь в примере «Addition closure» Posit лучше: нет NaN и Overflow ( и Exact больше)?
Думаю при сложении нет NaN, Overflow
Таким образом с одной стороны экономится память и затратность обращения к ней (в большинстве случаев можно использовать posit меньшей размерности чем float), с другой – увеличивается точность за счет использования аккамуляторов не приводящих к потере точности.
Таким образом с одной стороны экономится память и затратность обращения к ней (в большинстве случаев можно использовать posit меньшей размерности чем float), с другой – увеличивается точность за счет использования аккамуляторов не приводящих к потере точности.Это уже пробовали. 8087 называется. Сколько глюков это порождает — и не сосчитать…
В результате, например, программа, рассчитывающая сворачиваемость белков — теряли или дублировала варианты при запуске параллельно на нескольких машинах.
Появление x86-64 и переход на SSE2, где все числа всегда одного размера — это было просто «щастя» какое-то!
То есть выигрыш от более точных вычислений никого не волновал, но вот неустойчивость результата (на одинаковых входных данных получаем разные выходные) — напрягала очень сильно.
Это дезинформация. У х87 есть специальная настройка (поле PC в регистре управления) для понижения точности регистровых вычислений хоть до double хоть до single. Сделано как раз для того чтобы было проще делать детерминированные алгоритмы вычислений.
То что ряд компиляторов не умеет этим пользоваться — проблема компиляторов а не проца.
Появление x86-64 и переход на SSE2, где все числа всегда одного размера — это было просто «щастя» какое-то!
Соответственно, это тут тоже ни при чём.
То что ряд компиляторов не умеет этим пользоваться — проблема компиляторов а не проца.Если подавляющее большинство компиляторов чем-то не умеют пользоваться — это становится проблемами проца.
Беда в том, что достаточно быстро менять эти настройки можно только на оригинальном 8087. А на каким-нибудь Pentium'е это весь пайплайн сбивает. Потому и не дёргает их никто перед каждым вычислением.
Да и не было это так задумано. Идея была, что программисты будут всё рачками писать — отсюда и стековая архитектура и прочее… А они — не захотели, редиски. И сейчас не хотят.
Соотвественно, если рассматривать реальные программы, а не чьи-то влажные мечты, то SSE решил эту проблему.Появление x86-64 и переход на SSE2, где все числа всегда одного размера — это было просто «щастя» какое-то!Соответственно, это тут тоже ни при чём.
А на каким-нибудь Pentium'е это весь пайплайн сбивает. Потому и не дёргает их никто перед каждым вычислением.
Перед каждым и не надо, надо в начале работы программы и иметь соглашение что библиотеки эту настройку не ломают — то есть если меняют (а обычно им незачем) то на выходе ставят обратно, а конвеер там всё равно сбивается из-за передачи управления.
На всякий случай также упомяну, что с точки зрения си будет корректным проведение всех вычислений с точностью double. То есть разница между
float a,b,c;
a = a*b; a = a*c;
и
float a,b,c;
a = a*b*c;
(то, что во втором случае будет точнее и возможно другой результат) заложена в стандарте языка и х87 тут опять ни при чём.
Да и не было это так задумано.
Эта настройка была сделана именно для поддержки языков в которых подразумевается не максимальная точность вычислений.
Перед каждым и не надо, надо в начале работы программы и иметь соглашение что библиотеки эту настройку не ломают — то есть если меняют (а обычно им незачем) то на выходе ставят обратно, а конвеер там всё равно сбивается из-за передачи управления.Как вы это себе представляете, интересно? Во многих программах есть и 32-битные float и 64-битные double. Вы ж не можете в 8087 включить сразу два режима.
На всякий случай также упомяну, что с точки зрения си будет корректным проведение всех вычислений с точностью double.Только в режиме без оптимизаций. Если же вы перемножаете два числа и кладёте их в третье — в этот момент дожно произойти преобразование во float.
(то, что во втором случае будет точнее и возможно другой результат) заложена в стандарте языка и х87 тут опять ни при чём.Стандарт это разрешает, но не требует. И большинство процессоров всё считают в одинарной точности.
Эта настройка была сделана именно для поддержки языков в которых подразумевается не максимальная точность вычислений.Эта настройка была сделана, когда 8086 считался временной заглушкой до прихода iAPX 423. На 8087 (оригинальном) смена режимов была дёшева, а о XXI веке никто не думал…
8087 — это весьма продуманный дизайн для тех возможностей, которые были у его создателей. Но когда от него удалось отказаться в пользу SSE — разработчики «вздохнули с облегчением»…
Ну у вас какая-то странная логика. Только давайте если отвечать то отвечать на сообщение целиком а не на вырванные из контекста фразы.
Если кому-то хочется иметь один формат дробных чисел, который не придётся никуда конвертировать — ставим PC=double и используем переменные типа double. И никаких проблем с оптимизациями итд.
Если хотите использовать два типа — используйте два, но урезание точности до float при записи в float-переменную придётся указывать явно, да (например FSTP+FLD если хотим в итоге оставить в регистре).
SSE2 или что там ещё с "одним типом" никак не решает проблему конвертации, оно просто запрещает использовать второй формат. Но вас никто его и раньше использовать не заставлял. Так что если не хотите лишних конвертаций — не используйте два формата, при том "не использовать" можно добровольно даже при наличии технической возможности, а не только по причине отсутствия поддержки второго.
Только давайте если отвечать то отвечать на сообщение целиком а не на вырванные из контекста фразы.Давайте. Это к вам, в первую очередь, как бы, относится.
С чего всё началось? С моего утверждения:
Проблемой было то, что люди не программируют программы на ассемблере, а подавляющее большинство языков программирования «не умеют» в два представления. Соответственно в зависимости от того что, куда и когда компилятор решит посчитать на регистрах (в extended точности), а когда что-то попадёт в память результаты получались разные. Причём сильно разные.По моему всё достаточно понятно написано: никого не волнуют теоретические возможности, которые вы можете получить, программируя в машинных кодах. Людей интересует Fortran и С++, в первую очередь. И с реальными компиляторами, а не с чем-то возможным, но несуществующим…
Если хотите использовать два типа — используйте два, но урезание точности до float при записи в float-переменную придётся указывать явно, да (например FSTP+FLD если хотим в итоге оставить в регистре).Ну и кто это будет делать? Физик или химик, для которого слово «регистр» если что-то и значит, то либо что-то из музыки, либо из бухгалтерии?
SSE2 или что там ещё с «одним типом» никак не решает проблему конвертации, оно просто запрещает использовать второй формат.Как это? Инструкции cvtpd2ps и cvtps2pd никто не отменял. Оно просто делает использование двух форматов предсказуемым без выхода на уровень ассемблера и регистров.
С чего, собственно, вся дискуссия и началась. Теоретически 8087 — предоставляет кучу возможностей. Практически — пользоваться ими сложно, в некоторых случаях и невозможно.
Реализация технических возможностей — дело компилятора. В частности корректное приведение типа к float, будь оно задано в коде явно или неявно — его дело. Поэтому программистам на языках высокого уровня знать все эти детали как раз не нужно.
И лень или что там ещё у авторов компиляторов, которая помешала всё это правильно реализовать никак нельзя ставить в минус самому процу.
Подозреваю, что будь х87 система команд анонсирована сейчас, современные компиляторописатели к ней бы сделали бы поддержку всего что нужно. Но она была анонсирована в 80х годах, а тогда никто над этими вопросами не заморачивался, фронтэнд к калькулятору сделали, считает правильно и ладно (и вобщем-то для большинства применений и сейчас этого достаточно).
Инструкции cvtpd2ps и cvtps2pd никто не отменял.
У меня из предыдущих комментов сложилось впечатление что там всего один тип поддерживается вообще. Видимо неверное.
Ну, значит сделали более удобное приведение типов чем было в х87. С тем же успехом могли добавить инструкцию "обрезать точность ST(0) до single/double" для реализации приведений без трогающих оперативную память FSTP+FLD и тогда бы оно в плане обсуждаемой темы не уступало SSE.
И лень или что там ещё у авторов компиляторов, которая помешала всё это правильно реализовать никак нельзя ставить в минус самому процу.Вот как только вы «правильный» компилятор, который делает всё «как надо» — можете рассказывать сказки.
А пока реальные компиляторы этого делать не умеют — это недостаток процессора.
Подозреваю, что будь х87 система команд анонсирована сейчас, современные компиляторописатели к ней бы сделали бы поддержку всего что нужно.Нет. Не сделали бы. FLDCW сбрасывает конвеер, в результате какой-нибудь FMUL начинает занимает не «в среднем один такт», а «в среднем семь-восемь тактов».
Никакому компилятору подобное бы не простили — как не простили компилятору Ада пользователи iAPX432.
Да, можно пытаться как-то сдледить за control flow и прочее… но написать
-mfpu=sse проще и надёжнее.У меня из предыдущих комментов сложилось впечатление что там всего один тип поддерживается вообще. Видимо неверное.Нет, конечно. Один 32-битный тип — это первая версия SSE. Её только некоторые игры используют. SSE2 — это «второй набор» инструкций для работы с 64-битными числами. Вот его уже можно вместо 8087 использовать.
С тем же успехом могли добавить инструкцию «обрезать точность ST(0) до single/double» для реализации приведений без трогающих оперативную память FSTP+FLD и тогда бы оно в плане обсуждаемой темы не уступало SSE.Всё равно уступало бы. Стековая архитектура катастрофически плохо совместима с существующими компиляторами.
Но это уже эффект второго порядка, да.
Это имеет критическое значение, в частности, при реализации FFT произвольного размера (бо́льших 32K) алгоритмом Блюстейна.Ну и как вы собираетесь его делать на каких-нибудь ARM'ах?
В том-то и дело, что 80-битная точность — это удел маргиналов. Статья вообще говорит о том, что и 64 бита — не так уж нужны.
Но если реально сильно приспичило — есть же __float128. Да, он «ручками» реализован, но за счёт того, что x87 сейчас реализуют «по остаточному принципу»… может и быстрее 80-битного типа оказаться…
Ну и как вы собираетесь его делать на каких-нибудь ARM'ах?А причём тут ARM-ы? Платформ всяких разных множество существует, в том числе и без поддержки плавающей точки вообще в принципе.
В том-то и дело, что 80-битная точность — это удел маргиналов. Статья вообще говорит о том, что и 64 бита — не так уж нужны.Ну это только до тех пор, пока вы лично не столкнётесь с ситуацией, когда точности double катастрофически не хватает. Ещё один пример навскидку — аппроксимация методом наименьших квадратов.
Но если реально сильно приспичило — есть же __float128. Да, он «ручками» реализован, но за счёт того, что x87 сейчас реализуют «по остаточному принципу»… может и быстрее 80-битного типа оказаться…А я слышал (источник не приведу, в какой-то из специализированной литературе по низкоуровневой оптимизации), что и FPU, и SSE используют один и тот же модуль вычислений. Как минимум эксперимент на доступных PC показал, что FPU и SSE не распараллеливаются, а в единичных операциях SSE ничуть не быстрее.
Ну это только до тех пор, пока вы лично не столкнётесь с ситуацией, когда точности double катастрофически не хватает.Если я лично столкнусь, то мне придётся думать как без этого выкрутится.
А причём тут ARM-ы? Платформ всяких разных множество существует, в том числе и без поддержки плавающей точки вообще в принципе.Проблема в том, что подавляющее большинство платформ, в том числе практически все платформы на которых делают как «игрушечные» игровые рассчёты, так и большие, за миллионы долларов — поддерживают 32-битные флоаты и 64-битные. А никаких 80-битных — не поддерживают.
Вряд ли такая ситуация могла сохраняться долгие годы, если бы 80-битные были активно востребованы.
А я слышал (источник не приведу, в какой-то из специализированной литературе по низкоуровневой оптимизации), что и FPU, и SSE используют один и тот же модуль вычислений.Кто там такое сказал? О каком процессоре идёт речь? И вы вообще точно осознаёте о чём говорите?
Они не могут использовать один модуль как минимум потому, что для 8087 вам нужно поддерживать 80-битную точность, а для SSE этого не нужно.
Типичная ситуация такая: есть умножитель, который может выполняться вычисления с полной 80-битной точностью. Один. И есть 2-4 умножителя для 64-битных чисел и 4-8 для 32-битных. Можете у Agnerа посмотреть подробнее. При этом да, иногда умножитель на 80 бит используется как один из умножителей на 64-бит (это может быть проще, чем делать их совсем независимыми), но, как правило, умножителей на 64-бита больше.
То есть получить выигрыш от использования смеси 8087 и SSE кода не получится — а вот гемор получить можно изрядный. Проще всего про 8087 просто забыть и всё.
а в единичных операциях SSE ничуть не быстрее.А вот У Агнера — совсем другие результаты. Скажем PMUL на Ryzen'е имеет задержку 5 тактов и может исполнять 5 инструкций в параллель (то есть можно «снимать» один результат раз в такт), а MULSD — задержка 4 такта и 8 инструкций (то есть можно «снимать» два результата за такт). Ускорение — более, чем в два раза.
Понятно, что это всё микробенчмарки и в реальных вычислениях всё может быть не совсем так… но факт остаётся фактом: все известные мне люди, которые умеют что-то оптимизировать в один голос говорят: перед тем, как заниматься оптимизациями — убедитесь что у вас в программе нет 8087 инструкций. То есть вот это — абсолютно железное правило.
Никого не волнует что у вас там происходит в программе, если она не использует SSE. Корректность — да, важна, конечно, Legacy, всё такое. Скорость… ну что получится без большого влияния на SSE впихнуть — то и будет.
То есть получить выигрыш от использования смеси 8087 и SSE кода не получится — а вот гемор получить можно изрядный. Проще всего про 8087 просто забыть и всё.Искренне не понимаю, откуда такое неприятие — точно также и про SSE можно забыть, потому что AVX есть. Я в частности успешно применял смесь FPU и SSE для оптимизации double-double арифметики — получив двойное ускорение на сложение и десятикратное на умножение (по сравнению с си++ кодом отсюда).
все известные мне люди, которые умеют что-то оптимизировать в один голос говорят: перед тем, как заниматься оптимизациями — убедитесь что у вас в программе нет 8087 инструкций. То есть вот это — абсолютно железное правило.Я тоже люблю и (возможно) умею оптимизировать — одна их моих ассемблерных оптимизаций получила 40-кратное по сравнению с си, но вашей точки зрения категорически не разделяю. Вроде бы наоборот — максимальное оптимизация и достигается использованием всех доступных команд и ресурсов.
Искренне не понимаю, откуда такое неприятие — точно также и про SSE можно забыть, потому что AVX есть.Более того, не можно, а нужно.
Смесь AVX и SSE кода очень сильно бьёт по производительности. Тут можно почитать.
Я в частности успешно применял смесь FPU и SSE для оптимизации double-double арифметики — получив двойное ускорение на сложение и десятикратное на умножение (по сравнению с си++ кодом отсюда).На каком процессоре? И как, собственно, вы применяли FPU? Пересылки между FPU и SSE очень дорого стоят, плохо верится в десятикратное ускорение.
Вроде бы наоборот — максимальное оптимизация и достигается использованием всех доступных команд и ресурсов.Вот только 8087 и SSE — не есть дополнительные ресурсы. Как вы сами же написали чуть выше: Как минимум эксперимент на доступных PC показал, что FPU и SSE не распараллеливаются.
Собственно основная проблема со всеми этими хитрыми трюками — в том, что вам постоянно нужно следить за появлением новых CPU и «подкручивать» код под них.
А поскольку о ваших интересах ни разработчики процессоров, ни кто-либо ещё особо не думают…
И как, собственно, вы применяли FPU?
add_dd_dd PROC uses esi edi dd1:PTR DWORD, dd2:PTR DWORD, dd3:PTR DWORD;[+]
local s1s:XMMWORD, s2s:XMMWORD
mov esi, dd1
mov edi, dd2
mov eax, dd3
movupd xmm0, XMMWORD PTR [esi];=a
movupd xmm1, XMMWORD PTR [edi];=b
movapd xmm3, xmm0;=sum
addpd xmm3, xmm1
movapd xmm4, xmm3;=bv
subpd xmm4, xmm0
movapd xmm5, xmm3;=av
subpd xmm5, xmm4
subpd xmm1,xmm4;=br
subpd xmm0,xmm5;=ar
addpd xmm0, xmm1;=err
movupd XMMWORD PTR s1s, xmm3
movupd XMMWORD PTR s2s, xmm0
fld REAL8 PTR s1s
fadd REAL8 PTR s1s+8
fadd REAL8 PTR s2s
fadd REAL8 PTR s2s+8
fstp REAL8 PTR [eax]
fld REAL8 PTR [eax]
fsubr REAL8 PTR s1s
fadd REAL8 PTR s1s+8
fadd REAL8 PTR s2s
fadd REAL8 PTR s2s+8
fstp REAL8 PTR [eax+8]
ret
add_dd_dd ENDPmul_dd_dd PROC uses esi edi dd1:PTR DWORD, dd2:PTR DWORD, dd3:PTR DWORD
local t1:XMMWORD, t2:XMMWORD
mov esi, dd1
mov edi, dd2
mov eax, dd3
movupd xmm7, XMMWORD PTR _split64;=split
;
; product - x
fld REAL8 PTR [esi]
movsd xmm0, REAL8 PTR [esi];=a
fld REAL8 PTR [esi+8]; |a|a`|
movhpd xmm0, REAL8 PTR [edi];загрузка в старший разряд ;=b
fmul REAL8 PTR [edi]; |a|a`*b|
movsd xmm3, xmm0
fld REAL8 PTR [edi+8]; |a|a`*b|b`|
mulsd xmm3, REAL8 PTR [edi];=x=a*b
fmulp st(2),st(0); |a*b`|a`*b|
movsd REAL8 PTR [eax], xmm3
faddp st(1),st(0); |a*b`+a`*b|
; split
movapd xmm1,xmm0;=a, b
mulpd xmm1, xmm7;=ca, cb
movapd xmm2, xmm1;=ca, cb
subpd xmm2, xmm0;=abig, bbig
subpd xmm1, xmm2;=ahi, bhi
subpd xmm0, xmm1;=alo, blo
; перестановка чисел между двумя регистрами
movapd xmm2, xmm1;= ahi, bhi
unpcklpd xmm2, xmm0;= ahi, alo
unpckhpd xmm1, xmm0;= bhi, blo
; парное умножение
movapd xmm0, xmm2;= ahi, alo
mulpd xmm0, xmm1;= ahi*bhi, alo*blo
shufpd xmm1, xmm1, 1; перестановка чисел ;= blo, bhi
mulpd xmm1, xmm2;= ahi*blo, alo*bhi
; финальное вычитание
subsd xmm3, xmm0;err1
subsd xmm3, xmm1;err2
shufpd xmm0, xmm0, 1
shufpd xmm1, xmm1, 1
subsd xmm3, xmm1;err3
subsd xmm0, xmm3;y
movsd REAL8 PTR [eax+8], xmm0
; сложение с FPU
fadd REAL8 PTR [eax+8]
fstp REAL8 PTR [eax+8]
ret
mul_dd_dd ENDPЧестно говоря не верится. Всё что я тут вижу — «закат солнца вручную», когда вместо HADDPD/HSUBPD/etc используется, зачем-то, 8087.
Которые можно было бы точно так же вычислить не гоняя данные через память и не задействуя 8087.
когда вместо HADDPD/HSUBPD/etc используется, зачем-то, 8087Очевидно, зачем — чтобы не терять точность.
Которые можно было бы точно так же вычислить не гоняя данные через память и не задействуя 8087.Можно было — только по времени подольше выйдет. Я же ссылку на оригинальный код зачем давал? К чему столько вычислений, подробно расписано например вот в этой работе — и в первую очередь стоит обратить внимание на разницу в реализации Two-Sum и Quick-Two-Sum, в котором первое слагаемое должно быть больше второго по модулю.
Впрочем, вполне возможно, что и я, и авторы оригинального кода что-то упустили — поэтому буду безмерно благодарен за конкретную более быструю реализацию.
Честно говоря не верится
double z = 1e+200;
double zzz = sqrt(z*z + z*z);
cout << "default: " << zzz << endl;
__asm
{
fld z
fmul st,st
fadd st,st
fsqrt
fstp zzz
}
cout << "FPU: " << zzz << endl;Результат:
default: inf
FPU: 1.41421e+200Осталось понять только почему на HPC об этом никто не плачет и спокойно на GPU всё считают, которые в дополнительную точность не умеют.
Осталось понять только почему на HPC об этом никто не плачет и спокойно на GPU всё считаютОчевидно потому, что с такими задачами и алгоритмами не сталкивались — я же выше два конкретных примера привёл, и изучать тему тоже не от нечего делать начал. Требования к точности FFT принципиально отличаются в зависимости от того, используется оно для визуализации спектра или же для кепстральных преобразований. Для визуализации поворота вектора в 2D и single может быть избыточен.
Ваша "оптимизация" кода заключается лишь в одном: вы использовали аппаратную работу с 80-битными числами вместо программной работы со 128-битными.
Если же у вас есть необходимость работы только с 32- и 64-битными числами, то про FPU действительно стоит забыть.
Ваша «оптимизация» кода заключается лишь в одном: вы использовали аппаратную работу с 80-битными числами вместо программной работы со 128-битными.Было бы неплохо кавычки в слове «оптимизация» обосновать пруфом на бенчмарк, в котором программная реализация 128-бит флоат будет работать быстрее. Таки концепция double-double не моя, моя лишь ассемблерная реализация с FPU.
Если же у вас есть необходимость работы только с 32- и 64-битными числами, то про FPU действительно стоит забыть.Пример выше с переполнением до inf по-прежнему не аргумент? И входные, и выходные данные укладываются в 64 бит. Тут же половина комментариев по профиту posit-а упирает на то, что он в inf не вываливается.
Тут же половина комментариев по профиту posit-а упирает на то, что он в inf не вываливается.Ну дык собственно тот факт, что народ не ломанулся срочно переходить везде и всюду на posit'ы как раз и говорит о том, что само по себе это преимущество — важно, пожалуй, только для FP16 или гипотетического FP8.
Преимущество-то может и есть, но без аппаратного ускорения практический интерес к ним невелик.
народ не ломанулся срочно переходить везде и всюду на posit'ыБолее вероятная причина этого в том, что IEEE 754 — это уже устоявшийся, хорошо изученный аппаратно реализованный формат. А альтернативные реализации float-ов были ещё с самого начала его возникновения — и этот был выбран как наиболее удачный из предложенных на тот момент.
Начнём с того, что ваша ассемблерная реализация не является эквивалентом double-double, потому что вы работаете с 80-битными числами, а не 128-битными. Тут нечего сравнивать.
Но вы пошли на это осознанно, потому что точности 80-битных чисел для ваших задач вполне достаточно.
Пример выше с переполнением до inf по-прежнему не аргумент?
Нет, не аргумент. Вы упомянули задачу вычисления FFT. Но FFT — устойчивое преобразование, в нём не возникает больших чисел. Поэтому преимущество от posit-ов вы можете получить за счёт того, что у них больше бит отводится мантиссе для рабочих диапазонов значений чисел.
Нет, не аргумент. Вы упомянули задачу вычисления FFT. Но FFT — устойчивое преобразование, в нём не возникает больших чисел.А ещё я упоминал кепстр, которой считается как FFT от логарифма FFT — поэтому при обратном преобразовании небольшие погрешности как раз и могут приводить к очень большим числам.
Поэтому преимущество от posit-ов вы можете получить за счёт того, что у них больше бит отводится мантиссе для рабочих диапазонов значений чисел.Вполне может быть — возможное преимущество от posit-ов я не отрицаю.
Логарифм переводит относительные погрешности в абсолютные. То есть абсолютная разница между 100 и 101 будет такая же, что и между 1 и 1.01. Если при дискретизации логарифма вас устраивает подобное, то можно смело использовать fixed point и для логарифма.
Если же вы хотите остаться в рамках исходных погрешностей, то придётся посчитать погрешности.
Пусть исходный сигнал длиной N лежит в диапазоне [0, 1] с погрешностью eps. Коэффициенты ПФ будут лежать в диапазоне [0, N], при этом для восстановления сигнала по ним их погрешность должна быть не больше (eps/N). Нормируем коэффициенты ПФ в диапазон [0, 1], получим погрешность E = (eps/N^2).
Значения логарифма будут в диапазоне от (ln E) до нуля, а требуемая погрешность равна ln(1 + E) < E. Если логарифм нормировать в [-1, 0], то получим погрешность -E / (ln E).
Далее если второй раз требуется посчитать ПФ, то получим итоговую абсолютную погрешность:
E / (N^2 ln E) = eps / (N^4 (ln eps — 2 ln N)).
Если eps = 2^(-24), N = 32768, то имеем что-то близкое к 2^(-90). Вывод: 128-битного fixed point хватит с запасом, а вот 80-bit float будет маловато.
Я не совсем понял, почему у вас погрешность получилась обратно пропорционально N — в то время как на практике наоборот, чем больше N, тем больше погрешность.
Это требования к точности промежуточных вычислений, чтобы при использовании обратимых преобразований получить нулевую разницу между результатом и исходными данными. Чем больше N, тем выше должна быть точность промежуточных вычислений.
ваша ассемблерная реализация не является эквивалентом double-double, потому что вы работаете с 80-битными числами, а не 128-битными. Тут нечего сравнивать.double-double — это невычисленная сумма двух 64-битных чисел, а не одно 128-битное. Какой формат используется для промежуточных вычислений — не имеет значения, пока они укладываются в заданные погрешности.
Понятно, что это всё микробенчмарки и в реальных вычислениях всё может быть не совсем такВ реальных вычислениях всё совершенно точно будет не так — потому что данные для вычислений из памяти надо сначала считать, а результаты — записать обратно. И один-единственный промах кэша ценой в 1000 тактов аннулирует все микро-оптимизации на корню. Достаточно взять чуть более сложный, чем линейное суммирование элементов массива алгоритм, и разницы в скорости между FPU-only и SSE/AVX уже не будет.
Достаточно взять чуть более сложный, чем линейное суммирование элементов массива алгоритм, и разницы в скорости между FPU-only и SSE/AVX уже не будет.Если бы это было правдой, то ни в каких SSE1/2/3/4/etc не было бы смысла.
Думаю что истина, как обычно, посередине: двукратный выигрыш SSE не даст, но процентов 50 — возможно.
За исключением случаев, где, действительно, 80-битные числа могут что-то дать… но если бы они реально часто стреляли, то процессоров, где есть не только 32-битные и 64-битные числа было бы явно больше…
Появление x86-64 и переход на SSE2, где все числа всегда одного размера — это было просто «щастя» какое-то!
До появления операции FMA, в которой результаты промежуточных вычислений не огрубляются. Захотел компилятор оптимизировать код и использовать новую инструкцию — получите другой результат.
точность posit хуже на сложении маленьких чисел.Преведенные вами графики построены по двум разным системам исчестления, так что одинаковые координаты не соответствуют одним и тем же числам.
Если под маленькими числами подразумевать числа со значением экспоненты близкой к нулю (скажем от 1e-10 до 1e10) то позиты как раз таки предлагают гораздо более высокую точность в этом диапазоне.
Заголовок статьи совсем не совпадает со смыслом.
Никакого "избавления" нет, только альтернативное представление для хранения, которое лучше в одних случаях и хуже в других.
Если уполовинить количество битов, можно не только уменьшить объёмы кэша, памяти и хранилища для этих чисел, но и серьёзно уменьшить ширину канала, необходимого для передачи их на процессор и обратно. Это главная причина, по которой арифметика на базе posit, по его мнению, даст от двойного до четверного ускорения расчётов по сравнению с числами с плавающей запятой от IEEE.
Однако эта арифметика сама по себе более медленная. Причём неизвестно насколько, так как нет её реализации в железе. Скажем, если они будет в 5 раз медленнее, то толку от такого перехода? И, честно говоря, я сомневаюсь, что на самом деле, скажем, Float64 медленнее Float32 настолько. В своё время тестировал свои численные алгоритмы на разную точность и падение производительности было процентов 20, точно не больше 50. А хранить данные можно в меньшей точности, с сжатием, с преобразованиями в другой формат и так далее.
И, честно говоря, я сомневаюсь, что на самом деле, скажем, Float64 медленнее Float32 настолько.Двухкратная разница как минимум (Tesla V100), а может быть и разница раз в так в 30 (Quadro RTX 8000: 16312 GFLOPS vs 509.8 GFLOPS). То что вы получали разницу в 20% обозначает просто что плавучка не была «бутылочным горлышком» ваших алгоритмов. Может они в память упирались или просто там было много уровней индирекции.
Насчёт же огромной разницы у ускорителей на 32- и 64-битной точности. Конкретно на RTX 8000 стоит ровно в 32 раза меньше F64 блоков чем F32. У V100, как можете догадаться, первых блоков всего в два раза меньше. Поэтому разница не удивительна, так как фактически вычисления той и другой точности производятся на разных устройствах. Так что, давайте таки сравнивать одно и тоже, например CPU. Или Xeon Phi.
Или Xeon Phi.Xeon Phi используется на горсточке суперкомпьютеров. А вот GPU — на большинстве. При этом эти же рассчёты, зачастую, люди делают на своих домашних GPU, если нет доступа к суперкомпьютерам.
Интересно, а в posit возможны оптимизации вычислений при помощи магических констант типа «Магическая константа» 0x5f3759df? Интересно было бы почитать об этом.
Новый подход может помочь нам избавиться от вычислений с плавающей запятой