231 The prime factorisation of binomial coefficients

```import Number (primesUpTo)
count :: (Integral a) => a -> a -> a
count n p= sum.takeWhile (>).tail.iterate (`div` p) \$ n
p231 m n = sum.map f.primesUpTo \$ m
where f p = p * ( count m p - count n p - count (m - n) p )
main :: IO ()
main = print \$ p231 (2*10^7) (15*10^6)
```

C++とかで実装したら速いんだろうな．

#### 追記

C++で実装してみた．

```#include <iostream>
#include <vector>
#include <cmath>
#define FOR(x, xs) for (__typeof( (xs).begin() ) x = (xs).begin(); x != (xs).end(); x++)
using namespace std;
static const int m = 20000000;
static const int n = 15000000;
static const int u = sqrt(m);
int count (int q, int p) {
int c = ;
for ( ; q > ; c += (q /= p) );
return c;
}
int main() {
vector<bool> isPrime( m + 1);
FOR (q, isPrime) *q = true;
isPrime[] = false;
isPrime[1] = false;
for (int q = 2; q <= u; q++)
if ( isPrime[q] )
for (int s = q * q; s <= m; s += q)
isPrime[s] = false;
long long s = ;
for (int q = 2; q <= m; q++)
if ( isPrime[q] ) s += q * ( count(m , q) - count(n, q) - count(m - n, q));
cout << s <<endl;
}
```