Остаток от деления больших биномиальных коэффициентов

[malphunction]

Здравствуйте!

Помогите разобраться с алгоритмом.

Пытаюсь решить олимпиадную задачку по программированию, вот тут обсуждение разных решений ;
задачу удалось свести к вычислению ${n \choose k} \pmod {p^r}$, где $p$ -- простое число.

Нашел по ссылке с википедии публикацию:Andrew Granville, 1997, «Arithmetic Properties of Binomial Coefficients I: Binomial coefficients modulo prime powers». Canadian Mathematical Society Conference Proceedings 20: 253-275. MR1483922.. В публикации рассматриваются два алгоритма, один -- "старый" (Davis and Webb), и усовершенствованный новый. Судя по асимптотическим оценкам, для моих целей вполне подходит "старый" алгоритм из той статьи. Соответствующая теорема приведена на 2-ой странице этой PDF-ки.

Места под теорему нет, см. следующий ответ (или pdf-ку).

Для того, чтобы разобраться с этим алгоритмом, я попытался вычислить ${70 \choose 50} \pmod {3^3}$ (равно $9$, кстати), и тут же столкнулся с неправильным результатом. Вот что получилось:

Разложим $n, k, r$ по степеням $p$:

$$n = 70 = 2\cdot3^3 + 1\cdot3^2 + 2\cdot3^1 + 1\cdot3^0, (n_3 = 2, n_2 = 1, n_1 = 2, n_0 = 1)$$

$$k = 50 = 1\cdot3^3 + 2\cdot3^2 + 1\cdot3^1 + 2\cdot3^0, (k_3 = 1, k_2 = 2, k_1 = 1, k_0 = 2)$$

$$r = 20 = 0\cdot3^3 + 2\cdot3^2 + 0\cdot3^1 + 2\cdot3^0, (r_3 = 1, r_2 = 2, r_1 = 1, r_0 = 2)$$

Вычислим $N, M, R$:

$$N_0 = 16, M_0 = 23, R_0 = 20$$

$$N_1 = 23, M_1 = 16, R_1 = 6$$

$$N_2 = 9, M_2 = 5, R_2 = 2$$

$$N_3 = 2, M_3 = 1, R_3 = 0$$

Теперь посчитаем ($(n!)_p$ обозначает произведение всех $1..n$ без степеней $p$, т.е. $(n!)_p = n! / {[k/p]! p^{[k/p]}}$):

$$ {(N_0!)_p \over { (M_0!)_p (R_0!)_p }} \cdot {(N_1!)_p \over { (M_1!)_p (R_1!)_p }} \cdot {(N_2!)_p \over { (M_2!)_p (R_2!)_p }} \cdot {(N_3!)_p \over { (M_3!)_p (R_3!)_p }} $$

Подставляем значения:

$$ {(16!)_3 \over { (23!)_3 (20!)_3 }} \cdot {(23!)_3 \over { (16!)_3 (6!)_3}} \cdot {(9!)_3 \over { (5!)_3 (2!)_3 }} \cdot {(2!)_3 \over { (1!)_3 (0!)_3 }} $$

Сокращаем:

