2014 dxdy logo

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

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




Начать новую тему Ответить на тему На страницу 1, 2  След.
 
 Остаток от деления больших биномиальных коэффициентов
Сообщение16.03.2011, 08:07 


16/03/11
31
Решаю олимпиадную задачку по программированию, для её решения нужно быстро считать модули биномиальных коэффициентов.

Вот исходные данные:
* Нужно подсчитать $C^n_k (mod\ p)$, где
* $0 < k < n < 10^5$
* $0 < p < 10^5$
* Нужно подсчитать 500 коэффициентов за 3 секунды (500 - это количество входных данных), а лучше -- быстрее, т.к. выполняются другие операции
* Для справки: за секунду выполняется $10^9$ сложений, $10^8$ умножений и $10^7$ делений.

Обращаю ваше внимание, что не нужно считать весь коэффициент, нужно посчитать только остаток от деления. Это как-то должно ускорить вычисления, но как??

Пока придумал такой алгоритм: рассмотрим способ вычисления $C^n_k = {{n \times \ldots \times (n - k + 1) } \over {k!}}$. Ясно, что в числителе и знаменателе одинаковое количество чисел; этим и воспользуемся для того, чтобы избавиться от знаменателя. Будем называть числа в числителе "верхними числами", а в знаменателе -- "нижними". Так вот, ясно, что для каждого нижнего числа есть одно или несколько кратных ему верхних.

Если совпадение однозначно, то удается очень быстро избавиться от знаменателя; останется только числитель. Затем нужно быстро пробежаться по нему, чтобы найти модуль этого произведения. Тогда алгоритму требуется два прохода по массиву верхних чисел и $2*10^5$ делений, и $10^5$ умножений, и он укладывается в те самые 3 секунды.

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

Вот и получается, что по такому алгоритму я могу считать биномиальные коэффициенты в общем случае за 1 секунду, что в 166 раз хуже, чем надо.

Похоже, нужно как-то опираться на тот факт, что вычисления идут по модулю, но как -- не знаю.

Подскажите, пожалуйста!

На всякий случай, исходная задача -- подсчитать количество перестановок $1..N$ с неподвижными точками по модулю $p$, где $0 < N < 10^9$, она сводится к вычислению биномиального коэффициента и субфакториала.

 Профиль  
                  
 
 
Сообщение16.03.2011, 08:58 
Заслуженный участник
Аватара пользователя


14/02/07
2648
Сперва посчитать, какая степень $p$ в числителе, а какая в знаменателе. Если в числителе больше, сразу вывод.

Если одинакова - использовать теорему Вильсона (легко гуглится). Таким образом для каждого биномиального коэффициента останется меньше $2p$ умножений (деления тоже можно заменить умножениями по теореме Вильсона).

Если надо будет, напишу подробней.

Впрочем, я не понял, как эти биномиальные коэффициенты по модулю связаны с перестановками, ну и ладно.

 Профиль  
                  
 
 Re:
Сообщение16.03.2011, 09:16 


16/03/11
31
Хорхе в сообщении #423442 писал(а):
Впрочем, я не понял, как эти биномиальные коэффициенты по модулю связаны с перестановками, ну и ладно.

Сама задачка -- тут: https://www.spoj.pl/problems/SEQN/
Пусть $p_n(k)$ -- количество перестановок размера $n$ с $k$ неподвижными точками. Я смог вывести, что $p_n(k) = C_n^k \times !(n - k)$, где $!n$ -- субфакториал от $n$. Для субфакториала я смог найти эффективное вычисление по модулю, а $C_n^k$ от больших чисел приводится к меньшим по модулю. Засада в том, что даже с маленькими числами они долго считаются.

Про сравнение модулей в числителе/знаменателе и теорему Вильсона -- спасибо, сейчас поищу...

-- Ср мар 16, 2011 16:35:50 --

Хорхе в сообщении #423442 писал(а):
Сперва посчитать, какая степень $p$ в числителе, а какая в знаменателе. Если в числителе больше, сразу вывод.

Если одинакова - использовать теорему Вильсона (легко гуглится). Таким образом для каждого биномиального коэффициента останется меньше $2p$ умножений (деления тоже можно заменить умножениями по теореме Вильсона).


Хм, к сожалению, так и не понял, как это использовать. Поясните, пожалуйста.

Чтобы посчитать степень в числителе/знаменателе, мне нужно будет факторизовать числа оттуда. По-моему, это займет больше времени, чем 3 секунды. Или существуют быстрые способы факторизации?

