### PROBLEM LINK:

**Author:** Anton Lunyov

**Tester:** Hiroto Sekido

**Editorialist:** Anton Lunyov

### DIFFICULTY:

HARD

### PREREQUISITES:

Binomial coefficients modulo prime powers, 128-bit arithmetic

### PROBLEM:

For the given **n** and **R** such that **0 ≤ R < 2 ^{n}** you need to find the smallest

**K**such that

**C(2**,

^{n}− 1, K) mod 2^{n}= Rwhere

**C(N, K) = N! / K! / (N − K)!**, or report that such

**K**does not exist. Note, that

**n**could be up to

**120**.

### EXPLANATION:

At first note that the problem requires to code 128-bit arithmetic. This can be made by representing each number as a pair of 64-bit values (lowest 64 bits and highest 64 bits) and then some routine job should be made to code all necessary operations. The only non-trivial one is multiplication. It can be made using only six 64-bit multiplications (refer to the tester’s solution). Actually only five is enough (refer to the author’s solution).

The first thing to notice is:

**C(2 ^{n} − 1, K) = fact2(2^{n} − 1) / fact2(K) / fact2(2^{n} − 1 − K) * C(2^{n−1} − 1, K div 2)**, (1)

where

**fact2(N)**is the product of all odd numbers up to

**N**. (If some one want to see the proof ask me in comment I and will provide it when I will have a free time.)

It immediately follows from here that **C(2 ^{n} − 1, K)** is always odd. So if

**R**is even then the answer is

**−1**.

We will prove that for each odd

**R**in the range

**{0, 1, 2, …, 2**,

^{n}− 1}there exists unique

**K**in the range

**{0, 1, 2, …, 2**such that

^{n−1}− 1}**C(2**.

^{n}− 1, K) mod 2^{n}= RAt first we infer the following congruence from (1):

**C(2 ^{n} − 1, 2 * K + X) ≡ (−1)^{K + X} * C(2^{n−1} − 1, K) (mod 2^{n})**, where

**X**is

**0**or

**1**.

(Again if some one want to see the proof ask me in comment I and will provide it when I will have a free time.)

Also we will need the following property:

**C(2 ^{n} − 1, 2 * K) + C(2^{n} − 1, 2 * K + 1) ≡ 2^{n} (mod 2^{n+1})**, for

**n ≥ 2**and

**0 ≤ S < 2**.

^{n−1}Then after some non-trivial analysis we get the following simple iterative routine that will find the answer:

```
INVBINCF(n, R)
K = (R mod 4) div 2
for j = 3 to n do
par = K mod 2
if C(2^(j-1) - 1, K) mod 2^j != R mod 2^j
K = K xor 1
K = 2 * K + par
return K
```

The proof in Russian is contained here (pages 57-58) and contains a lot of formulas. If you really want to see the proof in English, ask me in comment and maybe some day I will translate it

Thus, the problem boils down to effective calculation of **C(2 ^{n} − 1, K)** modulo powers of two.

And here we ask for help of British mathematician Andrew Granville, the number theory expert

We refer to Theorem 2 of his article about binomial coefficients modulo prime powers (see end of the page 3).

For

**p = 2**it sounds as follows:

**fact2(2 * u) ≡ sgn * Product[ fact2(2 * j)**, for

^{b(r, j, u)}: 1 ≤ j ≤ r] (mod 2^{2 * r + 1})**r ≥ 2**, (2)

where

**b(r, j, u) = u / j * Product[ (u**,

^{2}− i^{2}) / (j^{2}− i^{2}) : 1 ≤ i ≤ r, i ≠ j]and

**sgn**is

**1**or

**−1**and can be found by comparing residues of LHS and RHS modulo

**4**.

Clearly

**b(r, j, u)**could be calculated modulo

**2**in view of the property:

^{2 * r − 1}**a**for

^{2n − 2}= 1 (mod 2^{n})**n ≥ 3**and odd

**a**.

The proof of (2) is based on 2-adic exponent and logarithm and is like magic for me.

But straightforward implementation of **C(2 ^{n} − 1, K)** using formulas (1) and (2), together with

`INVBINCF`

routine above leads to quite slow solution.Note, that actually we can maintain value **C(2 ^{j} − 1, K) mod 2^{n}** in the loop of

`INVBINCF`

and recalculate it using formula (1). Namely, if **K = K xor 1**is performed we actually replace

**K**by adjacent integer and can update

**C(2**straightforwardly using one multiplication and one division.

^{j}− 1, K)While when we do

**K = 2 * K + par**we should transform

**C(2**to

^{j}− 1, K)**C(2**and this exactly can be made by formula (1) by calculating three values of

^{j + 1}− 1, 2 * K + par)**fact2(X)**by formula (2). Also we can maintain numerator and denominator of

**C(2**separately and get rid of divisions by rewriting check in

^{j}− 1, K)`INVBINCF`

in obvious way. Refer to routine `invbin`

in author’s solution for clear implementation of this.So `INVBINCF`

is ready to be used but we need to figure it out how to calculate **fact2(X) mod 2 ^{n}** efficiently.

At first we should calculate **b(j, r, u) mod 2 ^{n−2}** efficiently. Here main complication is that denominator could have even factors and we can’t directly apply technique of inverse residues. To handle this we will use extended representation of integers in the form