$$ {(9!)_3 \over { (20!)_3 (6!)_p (5!)_3 } $$

Дальше считать не имеет смысла: значение в знаменателе значительно больше числителя, получается дробное число. А должно получиться $9$.

Помогите разобраться с этой теоремой и алгоритмом. Что я не так делаю?? В самой pdf-ке есть ссылка на Devis and Webb с доказательством, но я найти исходную статью не смог.

-- Вс апр 24, 2011 14:07:05 --

Вот теорема:
Дана степень простого числа $p^q$ и числа $n = m + r$.

Разложим $n$ по степеням $p$: $n = n_0 + n_1p + n_2p^2 + \ldots + n_dp^d$, и пусть $N_j$ --- наименьший положительный остаток (least positive residue) от $[n/p^j] \pmod {p^q}$. Аналогично введем $m_j, M_j, r_j, R_j$. Пусть $e_j$ --- количество индексов $i \ge j$, для которых $n_j < m_j$.

Тогда:

$$ {1 \over {p^{e_0}}} {n \choose m} \equiv (\pm 1)^{e_{q-1}} {(N_0!)_p \over { (M_0!)_p (R_0!)_p }} {(N_1!)_p \over { (M_1!)_p (R_1!)_p }} \cdots {(N_d!)_p \over { (M_d!)_p (R_d!)_p }} \pmod {p^q} $$

где $(\pm 1)$ равно $(-1)$, кроме случая $p = 2, q \ge 3$.

-- Вс апр 24, 2011 14:35:40 --

Попробовал зайти с другой стороны, заменим знаменатель на числитель (найдем обратные элементы в кольце по модулю $p^k$):

${1 \over {(20!)_3}} \pmod {27} \equiv 20$
${1 \over {(6!)_3}} \pmod {27} \equiv 4$
${1 \over {(5!)_3}} \pmod {27} \equiv 4$
$(9!)_3 = 2240
$2240 * 20 * 4 * 4 = 716800 \equiv 4 \pmod {27}


Но всё равно не то, получается $4$. А должно быть $9$.

Не понимаю, где ошибка...

[malphunction]

malphunction писал(а):

Попробовал зайти с другой стороны, заменим знаменатель на числитель (найдем обратные элементы в кольце по модулю $p^k$):

${1 \over {(20!)_3}} \pmod {27} \equiv 20$
${1 \over {(6!)_3}} \pmod {27} \equiv 4$
${1 \over {(5!)_3}} \pmod {27} \equiv 4$
$(9!)_3 = 2240
$2240 * 20 * 4 * 4 = 716800 \equiv 4 \pmod {27}

Не учел, что $27$ -- это не простое число, поэтому алгоритм бинарного возведения в степень не подходим. Так что воспользуемся алгоритмом Евклида, получим:


${1 \over {(20!)_3}} \pmod {27} \equiv -13$
${1 \over {(6!)_3}} \pmod {27} \equiv -2$
${1 \over {(5!)_3}} \pmod {27} \equiv -2$
$(9!)_3 = 2240
$2240 * -13 * -2 * -2 = -116480 \equiv 25 \pmod {27}

Какая-то фигня, в общем :(

[Sonic86]

Спасибо за статью! Очень интересно!

malphunction писал(а):

Для того, чтобы разобраться с этим алгоритмом, я попытался вычислить ${70 \choose 50} \pmod {3^3}$ (равно $9$, кстати), и тут же столкнулся с неправильным результатом.
...
Вычислим $N, M, R$:
$N_0 = 16, M_0 = 23, R_0 = 20$
$N_1 = 23, M_1 = 16, R_1 = 6$

Вы обратите внимание, что $N_j = \left[ \frac{n}{p^j}\right]$ и для остальных величин так же, а значит $N_j$ и остальные должны убывать экспоненциально, а у Вас $N_0<N_1$. И даже $N_0=n=70, M_0=m=50$. Ну и ошибка ползет дальше...

-- Вс апр 24, 2011 11:46:40 --

(Оффтоп)

Вообще как бы нехорошо, старая тема тут:
topic43247-15.html

[malphunction]

Sonic86 в сообщении #438191 писал(а):

Вы обратите внимание, что $N_j = \left[ \frac{n}{p^j}\right]$ и для остальных величин так же, а значит $N_j$ и остальные должны убывать экспоненциально, а у Вас $N_0<N_1$. И даже $N_0=n=70, M_0=m=50$.

Так ведь, насколько я понял из статьи, $N_j = \left[ \frac{n}{p^j}\right] \pmod{p^q}$. Поэтому и получается такие числа: $16, 23, 7, 2$ (для $N_j$), а не $70, 23, 7, 2$.

(Оффтоп)

Sonic86 в сообщении #438191 писал(а):

Вообще как бы нехорошо, старая тема тут:
topic43247-15.html

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

Кстати, спасибо за участие :)

-- Вс апр 24, 2011 16:46:58 --

Кстати, обратите внимание на множитель в формулировке теоремы: $1 \over p^{e_0}$; его можно перенести в правую часть равенства.

Теперь смотрите: у меня в вычислении получилось $25$, а ${25 \cdot 3^2} \pmod {3^3} \equiv 9$, т.е. тому, что нужно. Но тогда получается, что $e_0 = 2$. Но, как сказано в формулировке теоремы, $e_j$ есть количество индексов $i \ge j$, для которых $n_i < m_i$ ("Let $e_j$ be the number of indicies $i \ge j$ for which $n_i < m_i$"), и тогда $e_0 = 0$ (разложение $n = 2, 1, 2, 1$, а $m = 1, 2, 1, 2$, в порядке убывания степени). Или я неправильно перевел?

[Sonic86]

Блин, надо было сразу посмотреть, что $\gcd (9;27)=9$. У Вас значит все остальные вычисления, кроме вычисления $e_0$, верные.

malphunction писал(а):

Но тогда получается, что $e_0 = 2$. Но, как сказано в формулировке теоремы, $e_j$ есть количество индексов $i \ge j$, для которых $n_i < m_i$ ("Let $e_j$ be the number of indicies $i \ge j$ for which $n_i < m_i$"), и тогда $e_0 = 0$ (разложение $n = 2, 1, 2, 1$, а $m = 1, 2, 1, 2$, в порядке убывания степени). Или я неправильно перевел?

Перевели правильно, но действительно $e_0=2$, т.к. $n_1<m_1$ и $n_3<m_3$. Как Вы считали?
И еще, конечно же, $e_0 = \text{ord}_p(\binom{n}{m})$ - степень $p$ в биномиальном коэффициенте.

[malphunction]

А-а-а, имеется ввиду общее количество, а я, блин, минимальный индекс смотрю. Да, тогда всё сходится!!

Завтра утром попробую запрограммировать, сообщу о результатах. В теории, должно быть очень быстро; можно уложиться в нужные временные рамки...

[malphunction]

Вот тут у меня ошибка, но не критичная, как выяснилось ;)

