Пользователь
0,0
рейтинг
17 декабря 2015 в 12:13

Разработка → Вперед, на поиски палиндромов 3

После того, как вроде бы неплохой результат, полученный в предыдущей части, оказался лишь «локальным максимумом», я на некоторое время забросил задачку. Напомню условие:
«The decimal number 585 is 1001001001 in binary. It is palindromic in both bases. Find n-th palindromic number». Или, по-русски: «Десятичное число 585 в двоичной системе счисления выглядит как 1001001001. Оно является палиндромом в обеих системах счисления. Найдите n-й подобный палиндром».

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

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

Оригинальное объяснение и реализация используют префиксное дерево, но на мой взгляд, чуть более глубокое понимание математики алгоритма позволяет обойтись более простыми и быстрыми структурами. Думаю, что лучше всего разобрать алгоритм на примере.

Мы будем искать двоичные палиндромы среди десятичных палиндромов шириной 8. Исходных палиндромов ровно 9000: от 10000001 до 99999999.

Теперь возьмём 2 множества чисел:

  • А[90] = { 10000001, 11000011, 12000021, ..., 98000089, 99000099 }
  • B[100] = { 00000000, 00011000, 00022000, ..., 00988900, 00999900 }

Если описывать эти множества формально, то множество А является подмножеством десятичных палиндромов шириной 8, у которых средние 4 разряда — нули, а множество B состоит из 0, множества десятичных палиндромов шириной 2, умноженных на 103 и множества десятичных палиндромов шириной 4, умноженных на 102.

А если неформально, то множество А — это все возможные «края», а множество B — все возможные «серединки» десятичных палиндромов шириной 8. Очевидно, что множество всех десятичных палиндромов шириной 8 — это множество всех попарных сумм элементов множеств A и B.

for each (a in A) {
	for each (b in B) {
		palindrome = a + b;
	}
}

Далее, для краткости я буду использовать «a» вместо «элемент множества A», и «b» вместо «элемент множества B».

Теперь переведём элементы обоих множеств в двоичную систему счисления:

A
10000001 100110001001011010000001
11000011 101001111101100011001011
12000021 101101110001101100010101

98000089 101110101110101110011011001
99000099 101111001101001111100100011

B
00000000 000000
00011000 10101011111000
00022000 101010111110000

00988900 11110001011011100100
00999900 11110100000111011100

Я выделил по 6 верхних и нижних разрядов у всех a, и 6 нижних разрядов у всех b. Ширина 6 выбрана не случайно — это максимальная степень 2, не превышающая ширину «краёв» A = 102.

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

Например, для a = 10000001, все образованные от него десятичные палиндромы { 10000001, 10011001, 10022001, ..., 10988901, 10999901 } находятся в интервале [10000001, 11000011).

Это значит, что все десятичные палиндромы, образованные от a = 10000001 могут начинаться только со следующих 6 двоичных цифр:
100110 — первые двоичные цифры a = 10000001
100111
101000
101001 — первые двоичные цифры a = 11000011

Следовательно, чтобы быть двоичными палиндромами, все эти десятичные палиндромы могут заканчиваться только на обратную перестановку начальных двоичных цифр:
100110 -> 011001
100111 -> 111001
101000 -> 000101
101001 -> 100101

А учитывая, что сам а = 10000001 заканчивается на двоичные цифры 000001, то из всех возможных b нас интересуют только те, которые заканчиваются на двоичные цифры:
011001 — 000001 = 011000
111001 — 000001 = 111000
000101 — 000001 = 000100
100101 — 000001 = 100100

Только для таких b нужно проверять, является ли 10000001 + b двоичным палиндромом.