**2**, where

^{a}* b**b**is odd. So we could create a data structure that deals with this representation. Involving negative

**a**and inverse residue technique for odd part of the number, we can handle non-integer numbers as well. Thus, to calculate

**b(r, j, u)**we use the following plan:

- Before processing test data we precalculate values
**pseudo_fact[r][j] = 1 / j * Product[ 1 / (j**in extended format. They equal to inverse of denominators of^{2}− i^{2}) : 1 ≤ i ≤ r, i ≠ j]**b(r, j, u)**. - Inside
`fact2`

routine we calculate values**pref[k] = u * Product[ (u**and^{2}− i^{2}) : 1 ≤ i ≤ k]**suf[k] = Product[ (u**in two simple^{2}− i^{2}) : k ≤ i ≤ r]**O®**loops (again in extended format). - Then
**b(r, j, u) = pref[j - 1] * suf[j + 1] * pseudo_fact[r][j]**and after this, we can transform it to usual integer representation.

Thus, the complexity of precalc is **O(N ^{3})** (we calculate

**O(N**values

^{2})**pseudo_fact[r][j]**and each step requires finding inverse residue, which can be made in

**O(N)**time) and calculation of all values

**b(r, j, u)**for the given

**u**is only

**O(N)**. Actually, saving inverses for numbers up to

**N**in

**O(N**time we can speed up the precalc part to

^{2})**O(N**, but even

^{2})**O(N**precalc is many times faster then some tricky precalc below.

^{3})Now having values **b(r, j, u)** calculated we can find **fact2(u)** using **r** `PowerMod`

operations, where `PowerMod[a, b, n]`

returns **a ^{b} mod 2^{n}**. But even using exponentiation by squaring this will lead to

**O(N**complexity for the

^{2})`fact2`

routine and thus to **O(N**complexity for the

^{3})`INVBINCF`

routine, which is still too slow.To speed up this part we can do the following. Note that we need to calculate powers of small fact2-values in `fact2`

routine. So we can do some clever precalc for them. Namely, if we want to set up fast calculation of powers of some number **A** we can precalc the following its powers in advance: **powA[k][j] = A ^{j * 2H * k}** for all

**0 ≤ k < N/H**and

**0 ≤ j < 2**, where

^{H}**H**is some constant (divisor of

**N = 120**will be most convenient). Then, when we want to calculate

**A**for some

^{B}**B < 2**, we represent

^{N}**B**in

**2**-ary system as

^{H}**B = B**

_{0}+ B_{1}* 2^{H}+ … + B_{K}* 2^{H * K}and calculate

**A**as the product of known powers

^{B}**powA[0][B<sub>0</sub>] * powA[1][B<sub>1</sub>] * … * powA[K][B<sub>K</sub>]**in

**O(N / H)**time.

Thus precalculating such powers for each **fact2(2 * j)** (**1 ≤ j ≤ 60**) with constant **H = 15**, will speed up our solution in **H = 15** times to **O(N * N * N / H)**. But we can do even better. We can find prime factorization of **fact2(u)** and set up powers precalc only for odd primes up to **120** (there are only **29** of them). Then after calculation of **b(r, j, u)** we can do some simple **O(N)** loop to calculate degrees in factorization of **fact2(u)** and then calculate the answer as the product of prime powers in at most **29 * 8** multiplications. This trick also speeds up powers precalc in two times.

Overall complexity of the solution assuming constant time for each 128-bit operation is

**O(N ^{3} + π(N) * 2^{H} * N/H + T * N * (N + π(N) * N/H))**,

where

**N = 120**,

**H = 15**and

**π(N)**is the number of primes

**≤ N**. Here:

- term
**N**corresponds to precalc of^{3}**pseudo_fact[r][j]**; - term
**π(N) * 2**corresponds to precalc of prime powers;^{H}* N/H - term
**T * N * (N + π(N) * N/H)**means that for each of**T**tests we calculate**O(N)**values**fatc2(u)**in**O(N + π(N) * N/H)**time each.

You may also ask “Why did I annoy all by crappy 128-bit arithmetic?”. This is because I almost invented more elementary solution with **O(2 ^{N/3})** precalc and

**poly(N)**routine for each test. So I decided that staying under 64-bit arithmetic could be dangerous. Actually I spent the whole day to code the working 128-bit version having working 64-bit version. So this annoyed me too

### ALTERNATIVE SOLUTION:

@winger solution is very interesting and differs from described one.

It would be nice if he will describe his solution.

Solution by @mugurelionut is based on more elementary math than Granville formula. You can read about it here.

### AUTHOR’S AND TESTER’S SOLUTIONS:

Author’s solution can be found here.

Tester’s solution can be found here.

### RELATED PROBLEMS:

e-olimp - 322 - Binomial Coefficients 5

The current problem is extreme version of this one.

I’ve set this problem (for **n ≤ 40**) at Kharkiv Winter Programming School 2009 and editorial in Russian is available here (see page 57). So those who knew about this problem could have some advantage against others, but it seems that all who solved the problem was unaware of this link

The good problems to practice your skills on direct problems on calculating reduced factorials modulo prime powers:

Mathalon - 141 - Factorial trailing digits2

Mathalon - 144 - Binomial Coefficient 1 (this one is mine :P)