Алгоритм произвольной точности целочисленного квадратного корня?

9

Существуют ли какие-либо известные субквадратичные алгоритмы для вычисления минимального значения квадратного корня из nцелого бита?

Наивный алгоритм будет что-то вроде

def sqrt(x):
    r = 0
    i = x.bit_length() // 2
    while i >= 0:
        inc = (r << (i+1)) + (1 << (i*2))
        if inc <= x:
            x -= inc
            r += 1 << i
        i -= 1
    return r

Это требует O(n)итераций, каждая из которых включает в себя дополнения, которые являются O(n)временем, так что это O(n^2)общее время. Есть ли что-нибудь быстрее? Я знаю, что для случая умножения есть специальные алгоритмы, которые работают лучше, чем квадратичное время, но я не могу найти ничего для квадратных корней.

сурьма
источник
Мой ответ на что-то связанное может помочь cs.stackexchange.com/a/37338/12052 . Единственная проблема в том, что часть необходимого уравнения вам нужно будет найти эмпирически, чтобы настроить его точность.
Франческо Грамано
@FrancescoGramano: Извините, я не думаю, что это помогает.
Арьябхата
Кстати, это субквадратичное требование является частью большей проблемы? Потому что различие между простым квадратичным и сложным субквадратичным на практике может быть не таким уж большим. Или это просто теоретический интерес?
Арьябхата
@Aryabhata Извините, я не видел ваш комментарий ранее. Нет, это не часть большой проблемы, просто любопытство.
Сурьма

Ответы:

5

Вы можете использовать метод Ньютона или любой из ряда других методов для нахождения приближений к корням полинома .p(x)=x2c

Скорость сходимости для метода Ньютона будет квадратичной, что означает, что число правильных битов удваивается в каждой итерации. Это означает, что итераций метода Ньютона достаточно.O(lgn)

Каждая итерация метода Ньютона вычисляет

xj+1=xj(xj2c)/(2xj)=0.5xj+c2xj.

Битовая сложность умножения равна , чтобы умножить два разрядных целых числа (игнорируя факторы ). Битовая сложность для деления (до бит точности) одинакова. Следовательно, каждая итерация может быть вычислена в операциях . Умножая на итераций, мы находим, что общее время выполнения для вычисления квадратного корня до битов точности равно . Это субквадратичное.blglgbb O (nlgn)O(lgn)n O (n(lgn)2)O (blgb)blglgbbO (nlgn)O(lgn)nO (n(lgn)2)

Я думаю, что более тщательный анализ показывает, что это можно улучшить до времени выполнения (принимая во внимание, что нам нужно знать только каждый с точностью до бит точности, а не бит точности). Тем не менее, даже более базовый анализ уже показывает время выполнения, которое явно субквадратично.xjjnO (nlgn)xjjn

DW
источник
В двоичном коде также есть отличное начальное предположение, использующее тождество . Вместо того, чтобы вычислять журнал, можно приблизить к числу цифр в . Например, . войти 2 хx1/2=21/2log2xlog2xlog 2 101011 6xlog21010116
Ник Алджер
@DW: Но разве мы не ищем целочисленный квадратный корень? Если вы выполняете итерацию метода Ньютона, используя только целочисленную арифметику, то нам нужно дополнительное обоснование для утверждения , не так ли? В противном случае мы уже принимаем достаточно большую точность ... Извините, если я упускаю что-то очевидное. O(logn)
Арьябхата
@DW:«Скорость сходимости для метода Ньютона» не будет квадратичной, если , и я не знаю, что происходит для значений , которые не являются отрицательные реалы. Ваша оценка битовой сложности умножения является более жесткой, чем предполагает следующее замечание . Кроме того, нам «нужно знать каждый с точностью до« «битов точности». сc=0c2xj2j
@Aryabhata:Мы не совсем "ищем целочисленный квадратный корень"; мы ищем "пол квадратного корня". Вы правы в вопросе целочисленной арифметики, хотя те же битовые сложности имеют место для операций с плавающей запятой.
1
@RickyDemer, да, является частным случаем, потому что тогда корень имеет кратность 2, но когда , корень имеет кратность 1 , так метод Ньютона имеет имеет квадратичную сходимость. Я предполагаю, что никто не использовал бы метод Ньютона для вычисления квадратного корня из (потому что квадратный корень из нуля, очевидно, равен нулю). Так что ты пытаешься сказать? Является ли ваш комментарий тривиальным комментарием, который адресован путем добавления чего-то к моему ответу, который говорит «особый случай - квадратный корень из нуля», или здесь есть что-то более глубокое, чего я пропускаю? p ( x ) c > 0 c = 0c=0p(x)c>0c=0
DW
6

