You signed in with another tab or window. Reload to refresh your session.You signed out in another tab or window. Reload to refresh your session.You switched accounts on another tab or window. Reload to refresh your session.Dismiss alert
title = "Binomial Coefficients Modulo a Prime: Fermat's Theorem and the Non-Adjacent Selection Problem"
7
+
8
+
[header]
9
+
caption = ""
10
+
image = ""
11
+
12
+
+++
13
+
14
+
In the [previous post][efficient-impl], we implemented the closed form $F_{n,m} = \binom{n-m+1}{m}$ using Python's `math.factorial`, and with `scipy` and `sympy`. Here we cover the common competitive-programming case: computing the answer **modulo a large prime** $M$ (e.g. $M = 10^9+7$).
15
+
16
+
## Why modulo?
17
+
18
+
In counting problems, the result can be huge even for moderate input. Often the problem asks for the answer modulo a big prime so that it fits in a standard integer type. We could compute the full number and then take the remainder, but that forces expensive long-integer arithmetic. Computing **everything** modulo $M$ from the start is much faster.
19
+
20
+
## From binomials to modular inverses
21
+
22
+
We have $F_{n,m} = \binom{n-m+1}{m} = \frac{(n-m+1)!}{m!\,(n-2m+1)!}$. To compute this mod $M$, we need factorials mod $M$ and division mod $M$. Division mod $M$ is multiplication by the **modular inverse**: for prime $M$ and $0 < x < M$, the inverse of $x$ is $x^{M-2} \bmod M$ by [Fermat's little theorem](https://en.wikipedia.org/wiki/Fermat%27s_little_theorem). In Python we can use `pow(x, M - 2, M)`.
23
+
24
+
## Implementation
25
+
26
+
```python
27
+
import functools
28
+
29
+
M =10**9+7
30
+
31
+
deff_binom_mod(n, m):
32
+
assert n >=0and m >=0
33
+
34
+
if n +1<2*m:
35
+
return0
36
+
37
+
return binom_mod(n - m +1, m)
38
+
39
+
defbinom_mod(n, m):
40
+
assert0<= m <= n
41
+
42
+
return ((fact_mod(n) * inv_mod(fact_mod(m))) % M * inv_mod(fact_mod(n - m))) % M
43
+
44
+
@functools.lru_cache(maxsize=None)
45
+
deffact_mod(m):
46
+
if m <=1:
47
+
return1
48
+
49
+
return (m * fact_mod(m -1)) % M
50
+
51
+
definv_mod(x):
52
+
returnpow(x, M -2, M)
53
+
```
54
+
55
+
All operations stay in the ring of integers mod $M$. The only non-obvious part is modular division: we replace division by $d$ with multiplication by `inv_mod(d)` using Fermat's little theorem.
56
+
57
+
## Benchmarks
58
+
59
+
Compared to computing the full binomial and then taking the remainder, the modular version avoids long arithmetic and is much faster:
60
+
61
+
```python
62
+
fact_mod(10000) # for caching factorials
63
+
64
+
funcs = [f_binom_mod, f_binom, f_sci, f_sym]
65
+
66
+
test(10000, 1000, funcs)
67
+
test(10000, 2000, funcs)
68
+
test(10000, 3000, funcs)
69
+
```
70
+
71
+
Example output:
72
+
73
+
```
74
+
f(10000,1000): 450169549
75
+
f_binom_mod: 0.0000 sec, x 1.00
76
+
f_binom: 0.0073 sec, x 337.60
77
+
f_sci: 0.0011 sec, x 49.33
78
+
f_sym: 0.0076 sec, x 353.22
79
+
80
+
f(10000,2000): 75198348
81
+
f_binom_mod: 0.0000 sec, x 1.00
82
+
f_binom: 0.0063 sec, x 368.94
83
+
f_sci: 0.0026 sec, x 153.33
84
+
f_sym: 0.0053 sec, x 308.93
85
+
86
+
f(10000,3000): 679286557
87
+
f_binom_mod: 0.0000 sec, x 1.00
88
+
f_binom: 0.0060 sec, x 361.12
89
+
f_sci: 0.0056 sec, x 338.13
90
+
f_sym: 0.0053 sec, x 319.02
91
+
```
92
+
93
+
The same pattern—factorials mod $M$ plus Fermat-based inverses—works for any combinatorial formula that can be written in terms of factorials and binomials modulo a prime.
Copy file name to clipboardExpand all lines: content/post/efficient-implementation-non-adjacent-selection.md
+3-81Lines changed: 3 additions & 81 deletions
Display the source diff
Display the rich diff
Original file line number
Diff line number
Diff line change
@@ -15,7 +15,7 @@ In the [previous post][two-var-recursive], we derived the closed form for the no
15
15
16
16
$$ F_{n, m} = {n - m + 1 \choose m} $$
17
17
18
-
Now we discuss how to implement this efficiently in Python—from a simple factorial-based solution to library implementations and modular arithmetic for competitive programming.
18
+
Now we discuss how to implement this efficiently in Python—from a simple factorial-based solution to library implementations. For the common case of computing the answer **modulo a large prime** (e.g. in competitive programming), see the [next post][binom-mod].
19
19
20
20
## Fast Solutions Based on Binomials
21
21
@@ -142,86 +142,8 @@ You can play with running tests on different $n$ and $m$.
142
142
What I saw that actually there is no clear winner between the last 3 implementations.
143
143
Probably, the most of the time is spent on the long arithmetic computation.
144
144
145
-
## Modular Arithmetics
146
-
147
-
In questions where it is required to count some objects, not rarely the answer might be very big even on very small input.
148
-
In such case, typically it is asked to print the answer modulo some big prime integer, let's say, $M=1000^3+7$.
149
-
Since Python has built-in long arithmetics, we can apply modulo on the final result,
150
-
but executing the entire algorithm with long arithmetics while knowing that only small part of it is really important is very costly,
151
-
and of course, not that efficient.
152
-
153
-
Let's look, briefly, at very simple change we can do for `f_binom` function that will speed up the computation significantly:
154
-
155
-
```python
156
-
import functools
157
-
158
-
M =10**9+7
159
-
160
-
deff_binom_mod(n, m):
161
-
assert n >=0and m >=0
162
-
163
-
if n +1<2*m:
164
-
return0
165
-
166
-
return binom_mod(n - m +1, m)
167
-
168
-
defbinom_mod(n, m):
169
-
assert0<= m <= n
170
-
171
-
return ((fact_mod(n) * inv_mod(fact_mod(m))) % M * inv_mod(fact_mod(n - m))) % M
172
-
173
-
@functools.lru_cache(maxsize=None)
174
-
deffact_mod(m):
175
-
if m <=1:
176
-
return1
177
-
178
-
return (m * fact_mod(m -1)) % M
179
-
180
-
definv_mod(x):
181
-
returnpow(x, M -2, M)
182
-
```
183
-
184
-
As we can see, all the operations are computed modulo $M$.
185
-
The function `fact_mod` is recursive but uses Memoization.
186
-
The most tricky part is how to implement modular-division.
187
-
From [Fermat's little theorem](https://en.wikipedia.org/wiki/Fermat%27s_little_theorem),
188
-
we know that if $M$ is prime and $0 < x < M$, then $x^{-1} \equiv x^{M-2} \pmod M$.
189
-
This allows to compute the multiplicative inverse of $x$ using the Python's built-in function
Let's test the new approach against other implementations:
193
-
194
-
```python
195
-
fact_mod(10000) # for caching factorials
196
-
197
-
funcs = [f_binom_mod, f_binom, f_sci, f_sym]
198
-
199
-
test(10000, 1000, funcs)
200
-
test(10000, 2000, funcs)
201
-
test(10000, 3000, funcs)
202
-
```
203
-
204
-
It is not a surprise that taking the benefits of modular computations results in the huge speedup in running-time:
205
-
206
-
```
207
-
f(10000,1000): 450169549
208
-
f_binom_mod: 0.0000 sec, x 1.00
209
-
f_binom: 0.0073 sec, x 337.60
210
-
f_sci: 0.0011 sec, x 49.33
211
-
f_sym: 0.0076 sec, x 353.22
212
-
213
-
f(10000,2000): 75198348
214
-
f_binom_mod: 0.0000 sec, x 1.00
215
-
f_binom: 0.0063 sec, x 368.94
216
-
f_sci: 0.0026 sec, x 153.33
217
-
f_sym: 0.0053 sec, x 308.93
218
-
219
-
f(10000,3000): 679286557
220
-
f_binom_mod: 0.0000 sec, x 1.00
221
-
f_binom: 0.0060 sec, x 361.12
222
-
f_sci: 0.0056 sec, x 338.13
223
-
f_sym: 0.0053 sec, x 319.02
224
-
```
145
+
When the problem asks for the answer **modulo a large prime** (e.g. $10^9+7$), computing everything mod $M$ from the start is much faster than using long integers. We cover that in a [separate post][binom-mod]: binomial coefficients modulo a prime using Fermat's little theorem.
Copy file name to clipboardExpand all lines: content/post/two-var-recursive-func.md
+2-1Lines changed: 2 additions & 1 deletion
Display the source diff
Display the rich diff
Original file line number
Diff line number
Diff line change
@@ -81,9 +81,10 @@ Which actually equals to
81
81
82
82
$$ F\_{n, m} = {n - m + 1 \choose m} $$
83
83
84
-
In the [next post][efficient-impl], we discuss how to implement this closed form efficiently in Python—from a simple factorial-based solution to library implementations and modular arithmetic for competitive programming.
84
+
In the [next post][efficient-impl], we implement this closed form in Python (factorial-based and with scipy/sympy). For computing the answer modulo a large prime, see [Binomial Coefficients Modulo a Prime][binom-mod].
<p>In the <ahref="/post/efficient-implementation-non-adjacent-selection/">previous post</a>, we implemented the closed form $F_{n,m} = \binom{n-m+1}{m}$ using Python’s <code>math.factorial</code>, and with <code>scipy</code> and <code>sympy</code>. Here we cover the common competitive-programming case: computing the answer <strong>modulo a large prime</strong> $M$ (e.g. $M = 10^9+7$).</p>
985
+
<h2id="why-modulo">Why modulo?</h2>
986
+
<p>In counting problems, the result can be huge even for moderate input. Often the problem asks for the answer modulo a big prime so that it fits in a standard integer type. We could compute the full number and then take the remainder, but that forces expensive long-integer arithmetic. Computing <strong>everything</strong> modulo $M$ from the start is much faster.</p>
<p>In the <ahref="/post/two-var-recursive-func/">previous post</a>, we derived the closed form for the non-adjacent selection problem:</p>
985
1074
<p>$$ F_{n, m} = {n - m + 1 \choose m} $$</p>
986
-
<p>Now we discuss how to implement this efficiently in Python—from a simple factorial-based solution to library implementations and modular arithmetic for competitive programming.</p>
1075
+
<p>Now we discuss how to implement this efficiently in Python—from a simple factorial-based solution to library implementations. For the common case of computing the answer <strong>modulo a large prime</strong> (e.g. in competitive programming), see the <ahref="/post/binomial-modulo-prime/">next post</a>.</p>
987
1076
<h2id="fast-solutions-based-on-binomials">Fast Solutions Based on Binomials</h2>
988
1077
<p>We can reflect the closed form in very trivial Python code:</p>
989
-
<divclass="highlight"><pretabindex="0" style="color:#f8f8f2;background-color:#272822;-moz-tab-size:4;-o-tab-size:4;tab-size:4;"><codeclass="language-Python" data-lang="Python"><spanstyle="display:flex;"><span><spanstyle="color:#f92672">import</span> math
</span></span><spanstyle="display:flex;"><span><spanstyle="color:#66d9ef">assert</span> n <spanstyle="color:#f92672">>=</span><spanstyle="color:#ae81ff">0</span><spanstyle="color:#f92672">and</span> m <spanstyle="color:#f92672">>=</span><spanstyle="color:#ae81ff">0</span>
993
-
</span></span><spanstyle="display:flex;"><span>
994
-
</span></span><spanstyle="display:flex;"><span><spanstyle="color:#66d9ef">if</span> n <spanstyle="color:#f92672">+</span><spanstyle="color:#ae81ff">1</span><spanstyle="color:#f92672"><</span><spanstyle="color:#ae81ff">2</span><spanstyle="color:#f92672">*</span>m:
</span></span><spanstyle="display:flex;"><span><spanstyle="color:#66d9ef">return</span> binom(n <spanstyle="color:#f92672">-</span> m <spanstyle="color:#f92672">+</span><spanstyle="color:#ae81ff">1</span>, m)
</span></span><spanstyle="display:flex;"><span><spanstyle="color:#66d9ef">assert</span><spanstyle="color:#ae81ff">0</span><spanstyle="color:#f92672"><=</span> m <spanstyle="color:#f92672"><=</span> n
1001
-
</span></span><spanstyle="display:flex;"><span>
1002
-
</span></span><spanstyle="display:flex;"><span><spanstyle="color:#66d9ef">return</span> math<spanstyle="color:#f92672">.</span>factorial(n) <spanstyle="color:#f92672">//</span> math<spanstyle="color:#f92672">.</span>factorial(m) <spanstyle="color:#f92672">//</span> math<spanstyle="color:#f92672">.</span>factorial(n <spanstyle="color:#f92672">-</span> m)
1003
-
</span></span></code></pre></div><p>This implementation overperforms significantly the initial DP and memoization solutions from <ahref="/post/intro-to-dp/">Introduction to Dynamic Programming and Memoization</a>.
1004
-
A naive implementation of <code>math.factorial()</code> might make $n$ multiplications.
1005
-
This could still be faster than doing $\Theta(n)$ additions in DP approach.</p>
1006
1078
</div>
1007
1079
</a>
1008
1080
@@ -1035,7 +1107,7 @@ <h2 id="fast-solutions-based-on-binomials">Fast Solutions Based on Binomials</h2
0 commit comments