malphunction писал(а):

Вычислим $N, M, R$:
$N_2 = 9, M_2 = 5, R_2 = 2$

Должно быть $N_2 = 7$ (т.е. $[70 / {3^2}] \pmod {27}$).

Ну и все последующие вычисления неверны. Соответственно, после сокращения в числителе остается $7$, в знаменателе --- $5, 20, 6$. Вычисляем знаменатель ($7416253644800000$), обращаем его, получается: $2$ (проверяем: $(7416253644800000 \cdot 2) \pmod {27} = 1$). Перемножаем по модулю $27$ числитель и обращенный знаменатель, получаем $20$, умножаем на $p^{e_0} = 3^2$, получаем $18$ (по модулю $27$).

Далее, $e_{q-1} = 1$ (т.к. $n_3 > m_3$), так что $(-1)^{e_{q-1}} = -1$, так что ${{70} \choose {50}} \equiv {-18} \equiv 9 \pmod{27}$.

Продолжаю проверять программу...

[malphunction]

Не получается :(

Набросал программку для проверки, она нашла проблемы.

Посчитаем ${51 \choose 37} \pmod {3^3}$:

$n_i: [0, 2, 2, 1]$
$m_i: [1, 0, 1, 1]$
$r_i: [2, 1, 1, 0]$
$N_i: [24, 17, 5, 1]$
$M_i: [10, 12, 4, 1]$
$R_i: [14, 4, 1, 0]$

Произведение в знаменателе (p-факториал от $[10, 12, 4, 14, 4]$): $15840934100992000000$, обратное к нему по модулю $27$ равно $10$.

Числитель (p-факториал от $[24, 17, 5]$ ): $1144342679718285279232000000$.

Общее произведение по модулю: $16$; умножаем на $p^{e_0} = 3^2$, получается $9$ (по модулю $27$). А Wofram Alfa говорит, что C(51, 37) mod 27 = 21.

Что я делаю не так?

[Sonic86]

Ну раз $\gcd(27;21)=3$, а $\gcd(27;9)=3^2$, то неверно вычислена именно степень

[malphunction]

Да, опять я с индексами напутал...

А при чем здесь $gcd(27;21)$? Как это относится к задаче?

[Sonic86]

malphunction писал(а):

А при чем здесь $gcd(27;21)$? Как это относится к задаче?

Когда Вы считаете остаток от деления на $p^k$, Вы по формуле отдельно выделяете произведение обратимых элементов, состоящее из $M_j,N_j,R_j$ и отдельно степень числа $p$, которая не взаимно проста с $p^k$. Любой вычет из $\mathbb{Z}_{p^k}$ имеет вид $p^a \cdot r, \gcd (r,p)=1, r \in \mathbb{Z}_{p^{k-a}}^*$ (представление единственно). И значит, вместо сравнений полученных вычетов, я могу сравнивать пары $(a;r)$:
$9=3^2 \cdot 1$
$21=3 \cdot 7$
И проверив отдельно степени, я уже тогда могу сказать, что неправильно были вычислены степени (если бы они совпадали, я бы сказал, что неверно вычислено произведение из $M_j,N_j,R_j$).

-- Пн апр 25, 2011 11:22:55 --

Иными словами, я использую нечто вроде представления группы $<\mathbb{Z}_p, \cdot>$ в виде группы $\mathbb{Z}_p^*$ (только вот тут неточность будет из-за неединственности представления) и полугруппы $<\{ 1;p;...;p^k\},\cdot>$.
Типа того

[malphunction]

Всё, заборол этот алгоритм. Самая печаль была с $e_j$, метод "посчитать количество $n_i < m_i$" из публикации не проходит. Нужно считать количество переносов при сложении цифр $n_i$ и $m_i$: если есть перенос, увеличиваем $e_j$ (и переносим перенос на следующие цифры). $e_{q-1}$ тоже нужно считать аккуратно: нужно учесть переносы из цифр $[0..{q-2}]$.

-- Вт апр 26, 2011 14:45:53 --

Sonic86 писал(а):

Когда Вы считаете остаток от деления на $p^k$...

Спасибо за разъяснения; слаб я в алгебре, и понял лишь весьма поверхностно, но эти проверки помогли разобраться с алгоритмом :)

[Sonic86]

Поздравляю!
Да пожалуйста там кстати ничего сложного нет. Просто попрактиковаться в вычислении немного...

[maxal]

Вот здесь моя реализация вычисления биномиальных коэффициентов по модулю на PARI/GP:
http://www.cse.sc.edu/~maxal/gpscripts/