11

How to calculate binomial coefficient modulo 142857 for large n and r. Is there anything special about the 142857? If the question is modulo p where p is prime then we can use Lucas theorem but what should be done for 142857.

Jarod42
  • 203,559
  • 14
  • 181
  • 302
user1505986
  • 257
  • 4
  • 10
  • 1
    Note that the factorisation helps because you can use CRT and compute the coefficient modulo 11, 13, 27 and 37. – John Dvorak Oct 28 '12 at 05:40
  • 2
    Wikipedia links to a PDF about [Binary coefficients modulo prime powers](http://www.dms.umontreal.ca/~andrew/PDF/BinCoeff.pdf) – John Dvorak Oct 28 '12 at 05:48
  • I implemented that prime powers algo and it is giving correct answer when power=1 but when power!=1 it is giving wrong answer for some inputs and correct for some. something wrong with my code i guess. – user1505986 Oct 28 '12 at 10:29
  • For anyone who's interested: [here](https://github.com/niklasb/contest-algos/blob/master/number_theory.cpp#L162) is my implementation of a `binomod` function which was accepted during the contest. I used it together with the `crt_coprime` function in the same file. – Niklas B. Nov 09 '12 at 15:19

3 Answers3

14

You can actually calculate C(n,k) % m in O(n) time for arbitrary m.

The trick is to calculate n!, k! and (n-k)! as prime-power-vectors, subtract the two later from the first, and multiply the remainder modulo m. For C(10, 4) we do this:

10! = 2^8 * 3^4 * 5^2 * 7^1
 4! = 2^3 * 3^1
 6! = 2^4 * 3^2 * 5^1

Hence

C(10,4) = 2^1 * 3^1 * 5^1 * 7^1

We can calculate this easily mod m as there are no divisions. The trick is to calculate the decomposition of n! and friends in linear time. If we precompute the primes up to n, we can do this efficiently as follows: Clearly for each even number in the product 1*2*...*9*10 we get a factor of 2. For each fourth number, we get a second on and so forth. Hence the number of 2 factors in n! is n/2 + n/4 + n/8 + ... (where / is flooring). We do the same for the remaining primes, and because there are O(n/logn) primes less than n, and we do O(logn) work for each, the decomposition is linear.

In practice I would code it more implicit as follows:

func Binom(n, k, mod int) int {
    coef := 1
    sieve := make([]bool, n+1)
    for p := 2; p <= n; p++ {
        // If p is not sieved yet, it is a prime number
        if !sieve[p] {
            // Sieve of Eratosthenes
            for i := p*p; i <= n; i += p {
                sieve[i] = true
            }
            // Calculate influence of p on coef
            for pow := p; pow <= n; pow *= p {
                cnt := n/pow - k/pow - (n-k)/pow
                for j := 0; j < cnt; j++ {
                    coef *= p
                    coef %= mod
                }
            }
        }
    }
    return coef
}

This includes a Sieve of Eratosthenes, so the running time is nloglogn rather than n if the primes had been precalculated or sieved with a faster sieve.

Thomas Ahle
  • 30,774
  • 21
  • 92
  • 114
13

The algorithm is:

  • factorise the base into prime powers; 142857 = 3^3×11×13×37
  • compute the result modulo each prime power
  • combine the results using the Chinese Remainder Theorem.

To compute (n above k) mod p^q:

Source: http://www.dms.umontreal.ca/~andrew/PDF/BinCoeff.pdf, theorem 1

define (n!)_p as the product of numbers 1..n that are not divible by p

define n_j as n after erasing j least significant digits in base p

define r as n-k

define e_j as the number of carries when adding k+r, not counting the carries from j lowest digits, computing in base p

define s as 1 if p=2 & q>=3 and -1 otherwise

then (n above k) mod p^q := p^e_0 * s^e_(q-1) * concatenate(j=d..0)( (n_j!)_p / ((k_j!)_p*(r_j!)_p) ) with each term of the concatenation computing one base-p digit of the result, lowest j computing the least significant non-zero digits.

John Dvorak
  • 26,799
  • 13
  • 69
  • 83
  • I wonder how fast this is. Have you tried coding it? I'd like to compare it with my answer factorizing `n`. – Thomas Ahle Jul 01 '14 at 09:44
  • @ThomasAhle I haven't tried implementing it, but if you're going to compute factorial of a number larger than a few hundred, you are going to wait a while. – John Dvorak Jul 01 '14 at 10:02
  • The binomial is `product [(k+1)..n]/product [1..(n-k)]`, so you still need division, which you can't do modulo an arbitrary number. That's the whole point of this question, isn't it? – Thomas Ahle Jul 01 '14 at 10:07
3

What's special about 142857 is that 7 * 142857 = 999999 = 10^6 - 1. This is a factor that arises from Fermat's Little Theorem with a=10 and p=7, yielding the modular equivalence 10^7 == 10 (mod 7). That means you can work modulo 999999 for the most part and reduce to the final modulus by dividing by 7 at the end. The advantage of this is that modular division is very efficient in representation bases of the form 10^k for k=1,2,3,6. All you do in such cases is add together digit groups; this is a generalization of casting out nines.

This optimization only really makes sense if you have hardware base-10 multiplication. Which is really to say that it works well if you have to do this with paper and pencil. Since this problem recently appeared on an online contest, I imagine that's exactly the origin of the question.

eh9
  • 7,340
  • 20
  • 43
  • 2
    999999 is not really a more preferrable modulus to work with than 142857, so I don't see how this makes the problem easier to solve... You'll need Granville anyway – Niklas B. Nov 09 '12 at 14:48
  • 1
    @NiklasB. You can do the calculations by hand easier. I don't think I ever said this made a software implementation better. – eh9 Nov 09 '12 at 15:36
  • But the question was how the value could be computed? And in that context, there's nothing special about 142857 at all – Niklas B. Nov 09 '12 at 15:58
  • 1
    This question was asked: "Is there anything special about the 142857?" That's the question I answered. – eh9 Nov 09 '12 at 16:18