Одна из проблем метода Ньютона заключается в том, что он требует операции деления в каждой итерации, которая является самой медленной базовой целочисленной операцией.

Метод Ньютона для обратного квадратного корня, однако, не делает. Если - это число, для которого вы хотите найти , выполните итерацию:1x1x

ri+1=12ri(3xri2)

Это часто выражается как:

d i = 1 - w i x r i + 1 = r i + r i d i

wi=ri2
di=1wix
ri+1=ri+ridi2

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

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

Во-первых, давайте изменим масштаб :x

x=22ex

где мы хотели бы, чтобы было больше, но близко к . Если мы запустим вышеупомянутый алгоритм на вместо , мы найдем . Тогда . 1 x x r = 1x1xxr=1xx=2erx

Теперь давайте разделим на мантиссу и показатель степени:r

ri=2eiri

где - целое число Интуитивно , что представляет точность ответа. e iriei

Мы знаем, что метод Ньютона примерно удваивает количество точных значащих цифр. Таким образом, мы можем выбрать:

ei+1=2ei

С небольшой манипуляцией мы находим:

ei+1=2ei
wi=ri2
xi=x22eei+1
di=2ei+1wixi2ei+1
ri+1=2eiriridi2ei+1

На каждой итерации:

xrix2e+ei

В качестве примера, давайте попробуем вычислить квадратный корень из . Мы знаем, что ответ . Обратный квадратный корень равен , поэтому мы установим (это масштаб проблемы) и для нашего начального предположения выберем и . (То есть мы выбираем для нашей начальной оценки .)x=263231212231e=31r0=3e0=23412

Затем:

e1=4,r1=11
e2=8,r2=180
e3=16,r3=46338
e4=32,r4=3037000481

Мы можем решить, когда прекратить итерации, сравнив с ; если я правильно рассчитал, должно быть достаточно хорошим. Мы остановимся здесь и найдем:eieei>2e

2633037000481×263231+32=3037000481

Правильный целочисленный квадратный корень - , поэтому мы довольно близки. Мы можем сделать еще одну итерацию или оптимизированную последнюю итерацию, которая не удваивает . Детали оставлены в качестве упражнения.3037000499ei

Чтобы проанализировать сложность этого метода, обратите внимание, что умножение двух разрядных целых чисел требует операций. Тем не менее, мы упорядочили вещи так, чтобы . Таким образом, умножение для вычисления умножает два числа бит, чтобы получить число -бит, а два других умножения умножают два числа -бит, чтобы получить -битный номер.bO(blogb)ri<2eiwieiei+1ei+12ei+1

В каждом случае число операций на итерацию равно , и требуется итераций. Окончательное умножение порядка операций. Таким образом, общая сложность состоит из операций, которые субквадратичны по числу битов в . Это помечает все коробки.O(eilogei)O(loge)O(2elog2e)O(elog2e)x

Однако этот анализ скрывает важный принцип, который должен иметь в виду каждый, кто работает с большими целыми числами: поскольку умножение является суперлинейным по количеству битов, любые операции умножения должны выполняться только с целыми числами, которые имеют приблизительную величину текущей точности (и Я мог бы добавить, что вы должны попытаться умножить числа вместе, которые имеют схожий порядок). Использование целых чисел больше, чем это пустая трата усилий. Постоянные факторы имеют значение, а для больших целых они имеют большое значение.

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

Псевдоним
источник
Это отличные вещи. Однако один комментарий: не является ли битовая сложность деления асимптотически примерно такой же, как битовая сложность умножения? То есть вы говорите о чем-то, что дает постоянное улучшение фактора, а не асимптотическое улучшение, верно? Это не совсем понятно из твоего ответа.
DW
Вы говорите, что умножение двух битных целых чисел требует битовых операций. Я думаю, что правильный ответ - что-то вроде (верно?). Возможно, вы захотите указать, что вы игнорируете множители логарифмического логарифма (например, помещая тильду поверх своего большого О или чего-то еще). O ( b lg b ) O ( b lg b ( lg l g b ) O ( 1 ) )bO(blgb)O(blgb(lglgb)O(1))
DW
1
@DW:Нет, он говорит, что «умножение двух разрядных целых чисел требует операций». Слово "бит" появляется только один раз; иначе я бы уже указывал на это. O ( b log b )bO(blogb)
Это вопрос постоянных факторов, да. Лучшие алгоритмы деления больших целых чисел используют метод, очень похожий на весь алгоритм, такой как итерация Ньютона-Рафсона и удвоение эффективной точности на каждой итерации. Петля Ньютона-Рафсона внутри петли Ньютона-Рафсона накапливает постоянные множители! Рики Демер прав; Я думал в слове модель RAM. Я, наверное, должен был упомянуть об этом.
псевдоним