codeforces 711E. ZS and The Birthday Paradox 概率

已知一年365天找23个人有2个人在同一天生日的概率 > 50%

给出n,k ，表示现在一年有2^n天，找k个人，有2个人在同一天生日的概率，求出来的概率是a/b形式，化到最简形式，由于a,b可能非常大，对a,b分别%(10^6+3)

注意，这道题是先化到最简，再分别取模

首先，特判 k > 2^n 时，a = 1,b = 1

没有2个人同一天生日的概率为：

2^n * (2^n - 1) * ... * (2^n - k + 1) / 2^(nk)

所以a,b化简之前分别是：

a = 2^nk-n - (2ⁿ - 1) * (2ⁿ - 2) * ... * (2ⁿ - k + 1)

b = 2^nk-n

化简，就是要求gcd(a,b)

可以肯定，gcd(a,b)肯定是2^d这种形式

由于k < 2ⁿ，d实际上就等于(k-1)!分解素因子后2的个数

如果k >= P，化简取模后有a = b,(但是取模后不能再继续化简为a = b = 1)

如果k < P，log(k)求得d，则：

b = 2^nk-n-d，快速幂可求，由于n * k会超过long long,可以先 % phi(P)，欧拉定理

a = 2^nk-n-d - (2ⁿ - 1) * (2ⁿ - 2) * ... * (2ⁿ - k + 1) / 2^d

后面这个可以先暴力O(k)求，再乘以2^d的逆元

注意 a = (a + P) % P防止a < 0 的情况

代码：

                                            
  //File Name: cf711E.cpp
  //Created Time: 2017年01月06日 星期五 00时05分30秒
                                   
#include <bits/stdc++.h>
#define LL long long
using namespace std;
const int P = (int)1e6 + 3;
LL qp(LL x,LL y){
    LL res = 1;
    for(;y>0;y>>=1){
        if(y & 1) res = res * x % P;
        x = x * x % P;
    }
    return res;
}
LL cal_phi(LL n){
    LL res = n;
    for(int i=2;i*i<=n;++i){
        if(n % i == 0){
            res -= res / i;
            while(n % i == 0)
                n /= i;
        }
    }
    if(n > 1) res -= res / n;
    return res;
}
void solve(LL n,LL k,LL &a,LL &b){
    LL now = 1;
    for(LL j=1;j<=n;++j){
        now *= 2;
        if(now >= k) break;
    }
    if(now < k) a = 1,b = 1;
    else{
//        puts("fffff");
        LL d = 0;
        now = k - 1;
        while(now >= 2){
            d += now / 2;
            now /= 2;
        }
        LL phi = P - 1;
        b = qp(2,((n % phi) * (k % phi) - n % phi - d % phi + 2 * phi) % phi);
        if(k >= P) a = b;
        else{
            now = qp(2,n % phi);
            a = 1;
            for(LL i=1;i<k;++i)
                a = (now - i + P) % P * a % P;
            now = qp(2,d % phi);
            a = a * qp(now,P - 2) % P;
            a = (b - a + P) % P;
        }
    }
}
int main(){
    LL n,k;
    cin >> n >> k;
    LL a,b;
    solve(n,k,a,b);
    printf("%lld %lld
",a,b);
    return 0;
}