Используя данный подход, алгоритм для нахождения палиндромов в основаниях BASE0, BASE1 в интервале [1, N] может быть описан следующим образом:
Для каждой ширины W из [1, количество разрядов N в BASE0]
    Сгенерировать множества A и B, используя BASE0. Ширина «краёв» A WA = O(W/4), WA ≥ WB
    Перевести элементы A и B в BASE1
    Рассортировать элементы B по корзинам по последним цифрам в BASE1
    Для каждого a из множества A
        Для каждого X из множества возможных начальных цифр a + b
            Для каждого b, заканчивающегося на (X — конечные цифры a)
                Проверить, является ли a + b палиндромом в BASE1

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

  1. Множества A и B содержат O(N1/4) элементов.
  2. Для каждого a в среднем существует не более BASE0 вариантов X.
    (Мы изначально выбираем ширину интересующих нас начальных цифр в BASE1 так, чтобы получившееся из них число было не больше BASE0WA, а max(A) больше min(A) в BASE0 раз.)
  3. Для каждого X в среднем проверяется не более, чем BASE1 различных b.
    (X распределяется равномерно и почти не коррелирует с конечными цифрами a в BASE1. Значит каждая корзина с b выбирается равномерно, а всего таких корзин не более чем в BASE1 раз меньше, чем самих b.)
  4. Проверка на палиндром все равно занимает O(log(N)).

В общем, я предполагаю, что вычислительная сложность алгоритма O(N1/4 * log(N) * BASE0 * BASE1).

Первая же реализация этого алгоритма вполне предсказуемо была на пару порядков быстрее, а потратив еще немного времени на оптимизацию я довел время вычисления до чуть более 0.01 секунды, более чем в 1000 раз быстрее предыдыщего алгоритма. Этот результат наконец-то меня удовлетворил, но вполне естественно возникло желание проверить его на числах, бо́льших, чем 260. К этому времени я уже обнаружил интересующую меня последовательность в «Энциклопедии целочисленных последовательностей». В списке двойных палиндромов было 122 члена, самый большой из которых состоял из 131 бита. Это было впечатляюще, но и косвенно указывало на то, что никто пока не придумал алгоритм логарифмической сложности.

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

2131 / 260 = 271

271 * 1/4 < 218 = 262144

Учитывая 0.01 секунды, необходимых для предела в 260, выходило, что мой алгоритм должен был справиться с этой задачей примерно за 1 час! Я чувствовал какой-то подвох, но вызов был уже принят.

Продолжение следует.
Илья Никульшин @ilyanik
карма
5,0
рейтинг 0,0
Пользователь
Реклама помогает поддерживать и развивать наши сервисы

Подробнее
Спецпроект

Самое читаемое Разработка

