Квадратный корень на целочисленной арифметике

[B@R5uk]

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

Если кто-то встречал/писал что-то подобное, буду очень благодарен ссылочке/комментариям/идеям.

[AlexDem]

B@R5uk в сообщении #1019669 писал(а):

алгоритм вычисления квадратного корня, использующий исключительно целочисленную арифметику

Задача недоформулирована. Результат тоже целочисленным должен быть, дробные разряды отбрасываются? Или все вычисления вести над целыми, а в конце - упаковать? Такое ощущение, что если и есть такой алгоритм, то сильно нетривиальный.

[B@R5uk]

Ну, я поторопился создавать новую тему, поиск навёл на следующее:

uint32_t SquareRoot(uint32_t a_nInput)
{
    uint32_t op  = a_nInput;
    uint32_t res = 0;
    uint32_t one = 1uL << 30; // The second-to-top bit is set: use 1u << 14 for uint16_t type; use 1uL<<30 for uint32_t type


    // "one" starts at the highest power of four <= than the argument.
    while (one > op)
    {
        one >>= 2;
    }

    while (one != 0)
    {
        if (op >= res + one)
        {
            op -= res + one;
            res += one << 1;
        }
        res >>= 1;
        one >>= 2;
    }
    return res;
}

Это как раз тот алгоритм, что мне нужен. Более того, когда-то в детстве я этот же алгоритм изобрёл сам просто переложив в двоичную систему школьный метод вычисления квадратного корня в столбик.

AlexDem в сообщении #1019678 писал(а):

Задача недоформулирована. Результат тоже целочисленным должен быть, дробные разряды отбрасываются?

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

[AlexDem]

Мы в школе вроде в столбик корни не считали - не знал о таком. Интересно.

[slavav]

Кажется, метод Ньютона в целочисленной арифметике сходится к $\left\lfloor\sqrt{a}\right\rfloor$.

[whitefox]

Точно так, сходится с квадратичной скоростью. Почитать об этом можно в книге "Алгоритмические трюки для программистов" Генри Уоррена мл.

Применительно к квадратному корню, метод Ньютона ещё называют алгоритмом Герона.

[B@R5uk]

В методе Ньютона используется деление. А деление — это очень медленная операция, если нет аппаратной поддержки. Поэтому такой метод имеет серьёзные ограничения. Кроме того, зависимость времени выполнения от значения исходных данных в некоторых приложения может быть недопустима. Алгоритм, приведённый выше, свободен от первого недостатка полностью и освобождается от второго при некоторой доработке (надо избавится от первого цикла и слегка модернизировать содержимое второго).

[slavav]

Я сравнил метод Ньютона/Герона и ваш. Ваш лучше - он делает двоичный логарифм итераций, а Ньютон/Герон больше. Не смотря на квадратичную сходимость последнего. Думаю, что квадратичная сходимость проявляется когда вы близко к корню, А в целом случае у Ньютона/Герона нет "времени" чтобы обогнать ваш метод.
Если судить только по числу итераций, то ваш метод может оказаться эквивалентен двоичному поиску. Но он лучше и двоичного поиска, так как не умножает. Чистая победа.

[B@R5uk]

Это не мой метод. Сам не проверял, но говорят, что

Цитата:

...алгоритм был описан аж в 1945 в Von Neumann J. "First Draft of a Reaport on the EDVAC"

[Andrey A]

B@R5uk в сообщении #1019669 писал(а):

Интересует алгоритм вычисления квадратного корня, использующий исключительно целочисленную арифметику.

Логично предположить что и результат аппроксимации выражается простой дробью. Тогда трудно представить для фиксированного $m$ что-то быстрее последовательности $\frac{1}{0};\frac{P_1}{Q_1};...;\frac{P_{n+1}=KP_n-P_{n-1}}{Q_{n+1}=KQ_n-Q_{n-1}};...$, где пара $P_1,Q_1$ - первое решение уравнения Пелля, т.е. $P^2-mQ^2=1$, $K=2P_1$.
Для разложения $\sqrt{m} \approx a_1,a_2,...,a_k=\frac{P_1}{Q_1}$ есть доморощенный метод. В виде таблицы:

$\begin{matrix} \No & 1 & 2 & i+1\\ a & \left \lfloor \sqrt{m} \right \rfloor & \left \lfloor -\frac{2a_1}{z_1} \right \rfloor & \left \lfloor \frac{a_1+\left|b_i \right|}{\left|z_i \right|} \right \rfloor\\ b & a_1 & a_2z_1+a_1 & a_{i+1}z_i+b_i\\ z & a_1^2-m & a_2(b_1+b_2)+1 & a_{i+1}(b_i+b_{i+1})+z_{i-1} \end{matrix}$

Тут $a_i$ - собственно знаки непрерывной дроби, $b_i$ - некая вспомогательная переменная, $z_i$ - остатки вида $p_i^2-mq_i^2$, для которых имеет место рекуррентная зависимость. Вычисления ведутся до $z_k=1$, и простая дробь восстанавливается от последнего знака. Умножение присутствует, но только для целых чисел $<2\sqrt{m}$. Можно впрочем сократить количество вычислений и до полупериода, хотя тут тонкости начинаются. Думаю, Вам это не подойдет, но интересно бы оценить сложность алгоритма.

[venco]

B@R5uk в сообщении #1019669 писал(а):

Интересует алгоритм вычисления квадратного корня, использующий исключительно целочисленную арифметику.

А зачем обязательно целочисленная арифметика? На современных Интеловских процессорах встроенный корень быстрее, ну и скорость практически фиксированная.

[B@R5uk]

Однокристальные компьютеры как правило не имеют сопроцессора. Так что практических приложений много.

[ArtDen]

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

[Andrey A]

Andrey A в сообщении #1020357 писал(а):

... последовательность $\frac{1}{0};\frac{P_1}{Q_1};...;\frac{P_{n+1}=KP_n-P_{n-1}}{Q_{n+1}=KQ_n-Q_{n-1}};...$, где пара $P_1,Q_1$ - первое решение уравнения Пелля

P.S.
Еще быстрее такая: $\frac{P_1}{Q_1};...;\frac{P_{n+1}=2P_n^2-1}{Q_{n+1}=2P_nQ_n}$.

$\sqrt{2}\approx \frac{3}{2};\frac{17}{12};\frac{577}{408};\frac{665857}{470832};...$ Для всех пар выполняется $P^2-2Q^2=1$.

[=SSN=]

B@R5uk в сообщении #1020057 писал(а):

В методе Ньютона используется деление. А деление — это очень медленная операция, если нет аппаратной поддержки. Поэтому такой метод имеет серьёзные ограничения.

Вы просто не умеете его готовить.. ;)

По методу Ньютона нужно вычислять $\frac{1}{\sqrt{m}}$, после чего полученный результат нужно умножить на $m$.