2014 dxdy logo

Научный форум dxdy

Математика, Физика, Computer Science, Machine Learning, LaTeX, Механика и Техника, Химия,
Биология и Медицина, Экономика и Финансовая Математика, Гуманитарные науки




Начать новую тему Ответить на тему На страницу Пред.  1, 2, 3  След.
 
 Re: Квадратный корень на целочисленной арифметике
Сообщение28.05.2015, 11:58 
Аватара пользователя


26/05/12
1534
приходит весна?
Andrey A, в связи с вашим предложением применять цепные дроби для нахождения корня меня мучает вопрос: а какова средняя вычислительная сложность нахождения первого нетривиального решения уравнения Пелля для чисел в диапазоне $0...2^N-1$ ? А так же какова сложность проверки того факта, что эти самые нетривиальные решения вообще есть?

 Профиль  
                  
 
 Re: Квадратный корень на целочисленной арифметике
Сообщение28.05.2015, 17:06 
Заслуженный участник
Аватара пользователя


21/11/12
1870
Санкт-Петербург
B@R5uk в сообщении #1020671 писал(а):
... что эти самые нетривиальные решения вообще есть?

Тривиальное решение одно: $1^2-m\cdot 0^2=1$, а наличие бесконечной последовательности нетривиальных решений для любого $m\neq c^2$ - факт установленный. Доказательство следует из самой классической процедуры разложения квадратного радикала и конечного числа кв. вычетов $z_i=p_i^2-mq_i^2<m/2, где $\frac{p_i}{q_i}$ - $i$-ая подходящая дробь. Цепная дробь периодична, последовательность вычетов - тоже, а единица заложена уже в "нулевом" тривиальном решении. Вычислительная сложность - об этом как раз Вас хотелось поспрашивать. Я, простите, не математик. Знаю, что числа вида $(ab)^2\pm b$ и $(ab)^2\pm 2b$ имеют короткий период в 2 - 4 знака, но для общего случая статистики не встречал. Предполагаю, что процедура нахождения 1-го решения более затратная, чем хотелось бы, иначе использовалась бы мат. пакетами (хотя - кто знает, тема-то не популярная). Но в любом случае это разовая акция, а далее выигрыш может быть ощутимым, поскольку для заданной точности любой другой метод задействует бОльшие знаменатели, будь это степени десятки или что-то еще. Ну и предложенный метод разложения. Новая хау. Если это требует разъяснений - я готов, там такая тройная рекуррентная последовательность образовывается. Сложение и умножение относительно малых целых чисел, больше ничего. Excel сразу отвечает, возможно иллюзия. Получив последовательность $a_1,a_2,...,a_k$, собираем искомую простую дробь: $u_0=a_k;u_1=a_{k-1}a_k+1;...;u_i=a_{k-i}u_{i-1}+u_{i-2};...$ (два последних члена).
Для грубых приближений смысла нет огород городить, а для точных - вроде бы перспективно.

 Профиль  
                  
 
 Re: Квадратный корень на целочисленной арифметике
Сообщение29.05.2015, 10:20 
Заслуженный участник
Аватара пользователя


21/11/12
1870
Санкт-Петербург
Andrey A в сообщении #1020357 писал(а):
Вычисления ведутся до $z_k=1$, и простая дробь восстанавливается от последнего знака.

P.S. Лучше так: Вычисления ведутся до $z_k=\pm 1$, и простая дробь восстанавливается от последнего знака.
В сущности ведь используется тождество $(p^2-mq^2)^2=(p^2+mq^2)^2-m(2pq)^2$. Если $p^2-mq^2=z$ или $mq^2=p^2-z$, делая соответствующие подстановки получаем $(2p^2-z)^2-m(2pq)^2=z^2$. Для $z_k=-1$, а это не редкий случай, получаем первую дробь $\frac{P_1=2p_k^2+1}{Q_1=2p_kq_k}$, и далее по известной схеме. При больших $P_1,Q_1$ возведение в квадрат (назовем это "схема 2") может оказаться "слишком" быстрой процедурой. Тогда можно использовать схему из первого моего поста (схема 1). Перейти от 1-ой ко 2-ой можно всегда, но не обратно.

 Профиль  
                  
 
 Re: Квадратный корень на целочисленной арифметике
Сообщение02.06.2015, 12:15 
Заслуженный участник
Аватара пользователя


