Математический реверс-инжиниринг

[spamsink]

Или "разговор программиста с математиком о магии числ".
Послесловие к моему предыдущему посту, чтобы не вышло допытывания инженера у математика о непонятном, как xaxam упоминал.

Имеется такая довольно быстрая функция, вычисляющая квадратный корень с хорошей точностью:
Пусть x = a*2^2b, где 0.25 <= a < 1. Тогда sqrt(x) = sqrt(a)*2^b.
Для вычисления sqrt(a) положим
r = 0.365681 + 5a/8
и 3 раза выполним
r = (r + a/r)/2
Максимальная относительная погрешность результата объявляется равной 1.28₁₀-13.

Вопрос к математическим обратным инженерам - каким образом получены значения 0.365681 и 5/8?
(Чтобы гуглические обратные инженеры не мучились ищучи, откуда я это взял: https://pure.tue.nl/ws/files/3023549/674735.pdf)

В книге Computer Approximations by J. F. Hart et al., 1968 (цит. по pers. comm.) приводится такая формула:
Пусть x = a*2^b, где 0.5 <= a < 1.
Тогда sqrt(x) = sqrt(a)*2^[b/2]*{b mod 2 = 0: 1; b mod 2 = 1: sqrt(2)}.
Для вычисления sqrt(a) положим
r = -0.204405847*a²+0.890194099*a+0.313567057
и 2 раза выполним
r = (r + a/r)/2
Ее максимальная относительная погрешность утверждается равной 1.5₁₀-12.

Как получены приведенные значения коэффициентов?

(Символ ⏨ заменен на ₁₀ по просьбе трудящихся пользователей Win7.)

Comments

[lesamohval.livejournal.com]

Я ССЗБ(

[spamsink]

Win7?

[lesamohval.livejournal.com]

ога. но установлен универсальный набор шрифтов!любые иероглифы вижу

[spamsink]

Видимо, ни один шрифт из этого набора не стоит.

[lesamohval.livejournal.com]

что это?

[spamsink]

Это список фонтов, в которых есть использованный символ.
Я сейчас попробовал. У меня Symbola стоит, тут в соотв. строке он виден, а браузеры не показывают.
Похоже, что проблема в самом Win7. Ладно, исправлю пост.

[kcmamu.livejournal.com]

5/8 выбрано, чтобы умножать было легко (а оптимальное значение 1/3 + sqrt(2)/4).

[spamsink]

У меня в пятницу столько же получилось, но оно настолько отличается от 5/8, что я не поверил.

Вызывает удивление, что в обоих алгоритмах используется деление вместо аппроксимации обратного значения и умножения.

[janatem.livejournal.com]

Наверно константы можно найти чисто в лоб, записав соответствующий функционал невязки и (приблизительно) найдя его минимум. Количество итераций задается ad hoc (потом проверяется, что полученная погрешность удовл). Вот еще та знаменитая константа из Квейка (https://ru.wikipedia.org/wiki/%D0%91%D1%8B%D1%81%D1%82%D1%80%D1%8B%D0%B9_%D0%BE%D0%B1%D1%80%D0%B0%D1%82%D0%BD%D1%8B%D0%B9_%D0%BA%D0%B2%D0%B0%D0%B4%D1%80%D0%B0%D1%82%D0%BD%D1%8B%D0%B9_%D0%BA%D0%BE%D1%80%D0%B5%D0%BD%D1%8C), там объясняется, как она была найдена, но я что-то не вкурил.

[spamsink]

Понятно, как она была найдена, но непонятно, что за магическое число 0.512964. :)

[spamsink]

Количество итераций в случае кв. корня определяется желаемой точностью и точностью полученного приближения: каждая итерация увеличивает точность вдвое. Или можно идти назад: задав, сколько итераций мы готовы делать в конце, легко получить необходимую точность первоначального приближения.

[pappadeux.livejournal.com]

гггг

вспомнил о том же