Проверка что число - куб

[wrest]

guryev в сообщении #1331219 писал(а):

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

Это ме-е-е-е-дленно.

[iifat]

wrest в сообщении #1331218 писал(а):

Они должны быть взаимно-простыми?

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

-- 09.08.2018, 01:14 --

wrest в сообщении #1331222 писал(а):

Это ме-е-е-е-дленно

Не факт. Стоит проверить.

[wrest]

mihaild
Для n до 10^4 перебором нашлось n=9828 с максимальной "эффективностью": 31
Запустил перебор до 10^5, похоже это очень надолго, как вы нашли ваше 575757 - случайно набрали на клавиатуре?

-- 08.08.2018, 18:23 --

iifat в сообщении #1331224 писал(а):

Не факт. Стоит проверить.

У меня так примерно и работает сейчас. Я верю программистам PARI/GP и использую их встроенную функцию ispower которая как-то проверяет является ли число нужной степенью.

[mihaild]

Если $a = bc$, где $gcd(b, c) = 1$, то количество кубических вычетов по модулю $a$ равно произведению количеств кубических вычетов по модулям $b$ и $c$.

[Dmitriy40]

wrest
Ещё идейка по ускорению вычисления $\sqrt[3]{x}$. Ну таблицу кубов для всех интересующих $t$ уже посчитали заранее. Теперь считаем $a=\sqrt{x}$ (а он считается почти со скоростью деления), после чего двоичным поиском по таблице кубов ищем подходящее значение $t$ в диапазоне $(0,a=\sqrt{x})$. Для чисел $t\le10^4$ диапазон поиска миллион и достаточно всего 20 обращений в массив (или возведения в куб если памяти жалко). Это не так уж долго даже с умножениями без массива. Для $t\le10^n$ достаточно $5n$ возведений в куб.
И разумеется это делать лишь если проверка по модулю не отвергла $x$.
Вот будет забавно если ispower примерно так и действует ...

[wrest]

Dmitriy40 в сообщении #1331234 писал(а):

Ещё идейка по ускорению вычисления $\sqrt[3]{x}$.

Само вычисление корня кубического, напомню, не нужно. Вернее, нужно, но крайне редко. Нужна только проверка куб ли данное число, а вычисление-то так редко делается что хоть столбиком можно.

Dmitriy40 в сообщении #1331234 писал(а):

Вот будет забавно если ispower примерно так и действует ...

Я тоже надеюсь, что ispower действует максимально эффективно, но вряд ли использует таблицы предвычисленных кубов или "хеши".

-- 08.08.2018, 19:45 --

Dmitriy40 в сообщении #1331234 писал(а):

И разумеется это делать лишь если проверка по модулю не отвергла $x$.

Да, как проверку-то делать?
Вот есть у меня число 90909 и есть вектор v в котором перечислены все возможные остатки от деления кубов на это число, их 1365.
Теперь берем какое-нибудь (проверяемое) число, делим его на 90909 с остатком и проверяем не получился ли какой-нибудь из этих 1365?

Делаем вектор со всеми кубическими вычетами по модулю: 90909

v=vector(90909,i,i^3%90909)

Оставляем в векторе только неповторяющиеся остатки:

v=vecsort(v,,8)

Заводим функцию, проверяющую на потенциальную кубичность:

iscube(n)=my(d=90909);r=n%d;if(vecsearch(v,r),return(1),return(0))

Смотрим примерно 10^7 проверок

? for(i=10^6,10^7,iscube(i))
  ***   last result computed in 6,335 ms.
? for(i=10^6,10^7,ispower(i,3))
  ***   last result computed in 1,107 ms.

Что-то не так...
Повышаем ставки:

? for(i=10^9,10^9+10^7,iscube(i))
? ##
  ***   last result computed in 7,051 ms.
?  for(i=10^9,10^9+10^7,ispower(i,3))
? ##
  ***   last result computed in 875 ms.

Все равно фигня.
Конструкция без вызова функции, а роверки прямо по месту улучшает в два раза но все равно работает долго:

? for(i=10^9,10^9+10^7,vecsearch(v,i%90909))
? ##
 ***   last result computed in 2,090 ms.

[Dmitriy40]

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

(Оффтоп)