19/12/10
1546
slavav в сообщении #1020130 писал(а):
Я сравнил метод Ньютона/Герона и ваш. Ваш лучше - он делает двоичный логарифм итераций, а Ньютон/Герон больше. Не смотря на квадратичную сходимость последнего.

Вы уверены в правильности такой оценки?

Возьмём, например, 1024-битное число. По Вашей оценке, приведённый ТС алгоритм, выполнит порядка 1024 итераций.

В качестве первого приближения в алгоритме Герона возьмём произвольное 512-битное число. Пусть это приближение имеет только один верный бит. Квадратичная скорость сходимости означает, что после каждой итерации число верных битов удваивается, следовательно потребуется порядка 9 итераций чтобы получить 512 верных битов.

 Профиль  
                  
 
 Re: Квадратный корень на целочисленной арифметике
Сообщение02.06.2015, 14:17 
Заслуженный участник


26/05/14
981
whitefox в сообщении #1022716 писал(а):
slavav в сообщении #1020130 писал(а):
Я сравнил метод Ньютона/Герона и ваш. Ваш лучше - он делает двоичный логарифм итераций, а Ньютон/Герон больше. Не смотря на квадратичную сходимость последнего.

Вы уверены в правильности такой оценки?


Я просто проверил. Для чисел в районе $10^9$ алгоритм делает 18-19 итераций. Начальное значение выбиралось как $\frac{a}{2}$.

whitefox в сообщении #1022716 писал(а):
Возьмём, например, 1024-битное число. По Вашей оценке, приведённый алгоритм выполнит порядка 1024 итераций.

В качестве первого приближения в алгоритме Герона возьмём произвольное 512-битное число. Пусть это приближение имеет только один верный бит. Квадратичная скорость сходимости означает, что после каждой итерации число верных битов удваивается, следовательно потребуется порядка 8 итераций чтобы получить 512 верных битов.


Я не умею считать биты. Давайте оценим метод точно. Итерация:
$x_{i+1} = \frac{x_i + \frac{a}{x_i}}{2}$.

Обозначим $C_i = \frac{x_i}{\sqrt{a}}$. Получим:
$C_{i+1} = \frac{C_i + \frac{1}{C_i}}{2}$.

Ещё одно обозначение: относительная погрешность $\varepsilon_i = C_i - 1$. Получим:
$\varepsilon_{i+1} = \frac{{\varepsilon_i}^2}{2(1 + \varepsilon_i)}$.

Обычно, когда говорят "число верных битов удваивается" то имеют ввиду: $\varepsilon_{i+1} \leqslant {\varepsilon_i}^2$. Последнее утверждение гарантирует высокую скорость сходимости, если $\varepsilon_i$ заметно меньше единицы. А если $\varepsilon_i$ заметно больше единицы, то Герон уточняет только один бит. Это я и имел ввиду, говоря, что он медленно стартует. Он медленно стартует с неточного приближения.

 Профиль  
                  
 
 Re: Квадратный корень на целочисленной арифметике
Сообщение02.06.2015, 15:22 
Заслуженный участник
Аватара пользователя


19/12/10
1546
slavav в сообщении #1022761 писал(а):
Я просто проверил. Для чисел в районе $10^9$ алгоритм делает 18-19 итераций. Начальное значение выбиралось как $\frac{a}{2}$.

Ну что ж, давайте проверим.

Возьмём число 220491584. Это число было выбрано случайно, его битовый размер равен 28 бит.

В качестве первого приближения возьмем любое 14-битное число, например 16383 $(11 1111 1111 1111_2),$ то есть потребуется не более $\left(\lceil\log_2 14\rceil+1\right)=5$ итераций.

Первая итерация даст число 14920.
Вторая — 14849.
Третья — 14848.
Четвёртая — 14848, останавливаемся.

Итого потребовалось 4 итерации. Как у Вас получилось 18 - 19?

 Профиль  
                  
 
 Re: Квадратный корень на целочисленной арифметике
Сообщение02.06.2015, 15:39 
Заслуженный участник


26/05/14
981
whitefox в сообщении #1022786 писал(а):
slavav в сообщении #1022761 писал(а):
Начальное значение выбиралось как $\frac{a}{2}$.

 Профиль  
                  
 
 Re: Квадратный корень на целочисленной арифметике
