Комментарии 52
Вычислять квадратный корень может потребоваться и на процессоре, не имеющем не только его аппаратной реализации, но и вообще какой-либо поддержки операций с плавающей запятой -- скажем, на любом микроконтроллере с ядром Cortex-M3 или вообще на Ардуине. И вот там подобная магия может дать многократный прирост (особенно если писать на ассемблере). В общем, если нечто становится не особо актуальным для ПК, может быть вполне себе актуальным для других платформ.
А я на 8-битном микроконтроллере, при небольших значениях, иногда вычислял квадратный корень из целых чисел, вычитая из данного нечетные числа, начиная с единицы, пока разность не станет равной нулю. Их число и дает значение корня.
Например: корень из 9 равен 1- 3 -5 то есть три. И так далее.
Ну, естественно, когда квадрат не целый, то с округлением.
![a^b = \sum_{n=1}^{a^{\left [ \frac{b}{2} \right ]}}{a^{(\frac{b}{2})}-a^{\left [ \frac{b}{2} \right ]}+2n-1}](https://habrastorage.org/getpro/habr/upload_files/f3c/111/c36/f3c111c361c7768992eed550c6f47cdc.svg)
Лет 15 назад вывел такую формулу для разложения натуральной степени натурального числа в арифметическую прогрессию. Буду признателен, если кто-то ее узнает и подскажет как она называется.
Для четной степени получается вот что:
Тут как раз количество слагаемых равно корню из исходного числа:
Поиграть с формулой можно тут https://scastie.scala-lang.org/JHhFSE6XRB2r5R5eBgEjUA
Не узнал. Но если заменить 
на 
и поменять лево и право местами, то получим
Можно эту формулу ещё упростить, и даже вычислительную сложность понизить:
Логично в принципе.
Ну вы преобразовали формулу "в лоб" в сумму арифметической прогрессии. Ещё чуть-чуть упростить и до умножения дойдем))
С нечетными степенями интереснее, там слагаемые не с 1 начинаются.
Вычислять квадратный корень может потребоваться и на процессоре, не имеющем не только его аппаратной реализации, но и вообще какой-либо поддержки операций с плавающей запятой
Что-то мне подсказывает, что тут алгоритм влоб (с поглощением 2 битов исходного числа за раз и генерацией одного) может всё ещё быть быстрее. Конечно, первое приближение по таблице окажется мгновенным, зато потом умножать побитово или побайтово можно запариться. Исключение -- если процессор умеет умножать целые числа, особенно если 32-битные с 64-битным результатом (М3 как раз).
В бенчмарке измеряемые функции идут после тяжелой операции (rand). Из-за этого функция, которая идет первой (time(exact, EXACT)) выглядит чуть хуже по сравнению с другими функциями, чем на самом деле
Извините, но с чего бы ей выглядеть хуже? И когда LCG(https://en.wikipedia.org/wiki/Linear_congruential_generator) прямо таки стали "тяжелой операцией"(учитывая что в большинстве имплиментаций там и модуль по степени двойке)?
Спасибо, хорошая статья!
Буквально недавно на хабре спорил с одним человеком, который не верил, что вещественная арифметика может давать побитово разные результаты на процессорах Intel и AMD. А тут прямо конкретный пример программы по ссылкам.
Мы-то на практическом опыте на это натыкались ещё 20 лет назад. Но дампы в тетрадочку не записывали.
В общем случаи, не могут. Rsqrtss это совершенно отдельная инструкция, не прописанная ни в одном стандарте на плавучку. Которая дает аппроксимацию, а какой она будет это совершенно отдельный вопрос.
В стандартах вообще не прописываются конкретные инструкции процессора.
А суть вопроса в том, что одна и та же программа, будучи исполнена разными процессорами, может давать разные результаты в вещественной арифметике. Что бы ни было написано в стандартах, это просто эмпирический факт.
Нет это совершенно не эмпирический факт, когда последний раз плавучка отклонялась от стандарта, меняли целый выпуск процессоров. https://en.wikipedia.org/wiki/Pentium_FDIV_bug
И плавучка прописана в стандартах как минимум https://en.wikipedia.org/wiki/IEEE_754
и здесь обсуждение вопроса с совместимостью https://stackoverflow.com/questions/2234468/do-any-real-world-cpus-not-use-ieee-754
По вашей же собственной ссылке на Википедию написано, что хотя fdiv bug в большинстве случаев проявлялся малыми ошибками (подразумевая, что это, типа, дело житейское), но иногда ошибки могли достигать значительной величины. Поэтому и отзывали.
Вам приведён конкретный код, можете сами его скомпилировать и запустить.
Давайте вы мне покажете инструкции которые заявлены IEEE-754 complaint и которые дадут разные значения экспоненты и мантисы, и тогда мы продолжим разговор, что кто-то читал или нет. В данном случаи rsqrtss не является IEEE-754 complaint инструкцией, там даже в описании написано aproximate.Именно поэтому даже с агрессивной оптимизацией (-O3) компилятор ее не использовал, и пустился во все тяжкие лишь когда ему сказали "можешь ломать что хочешь" (-Ofast)
По-моему, тут каждый комментатор доказывает что-то своё. Лично мой тезис состоит в том, что процессоры Intel и AMD иногда дают побитово разные результаты на одинаковых программах с вещественными вычислениями, и это несомненный, эмпирически наблюдаемый факт. А как это соотносится с IEEE-754, меня вообще не очень волнует.
Так на чем основан данный эмпирический факт?
Инструкции rsqrtss, для которой черным по белому написано approximate, и соответственно сколь угодное отклонение от стандарта? И которую не один компилятор даже со стандартной агрессивной оптимизацией(-O3) не использует, именно из-за проблем с детерминированностью?
К примеру, я знаю кучу софта для научных вычислений, которые выдают одинаковые результаты на amd и intel, в которых очевидно любое отклонение от детерминированности привело бы к существенным изменениям результата, те же расчеты турбулентности.
PS: а так да, согласен, какие-то инструкции, могут давать разный результат, только если нужна воспроизводимость нормальные люди их не используют.
Софт для научных вычислений очень тщательно относится к вопросам оценки погрешности. И лишних разрядов обычно в результаты не выдаёт. Поэтому крайне маловероятно, чтобы битовое несоответствие результатов операций повлияло бы на расчёты - там по определению случайные погрешности данных и модели должны быть гораздо больше, чем обсуждаемая половина цены младшего разряда, иначе получится, что неправильно выбрана разрядность.
А стандарт IEEE-754 сам по себе – это вообще не верховная мера человеческого разума, не библия и не альфа и омега. На мейнфреймах, например, родная арифметика в принципе не имеет с IEEE-754 ничего общего.
Теория динамического хаоса началась с этого, что погрешность в младших разрядах приводила к совершенно разным прогнозам в простейших численных моделях погоды, которыми занимался Лоренц. Экспоненциальное разбегание сколь угодно близких траекторий, это широко известный постулат, для систем с динамическим хаосом, и турбулентности в частности. Так что да, в этих моделях даже одного флипа бита достаточно, чтобы получить абсолютно разные траектории результата.
Такую модель физически бессмысленно обсчитывать. Это значит, что истинные траектории, скорее всего, не могут быть представлены с помощью представимых в компьютере чисел, а находятся где-то между ними.
Дело не в точности, а в том что с какой точностью не считай, сколь угодно близкие результаты все равно будут разбегаться, экспоненциально. В широких массах это известно как "эффект бабочки". Расскажите про бессмысленность расчета тем, кто численно считает погоду.
https://en.wikipedia.org/wiki/Numerical_weather_prediction
"A more fundamental problem lies in the chaotic nature of the partial differential equations that govern the atmosphere. It is impossible to solve these equations exactly, and small errors grow with time (doubling about every five days). Present understanding is that this chaotic behavior limits accurate forecasts to about 14 days even with accurate input data and a flawless model."
“Small errors” здесь не означает цену младшего разряда.
А “эффект бабочки“ – это вообще из другой оперы. Его на самом деле, судя по всему, нет.
Да смоделируйте простейший аттрактор Лоренца, https://en.wikipedia.org/wiki/Lorenz_system и изменить младший разряд, и увидите как траектории разлетятся. Вы тут не со мной спорите, а с целой отраслью науки исследующей "динамический хаос" и его проявления с 60-х. Простейший доступный пример это турбулентность и погода.
"small errors" означает, что сколь угодно маленькая ошибка, при экспоненциальном разбегании, проявит себя довольно быстро.
Я 28 семестров изучал высшую математику так-то.
К тому моменту, как ошибка в младшем разряде проявит себя качественным изменением результата, всё вычисление теряет смысл.
А я практик, математик ваших не изучал, но если мы итеративную симуляцию, каждый новый шаг которой зависит от предыдущих расчетов, то сколь угодно малая погрешность приводит к абсолютно разным результам в конце. Вот наглядная демонстрация этого эффекта из-за различий в вычислениях с плавающей точкой между платформами. https://forum.unity.com/threads/soft-floating-points-calculations-and-determinism.1269275/ А вот решение этой проблемы с использованием софтверной эмуляции вычислений с плавающей точкой https://github.com/Kimbatt/unity-deterministic-physics И да, нам простым смертным, постоянно приходится работать с системами которые постоянно "разбегаются", ведь вообще не все физические движки гарантируют хоть какую-то детерминированность и расчеты физики у всех участников физики постоянно не сходятся друг с другом, но сетевые игры всё равно существуют, не смотря на "потерявшие весь смысл вычисления".
То, о чём вы пишете – это несколько другая проблема, чем обсуждавшаяся выше. У вас мы решаем вопрос воспроизводимости вычислительного процесса. Действительно, разные платформы аппаратно дают разные результаты, а софтверно их можно унифицировать.
Но это не значит, что один результат чем-то правильнее другого, и пунктуальное соблюдение IEEE-754 нам каким-то образом поможет с точки зрения научного физического моделирования. В реальности невозможно предсказать, как упадут мячики из вашего видео, поэтому оба разных результата одинаково незначимы с физической точки зрения.
Просто реальные мячики падают в полях вероятности :) а модель показывает одну из возможных секвенций :)
Но вообще не принято учитывать младший разряд, а для точности просто ещё раз считают и делают какой нибудь составной тип, для ученых весь вопрос выдает ли разный результат матлаб или их либа для чисел, если нет (что скорее всего), то и проблемы нет. А так то в теории можно и резонанс вызвать с вообще рандомными результатами и это от архитектуры и экземпляра железки зависит (тем кто говорит что нельзя - объясните, почему у меня видеокарта шебуршит когда в ллама считает строку, а когда в гпт2, то нет)
Смотрите работы по "Shadowing"
Частично это обходится расчетом статистических характеристик для ансамбля траекторий. Отдельную траекторию при этом не рассматривают.
При этом для диссипативных динамических систем рассматриваются характеристики "установившегося движения" на аттракторе, (не строго) ожидая что они в меньшей степени будут зависеть от указанной чувствительности чем частный "переходной процесс" установления режима.
"Истинные траектории" при этом могут не иметь физического смысла в отличии от "математически точного решения", которое можно апроксимировать с требуемой точностью (при необходимости) посредством "особо аккуратного" численного решения (и тут стандартной машинной точности уже не хватит).
Я писал игру, пошаговую стратегию. Полный обсчёт физики на клиентах позволяет обходиться только пересылкой нажатий клавиш в сетевой игре, что даёт мизерный трафик, а также другие плюшки. Физика простая, написана была на фрипаскале без всяких хитрых инструкций. В какой-то момент пришлось перейти с флоатов на программную реализацию чисел с фиксированной точкой. Просто потому что некоторые события в игре зависели от сравнений двух чисел, и на амд иногда получался другой результат сравнения, чем на интеловских процессорах, и симуляция расходилась. Бинарник был один и тот же, а дело было 17 лет назад.
Разработчики Factorio тоже сталкивались с подобным
Тут скорее вопрос в настройках компилятора. Нормализация (точнее ее отсутствие) влияет на точность, так как может сокращать число значащих разрядов. А ещё на разных архитектурах используются различные реализации математической библиотеки даже при одинаковых бинарниках, так как выбор реализации происходит при запуске бинарника. И если хотите одинакового результата везде, то нужно аккуратно настроить компилятор, статически ликовать некоторые библиотеки и лочить используемые расширения. Иначе никак. Особенно этим страдает интеловский компилятор, который на амд вообще подключает программную эмуляцию некоторых функций, которые при некоторых настройках оптимизации могут давать немного отличающиеся результаты. Плюс некоторые векторные функции имеют на разных процессорах разную разрядность. И если на какой-то платформе результат окажется более точным, то вы тоже получите расхождение.
Короче, там ничего сложного, но если вы действительно хотите бинарно идентичные результаты в математике, то нужно немного позаморачиваться.
В последний раз это произошло, когда https://randomascii.wordpress.com/2014/10/09/intel-underestimates-error-bounds-by-1-3-quintillion/
Вполне могут. Например, rsqrt инструкция точнее на AMD, что приводит к тому, что acosh, asin, asinh и atan2 фунции math.h точнее на AMD (не в glibc, в Intel Math Library, IML).
да, меня тоже про "22 значения из 2^32" улыбнуло)..
Было интересно почитать про математику лежащую в основе первого приближения. Но код замеров времени выполнения какой-то настораживающий. Где забота об изменении частоты, использование счетчиков производительности процессора, усреднение результатов, выбрасование аномалий из подсчета и прочее. Ведь недаром есть библиотеки для микробенчамрков типа этой https://github.com/google/benchmark . И "undefined behaviour" в коде немного напрягает. Понятно что исторически код был написан так, но при замере производительности вызов mempcy для преобразования "float" в "int" и обратно привел бы к такому же ассемблерному коду, но "undefined behaviour" бы не было.
а что за нуб на скрине с включенным cg_drawgun ??
p.s. простите не удержался
Так постоянно используем на DSP, MC, FPGA. В основном для вычисления модуля I/Q сигнала. В каких-то уже есть, так называемый, 10% обратный корень. Потом один или два приближения в зависимости от требований.
https://habr.com/ru/companies/infopulse/articles/336110/
А вот тут более подробно расписывается, как эта магия работает, и как её использовать для своих целей.
А я то думал что исходник игры запустят и заменят Хак на стандартный метод и протестят скорость
"Довольно странно, что первая аппроксимация стала вдвое быстрее. " - опечатка? Она же стала в 2 раза медленнее.
Полезен ли сегодня быстрый обратный квадратный корень из Quake III?