Но подумав ещё, думаю там всё же метод Ньютона, он быстрее должен сходиться. Проверка это не опровергает, для чисел до $2^{64}$ на x64 PARI ispower(x,3) выполняется за 700 тактов в среднем, что в начале, что в конце диапазона, что вполне похоже на полтора десятка делений.
Впрочем на битовое взвешивание с аккуратной реализацией даже больше похоже, как раз около 100 умножений по 5 тактов каждое плюс условие выбора бита. И удвоение времени при переходе границы $2^{64}$ тоже именно за это.

-- 08.08.2018, 20:07 --

wrest
Ну я же Вам процитировал практически готовый код (на С):

Dmitriy40 в сообщении #1331200 писал(а):

Там смысл в этих двух кусочках кода (привожу с сокращениями):

Используется синтаксис C++

//Заполнение массива флагов
const TInt CUBES_MOD = 575757;
std::vector<bool> cubesMod(CUBES_MOD, false);
for (t = 0; t < MAX_T; ++t) cubesMod[(t * t * t) % CUBES_MOD] = true;
//Использование массива флагов
if (!cubesMod[t_target % CUBES_MOD]) continue;//Такое t_target точно не является никаким кубом и извлекать корень не нужно

Т.е. словами, заполняем массив флагов так чтобы он по младшим цифрам числа указывал может ли куб какого-то числа иметь такие младшие цифры. Младшие цифры берутся не десятичные, а по модулю CUBES_MOD, при удачном выборе этого числа количество true в массиве будет сильно меньше количества false и для многих и многих t_target извлекать корень не понадобится.

Что тут непонятного? Объявляем массив длиной 575757, проходим по интересующим нас числам 0..MAX_T и расставляем в массиве флаг что число с таким остатком по модулю 575757 является кубом некоего числа. Это была инициализация.
При использовании вместо vecsearch просто обращаемся к массиву по индексу [x%CUBES_MOD] и если флаг сброшен, то не существует чисел, куб которых даёт такой остаток по данному модулю. Если же такие числа есть, то надо продолжить проверку и например тупо проверить ispower(x,3).

-- 08.08.2018, 20:27 --

Ну вот и PARI код для чисел до $10^4$ (лень ещё переменную вводить):

? MODS=575757;cubeMod=vector(MODS,i,false);for(i=0,10^4,cubeMod[i^3%MODS+1]=true);\\Инициализация
? for(i=10^10,10^10+10^7,if(cubeMod[i%MODS+1]==true,ispower(i,3)))
? ##
  ***   last result computed in 1,483 ms.
? for(i=10^10,10^10+10^7,if(true,ispower(i,3)))
? ##
  ***   last result computed in 905 ms.

Да, с массивом оказывается медленнее, значит ispower хорошо оптимизирована.

[mihaild]

iifat в сообщении #1331217 писал(а):

а в поле по модулю 575757 — в среднем 140

Это всё-таки не поле.
Ну и в среднем кубический корень имеет одно значение (т.к. возведение в куб - однозначная операция).

[wrest]

Dmitriy40 в сообщении #1331248 писал(а):

При использовании вместо vecsearch просто обращаемся к массиву по индексу [x%CUBES_MOD] и если флаг сброшен, то не существует чисел, куб которых даёт такой остаток по данному модулю.

А, то есть сначала все числа проверить на псевдокубовость и запомнить результат, а потом там рыться чтобы не делить с остатком...
Ну как-то оно не быстрее-то vecsearch вышло почему-то (а можете у себя с vecsearch) проверить?
Мне вот интересно, mihaild писал что у него проверка с 400 по 500 в пределах секунды. А у вас? Можете сравнить С++ а и PARI/GP?

Вот код на PARI/GP без предвычисления, для прохождения головоломки от 400 до 500

(Оффтоп)