Сообщение02.06.2015, 15:43 
Заслуженный участник
Аватара пользователя


19/12/10
1546
А зачем? Ведь квадратный корень содержит вдвое меньше цифр.
Неэффективное применение алгоритма ещё не говорит о его неэффективности.

 Профиль  
                  
 
 Re: Квадратный корень на целочисленной арифметике
Сообщение02.06.2015, 16:36 
Заслуженный участник


26/05/14
981
Пытаюсь оправдаться:
На языке высокого уровня трудно найти число бит в числе. Этому посвящена глава 5.3 в Алгоритмических трюках младшего Уоррена.
Там же в главе 11.1 описаны алгоритмы поиска корня на целых числах. Половина главы относится к быстрому получению первого приближения (после этого достаточно 5 делений для 32-битовых целых).

 Профиль  
                  
 
 Re: Квадратный корень на целочисленной арифметике
Сообщение16.06.2015, 21:28 


09/02/09
90
Novosibirsk
slavav в сообщении #1022813 писал(а):
На языке высокого уровня трудно найти число бит в числе.

Выделить старший ненулевой бит несложно, сложнее быстро посчитать в какой позиции он оказался.
У процессоров с архитектурой ARM, впрочем, есть очень полезная команда CLZ:
Цитата:
Count Leading Zeros returns the number of binary zero bits before the first binary one bit in a value.

 Профиль  
                  
 
 Re: Квадратный корень на целочисленной арифметике
Сообщение16.06.2015, 22:50 
Заслуженный участник


20/08/14
11065
Россия, Москва
У процессоров x86 начиная с 80386 есть такая же команда BSR. Не самая быстрая, но быстрее вычислений на С. Да и вообще почти на всех процессорных архитектурах есть. Вот в микроконтроллерах редкость.

 Профиль  
                  
 
 Re: Квадратный корень на целочисленной арифметике
Сообщение03.07.2015, 11:55 


24/05/09

2054
Квадратный корень.

1. Считываем входное число при osn = 100. В дальнейшем считаем что osn = 10.
2. Подбираем цифру, квадрат которой не превосходит старший разряд входного числа. Это есть первая цифра ответа.
3. Найдем разницу между старшим разрядом входного числа и квадратом найденной цифры. Ответ припишем слева ко второму
разряду входного числа. Получим некое число A.
4. Удвоим имеющуюся часть ответа, получим число a. Подберем наибольшую цифру x, такую что произведение чисел ax
(a - десятки, x - единицы) и x не превосходило А. Число x - вторая цифра ответа.
5. Число A - ax * x припишем слева к третьему разряду входного числа. Получим некое число B. Повторяем итерации
4-5 до тех пор пока не закончатся разряды во входном числе.

Более нагляднее рассмотрим пример на основе числа 138384.

1) Храним число как 84-83-13
2) 3*3 <=13. Текущий ответ: 3
3) 13-3*3 = 4. Припишем 4 к 83. Получаем А = 483.
4) a = 3*2 = 6. 6x * x <=483. x = 7. Текущий ответ: 37.
5)483-67*7 = 14. Припишем 14 к 84. Получаем B = 1484.
6). b = 37*2 = 74. 74x*x <=1484. x = 2. Текущий ответ = 372. Он же и окончательный.

http://cppalgo.blogspot.co.uk/2010/09/blog-post_11.html

 Профиль  
                  
 
 Re: Квадратный корень на целочисленной арифметике
Сообщение15.01.2020, 23:41 
Аватара пользователя


26/05/12
1534
приходит весна?
Снова возникла необходимость считать квадратный корень, но уже на платформе, где деление как операция присутствует и, кроме того, более-менее быстро выполняется. Теперь метод Ньютона для целочисленного корня уже не выглядит такой нелепостью, а больше похож на заманчивую перспективу.