Комментарии (38)

  • 0
    Предсказываю, что для 131 бита потребуется перебор 5^20 чисел (примерно 10^14).
    • 0
      Почему 520? Только первая цифра десятичного палиндрома нечетная, на остальные ограничений нет.
  • 0
    Я рассмотрел слегка другой алгоритм — рекурсивный.
    Когда мы выбрали K старших десятичных цифр палиндрома, мы почти всегда знаем K+1 его старших бит (обычно, почти 3*K, но для наших целей это неважно). Следовательно, знаем K+1 младший бит, а поскольку K младших цифр нам известно, то K+1-я цифра может принимать 5 разных значений — быть либо чётной, либо нечётной (в зависимости от ситуации). Так что число кандидатов увеличилось в 5 раз.
    В редких случаях, когда старшие K+1 бит неизвестны, нам придётся поставлять все 10 кандидатов на K+1-ю цифру, но половина из них сразу отсеется.
    Какой-то из пунктов 2 или 3 в Ваших рассуждениях вызывает сомнения — в нём (не знаю пока, в каком) должно быть (BASE0/BASE1)^(W/4) вариантов (или (BASE1/BASE0)^(W/4) — из текста непонятно, кто из них 2, а кто 10).
    • 0
      За 20 минут алгоритм нашёл числа до 10^25 (т.е. 80 первых членов последовательности). На поиск чисел до 10^40 ушло бы около 5 лет.
    • 0
      Я брал BASE0 > BASE1, но не думаю, что это принципиально.

      80 первых членов я нахожу за 3.5 секунды.
    • 0
      Это не строгое доказательство, а попытка объяснить реальные результаты :)

      Например, для ширины в 24 десятичных знака, и следовательно 900000 элементов А, проверено 3958921 палиндромов, в среднем по 4.4 различных b на каждый a.
      • 0
        проверено 3958921 палиндромов

        Действительно палиндромов? Или значений X? Палиндромов должно быть в 64 раза больше — ведь все числа из множества B в двоичной записи заканчиваются на 6 нулей. Или я чего-то совсем не понял.
        Какая там получается ширина края в двоичной записи — 19 или 20?
        • 0
          Ширина края — 19. Т.е. 1000000 палиндромов-серединок распределяются по 2^19 корзин (в среднем по 1.9 в каждой корзине).
          Палиндромы-края изменяются от 100000хххххххххххх000001 до 999999хххххххххххх9999999, и для каждого палиндрома-края проверяется всего несколько корзин.
          • 0
            К сожалению, заполнено только 2^13 корзин, по 122 в каждой…
            • 0
              Эээ..., а почему 2^13 корзин? :)

              Это же мой алгоритм, у меня 2^19? Откуда вы вообще взяли число 13?
              • 0
                Но вы же кладёте в них числа из множества B, которые в десятичной системе оканчиваются на 6 нулей, не так ли?
                • 0
                  А, из-за того, что 10 делится на 2?

                  Да, распределение по корзинам получается очень неравномерное. Но из-за равномерного выбора каждой из 2^19 корзин, в среднем выходит 1.9: на каждую корзину с ~122 элементами 63 пустых.
        • 0
          Например, первого палиндрома края я проверяю 5 корзин:

          100000_000000000000_000001 = 1010100101101000000_1011000111111000010100101011110110100000000000000000000001
          100001_000000000000_100001 = 1010100101101000100_0010101000100101111111111010011101111001011000011010100001
  • 0
    Ширина края — 19. Т.е. 1000000 палиндромов-серединок распределяются по 2^19 корзин (в среднем по 1.9 в каждой корзине).
    Палиндромы-края изменяются от 100000хххххххххххх000001 до 999999хххххххххххх9999999, и для каждого палиндрома-края проверяется всего несколько корзин.
  • 0
    Ну что же. В их последовательности пропущено 31-значное число :) Сможете найти?
  • 0
    Ещё они пропустили, по меньшей мере, два 39-значных числа (оба начинаются на 12...)…
  • 0
    Пропущенное 40-значное число: 1017421766189445102992015449816671247101
    То есть, их последнее число — не 122-е, а 126-е.
  • 0
    Ну, и ещё одно 40-значное число: 7114907950920173924554293710290597094117. На всё ушло 49 минут на одном ядре.
    • 0
      Круто! А как же 5 лет?
      • +1
        Чуток соптимизировал.
        Вот ещё:
        128 14327425216354951264146215945361252472341
        129 31636759764024794204540249742046795763613
        130 38090421176450233778487733205467112409083
        131 56545858306667087923432978076660385854565
        132 98801466348600079992129997000684366410889
        133 128795669673344381770077183443376966597821
        134 139035351443367699760067996763344153530931
        135 170341815153453197154451791354351518143071
        136 309612431907274418544445814472709134216903
        137 595943598626320807905509708023626895349595
        138 905357630732463833436634338364237036753509
        139 1816060344791285708869688075821974430606181
        140 3381578420704603001613161003064070248751833
        141 7342779513978827245484845427288793159772437
        142 9101547767757547725021205277457577677451019

        44-значных программа не нашла.
        Думаю, на этом остановлюсь (до очередного прогрессивного способа) — на следующий разряд ушли бы сутки счёта.
        • 0
          Мне кажется, вам надо принять эстафету «Вперед, на поиски палиндромов» :)

          Потому что ваш алгоритм уже быстрее O(N1/4)
          • 0
            До 40-значных включительно — O(N7/34), потом — O(N0.35).
            Но алгоритм отличается от вашего только шириной краёв в десятичной записи: я беру не 1/4, а 5/17 от длины числа (пока хватает памяти на серединки). Тогда в непустые корзины действительно попадает по одному числу. Ну, и края фильтруются по младшим битам.
            • 0
              Ну, с шириной краев и середины я тоже «играл», тем более, что потом всё равно кончается память. Но дальше вы по краям тоже ведь не просто проходите…
              • 0
                Память тратится только на середину. 10^8 чисел удержать можно (если хранить их как 128-битные). Края строятся динамически, их запоминать не нужно.
            • 0
              Мне кажется, что я придумал модификацию вашего алгоритма, позволяющую при том же количестве оперативной памяти увеличить количество десятичных разрядов N, для которых вычислительная мощность все еще равна O(N7/34) в 10/7 раз.

              Дело в том, что выбор верхних десятичных разрядов краёв влияет на четность всех цифр палиндрома, включая цифры серединок. А это значит, что можно разделить полный обход каждой конкретной ширины N на 2X подзадач. Каждая подзадача включает в себя прегенерацию 10Y*5X серединок (X нижних разрядов фиксированной четности) и обход только тех краёв, для которых Х нижних разрядов серединок имеют нужную четность.

              Благодаря свойствам алгоритма обхода краёв, ненужные ветки будут отсекаться относительно рано, поэтому суммарное время обхода всех краёв по всем подзадачам будет лишь незначительно больше, чем в случае, когда у нас есть достаточно оперативной памяти для всех 10Y + X серединок сразу.
        • 0
          Но краев-то 5*10^13 (для ширины 44). Только на их генерацию должно было уйти много часов, даже при скорости генерации 10^9 в секунду.
          • 0
            С учётом предварительной фильтрации — всего лишь 5^14. Ушло 3-4 часа.
            Если бы я не поленился и написал 192-битные числа, получилось бы ещё быстрее.
            • 0
              Так вы профильтровали 5*10^13 краев, или сразу сгенерировали 5^14?
              • 0
                Можно сказать, что сразу — фильтрация шла для каждой очередной цифры.
                • 0
                  Здорово!

                  Но, если я не ошибаюсь, ваше ускорение не применимо в случае взаимно простых оснований.
                • 0
                  А вы как-то оптимизировали/сортировали структуру, содержащую серединки?

                  Я только что посмотрел — на «случайный» доступ к этой огромной структуре тратится 70-80% всего времени, в разы больше чем на всю математику.
                  • 0
                    У меня основное время уходит на разворот длинных чисел в двоичной записи.
                    Структура, содержащая серединки, устроена просто:
                    struct{
                    UInt128 val;
                    int next;
                    }
                    где next — индекс в массиве следующего элемента в списке. next=0 — ячейка пуста, next=-1 — этот элемент последний. Первые 2^X элементов — корзины (начала списков), дальше — свободные места для продолжения списков. Когда массив оказался длиннее 2ГБ, пришлось переделать его на массив массивов.
                    Если захочу считать дальше, переделаю его на файл и буду хранить серединки, упорядоченные начиная с младшего бита. Понятие корзины при этом исчезнет, всё будет делаться сопоставлениями двух отсортированных списков. Но это потом.
                    • 0
                      А вы использовали этот прием:
                      graphics.stanford.edu/~seander/bithacks.html#ReverseParallel
                      ?
                      • 0
                        Нет, я даже инверсией по табличке не пользовался. Поскольку длина в двоичной записи не кратна 8, пришлось бы совмещать её со сдвигом длинного числа. Но было лень :(
        • +1
          143 557019370941199612258585852216991149073910755
          144 706400325289926993853434358399629982523004607
          145 903240486809073407096959690704370908684042309
          146 5318935032766104711807997081174016672305398135
          147 9335388324586156026843333486206516854238835339
          

          45 и 46-значные: 5:30 и 4:45 часа соответственно.

          Огромное спасибо за вашу прекрасную оптимизацию и конкуренцию!
        • +1
          Я обновил список на сайте «Энциклопедии целочисленных последовательностей». Сейчас он доступен здесь, но после утверждения правки он скорее всего заменит предыдущий список.

          Также, я указал ваше имя (Andrey Astrelin) — без вашей идеи я бы дошел лишь до 10^41.
        • 0
          47-значные, 9:45 часов:

          148 12328899531897059171731113717195079813599882321
          149 16422159001061376847917371974867316010095122461
          150 36755925874534219715185758151791243547852955763
          151 56243939994005432191600400619123450049993934265
          152 72224512737657344148643434684144375673721542227
          153 72830033748815722240681118604222751884733003827
          154 74953152103456169263148284136296165430125135947
          155 78737696079148631316169196161313684197069673787
          

          И это уже с первой робкой попыткой prefetch, даже с которой на хаотическое обращение к памяти уходит 75-80% времени.
        • 0
          48-значные, 0:32 часа:

          156 102735644379963218861031130168812369973446537201
          157 502186653128032879493361163394978230821356681205
          158 597817365870480462496846648694264084078563718795

          49-значные, 2:41 часа:

          159 1873893355166996611906735376091166996615533983781
          160 7109242970502610649115178715119460162050792429017
          161 9142687306518774468290258520928644778156037862419

          50-значных, 2:30 часа: нет

          Основные «ускорители»: 109 серединок, 8 потоков.

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