Нашел только теорему Вильсона про простоту числа, о том, что $n$ простое, если $(n - 1)! + 1 \equiv 0\pmod{n}$. Придумал только один способ применения: разложим $p$ на простые $p_i$, после чего можно будет заменить ${p_{max}}!$ на 1 (по модулю. Сюда же как-то ещё добавить Китайскую теорему об остатках, да?), и это несколько ускорит умножение.

 Профиль  
                  
 
 
Сообщение16.03.2011, 10:36 


16/03/11
31
Понял, как можно сравнительно быстро проверить вхождения $p$ в числитель/знаменатель. Для этого нужно факторизовать одно только $p$ (на простые множители $pp_i$ со степенями $ps_i$). После чего пробежаться по числам числителя/знаменателя и посчитать степени вхождения простых чисел $pp_i$ в эти числа. Получатся степени $up_i$ для числителя и $down_i$ для знаменателя. По этим массивам получаем кратность $pp_i$ и, соответственно, степерь $p$. Так?

Но тогда получается оценка: $10^5$ чисел в числителе/знаменателе, и около 1-17 простых чисел в разложении $p$. Т.е. нужно будет около $10^6 \ldots 10^7$ делений, а это уже на грани допустимого. Т.е. такой способ не очень подходит :(

 Профиль  
                  
 
 Re: Остаток от деления больших биномиальных коэффициентов
Сообщение16.03.2011, 13:55 
Заслуженный участник


08/04/08
8560
Насколько я помню, $C_n^k \mod p$ можно найти быстро как Хорхе говорит, только там уже была готовая формула, именная, я ее забыл
Что-то вроде $\binom{n}{k} \equiv \binom{n \mod p}{k \mod p} \binom{\left[ \frac{n}{p} \right]}{\left[ \frac{f}{p} \right]} \pmod p$. И доказывалась она так же: представляем $n=pq+r, k=pu+v$, подставляем, выделяем отдельно множители с $pq+1,...,pq+r;pu+1,...,pu+v$ - это дает один биномиальный коэффициент, а во втором с помощью теоремы Вильсона некратные $p$ множители дают 1, а в кратных сокращается $p$ и получается бином. коэффициент, только меньший...

-- Ср мар 16, 2011 16:56:19 --

Вот оно! формула Люка!
http://ru.wikipedia.org/wiki/%D2%E5%EE% ... B%FE%EA%E0

 Профиль  
                  
 
 Re: Остаток от деления больших биномиальных коэффициентов
Сообщение16.03.2011, 16:15 
Заслуженный участник


04/05/09
4586
Если $p$ простое, то я бы решал так:
Заранее вычислить обратные элементы всех чисел $[1...{10}^5]$ по модулю $p$. Это можно сделать за порядка ${10}^5$ операций.
А дальше вычислять результат, пользуясь умножением на обратный элемент вместо деления.

 Профиль  
                  
 
 Re: Остаток от деления больших биномиальных коэффициентов
Сообщение16.03.2011, 16:27 


16/03/11
31
Sonic86 в сообщении #423517 писал(а):
Вот оно! формула Люка!
http://ru.wikipedia.org/wiki/%D2%E5%EE% ... B%FE%EA%E0


Эх, я уже использовал её, чтобы привести $C_n^k$$0 < k < n < 10^9$) к произведениям $C_{n'}^{k'}$$0 < k' < n' < p < 10^5$).

Остается только возможность ещё раз использовать эту формулу для множителей $p$, но если $p$ -- простое и большое, то эта формула не помогает.

 Профиль  
                  
 
 Re:
Сообщение16.03.2011, 18:25 
Заслуженный участник


04/05/09
4586
malphunction в сообщении #423451 писал(а):
Но тогда получается оценка: $10^5$ чисел в числителе/знаменателе, и около 1-17 простых чисел в разложении $p$. Т.е. нужно будет около $10^6 \ldots 10^7$ делений, а это уже на грани допустимого.
Зачем так много делений?
Количество вхождений простого делителя в факториал вычисляется по простой формуле.

 Профиль  
                  
 
 Re: Остаток от деления больших биномиальных коэффициентов
Сообщение17.03.2011, 11:06 
Заслуженный участник


08/04/08
8560
venco писал(а):
Количество вхождений простого делителя в факториал вычисляется по простой формуле.

