ZS and The Birthday Paradox
题目链接：http://codeforces.com/contest/711/problem/E
数学题（Legendre's formula）
这题是以生日悖论（如果有23个或23个以上的人，那么至少有两个人的生日相同的概率要大于50%）为背景的数学题。公式本身不难推，主要是对大数的处理。
首先，我们需要找到分子和分母的公因数，这里用到了Legendre's formula（举个例子：10!的因数中质数3的个数为[10/3+10/(3^2)]）。因为若2^n-x与2^n有公因数，那么x与2^n有相同的公因数，所以求(2^n)(2^n-1)(2^n-2)*...*(2^n-k+1)与(2^n)^k的公因数，就转化为了求(k-1)!与(2^n)^k的公因数。
之后就是分别求分子和分母：
对于分母来说，只需要用快速幂（也可以用费马小定理）就可以很容易的求出，再乘上公因数的逆元即可；
而对于分子来说，稍微有点麻烦：
1.如果k>=M，即(2^n)(2^n-1)(2^n-2)*...*(2^n-k+1)中至少有连续的M个整数，
那么(2^n)(2^n-1)(2^n-2)*...*(2^n-k+1)一定为M的倍数，所以它被M求模后余0；
2.如果k<M，因为M很小，所以枚举一下，就可以求出分子，再乘上公因数的逆元即可。
总的时间复杂度为O（M+lgk+lgn）
（感谢游少半夜教我证明Orz）
证明：(A/B)modM=(A*(BmodM)')modM，其中B'为B在M下的逆元
令B=b1*b2*b3*...*bn，则((BmodM)')modM
=((b1*b2*b3*...*bn mod M)')modM
=(b1modM*b2modM*...*bn mod M）'modM
=(b1*b2*...*bn)modM*(b1modM*b2modM*...*bn mod M）'modM*(b1*b2*...*bn)'modM
=(b1modM*b2modM*...*bn mod M)*(b1modM*b2modM*...*bn mod M）'modM*(b1*b2*...*bn)'modM
=1*(b1*b2*...*bn)'modM=(b1*b2*...*bn)'modM
=b1*b1'modM*b2*b2'modM...bn*bn'mod M*(b1*b2*...*bn)'modM
=(b1*b2*...*bn)modM*(b1'modM*b2'modM...bn'mod M)*(b1*b2*...*bn)'modM
=b1'modM*b2'modM...bn'mod M
代码如下：
#include<cstdio>
#include<iostream>
#define M (long long)(1e6+3)
using namespace std;
typedef long long LL;
LL n,k,cnt,molecular,numerator,x,y;
LL exGCD(LL a,LL b){
if(b==){
x=,y=;
return a;
}
LL r=exGCD(b,a%b);
LL temp=x;
x=y;
y=temp-(a/b)*y;
return r;
}
LL mod(LL a,LL b){
LL base=a,temp=;
while(b){
if(b&)temp=(temp*base)%M;
base=(base*base)%M;
b>>=;
}
return temp;
}
int main(void){
cin>>n>>k;
if(n<=&&k>(((LL))<<n)){//总人数大于总天数
cout<<""<<" "<<""<<endl;
return ;
}
LL temp=;
while(k->=temp){//根据Legendre's formula,求出(k-1)!的因数中质数2的个数
cnt+=((k-)/temp);
temp<<=;
}
if(cnt){//若有公因数，则求出公因数的逆元
LL gcd=mod(,cnt);
exGCD(gcd,M);
x=(x+M)%M;
}else x=;//若没有公因数，则令x=1，来消除对后面计算的影响
temp=mod(,n);//计算2^n
numerator=(mod(temp,k-)*x)%M;//计算(2^n)^(k-1)
if(k>=M)molecular=;//若分子出现连续的M个整数，则分子一定为M的倍数
else{
molecular=;
for(LL i=;i<k;++i){
molecular=((temp-i+M)%M*molecular)%M;//计算分子
}
molecular=(molecular*x)%M;//除以公因数
}
molecular=(numerator-molecular+M)%M;
cout<<molecular<<" "<<numerator<<endl;
}