k(z)=my(n=0,z5=z^5,xmax=sqrtnint(z5,4)+1,ymax=xmax,t=0);for(x=1,xmax,for(y=x,ymax,n=x^4+y^4;if(n>z5,next(2));if(ispower(z5-n,3,&t),print("x=",x," y=",y," z=",z," t=",t))))
for(i=400,500,k(t)

[Dmitriy40]

wrest в сообщении #1331227 писал(а):

Для n до 10^4 перебором нашлось n=9828 с максимальной "эффективностью": 31

Ещё $n=15561$ даёт хорошую эффективность, всего $315$ вхождений, эффективность $49{,}4$. Предыдущее $n=14364$ давало эффективность лишь $32{,}6$. Код для поиска (выводит только при увеличении эффективности):

best=0;for(n=2,10^6,s=vector(n);for(i=0,#s-1,s[i^3%n+1]=1);w=vecsum(s);if(n/w>best,best=n/w;printf("n=%d / %d = %0.1f\n",n,w,best)))

-- 08.08.2018, 21:39 --

wrest в сообщении #1331262 писал(а):

(а можете у себя с vecsearch) проверить?

Могу:

? v=vector(90909,i,i^3%90909);v=vecsort(v,,8);iscube(n)=my(d=90909);r=n%d;if(vecsearch(v,r),return(1),return(0));
? for(i=10^6,10^7,iscube(i))
? ##
  ***   last result computed in 32,402 ms.
? for(i=10^6,10^7,ispower(i,3))
? ##
  ***   last result computed in 515 ms.
? for(i=10^9,10^9+10^7,iscube(i))
? ##
  ***   last result computed in 36,146 ms.
? for(i=10^9,10^9+10^7,ispower(i,3))
? ##
  ***   last result computed in 561 ms.

Почему у меня более 30с вместо 6-7с как у Вас я не знаю.
C++ код проверить не могу.

-- 08.08.2018, 21:44 --

wrest
Кстати, если будете делать с операцией % в индексе векторов PARI, имейте в виду что индекс векторов начинается с 1, а не с 0 как в С/С++, и программа может налететь на ошибку если получится нулевой остаток, именно поэтому я в коде делаю +1 индексу массива.

[lel0lel]

В тестах простоты для определения является ли число полным квадратом сначала проверяют является ли оно квадратичным вычетом по нескольким небольшим простым модулям (благо - символ Лежандра вычислять довольно просто). Если число не является полным квадратом, тогда вероятность того, что оно пройдет $n$ таких проверок $1/2^n$.

По аналогии можно поступать и с определением является ли число кубом. Отличие от предложенного выше в том, что не потребуется создавать предварительно рассчитанные массивы в памяти. Выбираем простые вида $6k+1$, например набор $\{7, 13, 19, 31\}$. Для данного $x$ и каждого простого $p$ из выбранного набора проверяем является ли $x$ кубическим вычетом по модулю $p$. Поскольку удобной процедуры (типа вычисления символа Лежандра с использованием Закона квадратичной взаимности), по-моему, не известно, то проверяем в лоб $x^{\frac{p-1}{3}}\equiv1\pmod{p}$. Если число не куб, то вероятность, что $n$ проверок будут пройдены равна $1/3^n$.

[Dmitriy40]

Мне кажется одно деление и одно обращение в массив, причём если уж выдаёт ложь, то точно не куб - быстрее множества умножений по модулю.

[lel0lel]

Это надо проверять. Не маленькая вероятность, что будет достаточно одной-двух проверок. Выше я не совсем аккуратно написал, что проверяем для всех простых из множества. Разумеется, если сравнение не выполнено для какого-то простого, то число кубом не является и его мучить не надо.

[wrest]

Dmitriy40
Кстати, у pari/gp открытый исходный код и можно же посмотреть где там ispower и как реализована, но я в этом не понимаю и найти не могу. Нашел файл с таки именем, но там текст из справки, которая выдается по ??ispower

[Dmitriy40]

И всё же метод с дополнительным массивом для проверки по остаткам имеет смысл. Написал примитивную (без особых оптимизаций) утилиту извлечения кубического корня для чисел до $2^{31}$ и сравнил её скорость без массива и с предварительной проверкой по массиву:

C:\TEMP>sqrt3.exe
n=15561 / 315 : 10780ms, 3994ms.

Это посчитаны кубические корни миллиарда чисел от одного до двух миллиардов. Т.е. 10 или 4 нс на число в среднем. Ускорения в 50 раз не видно, но в 2.5 раза есть. Думаю выигрыш съедается операцией деления, а 10-ти кратное битовое взвешивание с возведением в куб (всего 20 умножений) недостаточно долгое чтобы сильно перевесить деление.
PARI:

? for(i=10^9,2*10^9,ispower(i,3))
? ##
  ***   last result computed in 57,002 ms.

-- 09.08.2018, 00:32 --

Зато если резко замедлить вычисление корня (использовав int64 числа, x32 код работы с которыми после компилятора ужасен), то выигрыш от предпросмотра по массиву очень даже проявляется:

C:\TEMP>sqrt3.exe
n=15561 / 315 : 220726ms, 10593ms.

Условия те же, весь второй миллиард. Ускорение уже более 20 раз.