Угу, вот по такой:
$$\text{ord}_p(n!) = \sum\limits_{k=1}^N \left[ \frac{n}{p^k}\right]$$
$N = [\log _p n]$
Сразу же можно из нее вывести и формулу для $\text{ord}_p(C_n^k)$, если нужно.

 Профиль  
                  
 
 Re: Остаток от деления больших биномиальных коэффициентов
Сообщение20.03.2011, 08:05 


16/03/11
31
Всем спасибо за ответы! В результате у меня получается следующая схема:

  • Разложить $p$ на простые множители $p_i$.
  • Для каждого $p_i$ вычислить $c_i = C_n^k \pmod{p_i}$, вычисление такого бином.коэфф. сокращается по теореме Люка к вычислению $C_{n_{j}}^{k_{j}}$, где $n_j < p_i$.
  • Всё равно, возможно, что $p_i$ будет большим (или $p$ -- простым), тогда для вычисления модуля б.к. обращаем знаменатель (используя расширенный алгоритм Евклида, например).
  • Далее, из $c_i$ сооружаем целиком $C_n^k$, можно для этого использовать метод независимых остатков из "Конкретной математики" (стр.152). Правда, тут проблема в том, что при разложении $p$ могут получаться неединичные степени...

Сложность этого алгоритма не считал, но, кажется, она будет меньше, чем $10^9$ умножений.

Как вы считаете, всё ли правильно?

 Профиль  
                  
 
 
Сообщение20.03.2011, 15:58 
Заслуженный участник


08/04/08
8560
malphunction писал(а):
Всё равно, возможно, что $p_i$ будет большим (или $p$ -- простым), тогда для вычисления модуля б.к. обращаем знаменатель (используя расширенный алгоритм Евклида, например)

Я вот не знаю, если честно, но, по-моему, очевидного преимущества алгоритма Евклида по сравнению с простым делением нету. Обращение числа $j$ по модулю $p$ потребует $O(\ln p)$ операций, для $1 \leq j \leq k$ надо $O(k \ln p)$ операций, а так всего $k$ делений :roll:
malphunction писал(а):
Правда, тут проблема в том, что при разложении $p$ могут получаться неединичные степени...

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

Еще мелочь: $C_n^k = C_n^{n-k}$ - можно добиться $k \leq \frac{n}{2}$

-- Вс мар 20, 2011 19:03:35 --

Еще. Пусть мы считаем $C^k_n = \prod\limits_{p_j \leq n} p_j^{a_j}$.
Из формулы для $a_j$ получается, что надо $\frac{\ln n}{\ln p_j}$ операций выполнить, да потом перемножить $O(\pi (n)) = O(\frac{n}{\ln n})$ чисел. Это всего $\sum\limits_{p \leq n}\frac{\ln n}{\ln p_j} + O(\frac{n}{\ln n}) = O(\ln n \frac{n}{\ln n}) + O(\frac{n}{\ln n}) = O(n)$ операций, так что этот способ, видимо медленнее, чем просто $k$ умножений и $k$ делений, хотя точно не уверен...

 Профиль  
                  
 
 Re:
Сообщение21.03.2011, 03:21 


16/03/11
31
Sonic86 в сообщении #425099 писал(а):
Я вот не знаю, если честно, но, по-моему, очевидного преимущества алгоритма Евклида по сравнению с простым делением нету.


Нужно делать вычисления в целочисленной арифметике, оставаясь в пределах машинных слов ($n \in [0..(2^{32}-1)]$ в моём случае). Я пробовал, чисто для проверки, использовать всякие библиотеки длинной арифметики, типа NTL; там можно посчитать б.к. по обычной формуле, разделив громадный числитель на знаменатель, но получается слишком медленно.