Однако, в процессе реализации столкнулся с проблемой останова. Дело в том, что если пытаться завершить процедуру при совпадении двух последних приближений, то есть шанс, что она никогда не закончится: приближения будут прыгать туда-сюда на единичку. Такое происходит, если пытаться вычислить корень из чисел вида $N^2-1$. Результатом должно быть число N, но когда это число подставляется в рекуррентную формулу, то получается вот что: $${{x}_{k-1}}=N$$ $${{x}_{k}}=\left\lfloor \frac{1}{2}\left( {{x}_{k-1}}+\left\lfloor \frac{{{N}^{2}}-1}{{{x}_{k-1}}} \right\rfloor  \right) \right\rfloor =\left\lfloor \frac{1}{2}\left( N+\left\lfloor \frac{{{N}^{2}}-1}{N} \right\rfloor  \right) \right\rfloor =\left\lfloor \frac{1}{2}\left( N+N-1 \right) \right\rfloor =N-1$$ $${{x}_{k+1}}=\left\lfloor \frac{1}{2}\left( {{x}_{k}}+\left\lfloor \frac{{{N}^{2}}-1}{{{x}_{k}}} \right\rfloor  \right) \right\rfloor =\left\lfloor \frac{1}{2}\left( N-1+\left\lfloor \frac{{{N}^{2}}-1}{N-1} \right\rfloor  \right) \right\rfloor =\left\lfloor \frac{1}{2}\left( N-1+N+1 \right) \right\rfloor =N$$
Вот, думаю теперь, как побороть эту проблему.

 Профиль  
                  
 
 Re: Квадратный корень на целочисленной арифметике
Сообщение16.01.2020, 01:29 
Заслуженный участник


20/08/14
11065
Россия, Москва
B@R5uk
Сохраняйте значение предыдущей итерации и проверяйте условие $|x_{k}-x_{k-1}| \le 1$. Где-то видел что можно проверять смену знака производной (например что $x_k<x_{k-1}$), но это кажется зависит от выбора начального приближения.
Какое из двух значений принять за корень легко проверить обратным умножением, если получили больше исходного, значит другое.

Кстати не вижу особой проблемы в выборе начального приближения для не слишком больших чисел на ЯВУ, уж выполнить 15-20 сравнений всяко быстрее даже одного деления (на x86/x64).

-- 16.01.2020, 01:38 --

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

 Профиль  
                  
 
 Re: Квадратный корень на целочисленной арифметике
Сообщение16.01.2020, 01:48 
Аватара пользователя


26/05/12
1534
приходит весна?
На stackoverflow нашёл вариант решения:

Используется синтаксис Java
public static BigInteger sqrt(BigInteger x) {
    BigInteger div = BigInteger.ZERO.setBit(x.bitLength()/2);
    BigInteger div2 = div;
    // Loop until we hit the same value twice in a row, or wind up alternating.
    for(;;) {
        BigInteger y = div.add(x.divide(div)).shiftRight(1);
        if (y.equals(div) || y.equals(div2))
            return y;
        div2 = div;
        div = y;
    }
}
 


Я, правда, для обычных чисел хотел этот алгоритм приспособить. И я решил проверять условие $0\le A-B^2\le 2B$, где A — число, из которого вычисляется корень, а B — очередное приближение. Таким образом не надо делать две лишних итерации, но требуется дополнительное умножение на каждом шаге. В конце концов (long)Math.sqrt(argValue) оказалось быстрее.

Dmitriy40 в сообщении #1435417 писал(а):
Кстати не вижу особой проблемы в выборе начального приближения для не слишком больших чисел

Можно составить таблицу размером где-нибудь 64-256 чисел (хотя первая четверть не понадобится особо) и из неё выбирать первые несколько бит начального приближения. Делается это с помощью сдвигов как-то так:

Используется синтаксис Java
argBits = 64 - Long.numberOfLeadingZeros(argValue);
shiftSize = (argBits - refBits + 1) / 2;
tmpValue = refRoots[(int) (argValue >> (2 * shiftSize))] << shiftSize;
 


По идее начальный запас точности позволяет сократить количество итераций. Хотя установить в единицу первый бит результата — это уже хороший запас точности.

 Профиль  
                  
Показать сообщения за:  Поле сортировки  
Начать новую тему Ответить на тему  [ Сообщений: 40 ]  На страницу Пред.  1, 2, 3  След.

Модераторы: Karan, Toucan, PAV, maxal, Супермодераторы



Кто сейчас на конференции

Сейчас этот форум просматривают: нет зарегистрированных пользователей


Вы не можете начинать темы
Вы не можете отвечать на сообщения
Вы не можете редактировать свои сообщения
Вы не можете удалять свои сообщения
Вы не можете добавлять вложения

Найти:
Powered by phpBB © 2000, 2002, 2005, 2007 phpBB Group