Поэтому и хочется как-то исхитриться и остаться в машинной разрядной сетке. Если заменить деление умножением, то мы избавимся от потенциальных вещественных чисел (дробей), и при умножении у нас всегда будут целые числа. Рост произведения можно контролировать, деля его на $p$. То есть в принципе, б.к. можно посчитать где-то за $O(n \log n)$ умножений, но вот коэффициент при этом $O$ будет большой :(

Ну и в самом начале топика я предлагал такой способ: пройтись по числителю/знаменателю и посокращать, сколько можно. Проход требует линейного времени и очень быстр. Но всё равно, часть этого знаменателя остается, где-то $1/50$ чисел. Можно обращать именно эти числа, а не весь знаменатель.

За остальные подсказки -- спасибо, сегодня вечером ещё поразбираюсь...

-- Пн мар 21, 2011 10:25:35 --

venco в сообщении #423580 писал(а):
Заранее вычислить обратные элементы всех чисел $[1...{10}^5]$ по модулю $p$. Это можно сделать за порядка ${10}^5$ операций.


А как можно сделать за ${10}^5$ операций? Если использовать алгоритм бинарного возведения в степень (http://e-maxx.ru/algo/reverse_element, вроде бы он быстрее сходится, чем алг.Евклида), то на каждый элемент потребуется где-то $\log 32$ умножений и делений, т.е. около 10 операций. И самих элементов $10^5$, получается всего $10^6$ операций.

Как же сделать быстрее?

 Профиль  
                  
 
 
Сообщение21.03.2011, 09:19 
Заслуженный участник


08/04/08
8560
malphunction писал(а):
там можно посчитать б.к. по обычной формуле, разделив громадный числитель на знаменатель, но получается слишком медленно.

Ааа, понял! Я туплю, я просто забыл, что $C_n^k$ очень длинный получается...
А вот это вот:
malphunction писал(а):
Ну и в самом начале топика я предлагал такой способ: пройтись по числителю/знаменателю и посокращать, сколько можно. Проход требует линейного времени и очень быстр. Но всё равно, часть этого знаменателя остается, где-то $1/50$ чисел. Можно обращать именно эти числа, а не весь знаменатель.

и аналогичное сообщение в первом посте я не понял. Как это часть знаменателя остается? Она никогда не остается, биномиальный коэффициент-то целый! Это только если Вы на каком-то шаге вычисления не берете остаток от деления - тогда сразу возникает облом.
Но если сокращать и даже игнорировать то, что тогда придется использовать длинные целые, то будет то, что я вот тут написал:
Sonic86 писал(а):
Еще. Пусть мы считаем $C^k_n = \prod\limits_{p_j \leq n} p_j^{a_j}$.... $O(n)$ операций

+ хранить множество простых за $O(\frac{n}{\ln n})$ памяти.

Кстати, я еще вспомнил треугольник Паскаля. $C_n^k = C_{n-1}^k + C_{n-1}^{k-1}$. Никакого деления и умножения, только взятие остатка. Я не знаю, как машинно считается остаток, но если и его считать сложно (через деление), можно его заменить if-ом с одним вычитанием. Только надо с ним похимичить еще наверное...

(Оффтоп)

вообще тема интересная, меня вот интересовало как считать $C_{2n}^n \mod n$, но тогда компа у меня еще не было...


-- Пн мар 21, 2011 13:08:08 --

Наверное я с треугольником Паскаля загнул. Это надо завести массив из $k$ элементов и выполнить в нем $O(kn)=O(n^2)$ операций сложения.
А если попытаться использовать тождество $C_n^k = \sum\limits_{j=0}^k C_k^jC^{k-j}_{n-k}$, то тут будет и умножение, и $O(n)$ сложений, и рекурсивное вычисление...

 Профиль  
                  
 
 Re:
Сообщение21.03.2011, 11:26 


16/03/11
31
Хорхе в сообщении #423442 писал(а):
Если одинакова - использовать теорему Вильсона (легко гуглится). Таким образом для каждого биномиального коэффициента останется меньше $2p$ умножений (деления тоже можно заменить умножениями по теореме Вильсона).


Кто-нибудь, можете объяснить, как применить теорему Вильсона? Особенно для замены деления на умножения??
Хорхе, к сожалению, не отвечает :(

 Профиль  
                  
 
 Re: Остаток от деления больших биномиальных коэффициентов
Сообщение21.03.2011, 11:58 
Заслуженный участник


08/04/08
8560
malphunction писал(а):
Кто-нибудь, можете объяснить, как применить теорему Вильсона? Особенно для замены деления на умножения??
Хорхе, к сожалению, не отвечает :(

Я могу. Сначала используете формулу Люка и тогда $0 \leq k \leq n <p$, а потом поскольку $(p-1)! \equiv -1 \pmod p$ будет $C_n^k = \frac{n(n-1)...(n-k+1)}{k!} \equiv - \frac{n(n-1)...(n-k+1)(p-1)!}{k!} = $
$= - n(n-1)...(n-k+1) \cdot (k+1)(k+2)...(p-1) \pmod p$

Для модулей $p^a,a>1$ теорема вроде уже не работает :roll:

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

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



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

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


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

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