luogu P3726 [AH2017/HNOI2017]抛硬币

题面传送门

exLucas的第一道例题。

我们考虑我们要算的东西，显然是\(\sum\limits_{i=0}^{a}{\sum\limits_{j=0}^{b}{C_{a}^iC_b^j[i>j]}}\)。

为了快速计算这个东西，我们考虑一种特殊情况\(a=b\)，可以发现，除了平局，输赢其实是对称的，因此我们只要把平局减去，然后除二即可。

平局的方案数显然是\(\sum\limits_{i=0}^{n}{(C_{n}^i)^2}\)。这个东西可以写成\(\sum\limits_{i=1}^{n}{C_{n}^iC_{n}^{n-i}}\)，因此根据范德蒙特卷积，这个东西等于\(C_{2n}^n\)。

现在压力来到了快速计算组合数，一个非常棘手的事情是这个题的模数似乎不是质数，而且和模数的质因数分解还非常小因此一定不互质。

这时候就需要exLucas了，具体的，它只是将模数拆成质数幂次乘积然后CRT合并而已。

考虑\(a\ne b\)的情况，这种情况下似乎不是很好做，因为输赢并不对等。

那我们把赢的多出来的情况加上去然后除二不就好了吗。

同理根据上面的推法，可以发现答案要加上\(\sum\limits_{i=b+1}^{a-1}{C_{a+b}^i}\)。直接exLucas计算即可，有点卡常。

时间复杂度\(O(T(a-b)\log ^2n)\)

code:

#include<bits/stdc++.h>
#include<iostream>
#define I inline
#define ll long long
#define db double
#define lb long db
#define N (2000000+5)
#define M ((1<<16)+5)
#define K (1500+5)
#define mod 1000000007
#define Mod (mod-1)
#define eps (1e-5)
#define ull unsigned ll
#define it iterator
#define Gc() getchar() 
#define Me(x,y) memset(x,y,sizeof(x))
#define Mc(x,y) memcpy(x,y,sizeof(x))
#define d(x,y) ((k+1)*(x)+(y))
#define R(n) (1ll*rand()*rand()%(n)+1)
#define Pc(x) putchar(x)
#define LB lower_bound
#define UB upper_bound
#define PB push_back
using namespace std;
ll n,m,F[N];int k,p1,p2;
I ll mpow(ll x,ll y,int p){ll Ans=1;while(y) y&1&&(Ans=Ans*x%p),y>>=1,x=x*x%p;return Ans;}
I void EG(ll x,ll y,ll &a,ll &b){if(!y){a=1;b=0;return;}EG(y,x%y,b,a);b-=x/y*a;} 
I ll GI(ll x,ll p){ll y,z;EG(x,p,y,z);return (y%p+p)%p;}
I ll CRT(ll p1,ll S1,ll p2,ll S2){return (S1*p2%(p1*p2)*GI(p2,p1)+S2*p1%(p1*p2)*GI(p1,p2))%(p1*p2);}
I ll GG(ll n,int p){if(!n) return 0;return n/p+GG(n/p,p);}I ll GF(ll n,ll p,ll pk){if(!n) return 1;return GF(n/p,p,pk)*mpow(F[pk],n/pk,pk)%pk*F[n%pk]%pk;}
I ll calc(ll n,ll m,ll p,ll pk){ll B=GG(n,p)-GG(m,p)-GG(n-m,p);if(B>=k+1) return 0;return GF(n,p,pk)*GI(GF(m,p,pk),pk)%pk*GI(GF(n-m,p,pk),pk)%pk*mpow(p,B,pk)%pk;}
I void Solve(){
	int i,j;ll T1,T2;p1=2;p2=1;for(i=1;i<=k;i++) p1*=2,p2*=5;T1=mpow(2,n+m,p1);T2=mpow(2,n+m,p2);
	for(F[0]=i=1;i<=p1;i++) F[i]=F[i-1]*(i%2?i:1)%p1;if(n==m)T1-=calc(2*n,n,2,p1);else for(ll i=m+1;i<n;i++) T1+=calc(n+m,i,2,p1);T1=(T1%p1+p1)%p1;
	for(F[0]=i=1;i<=p2;i++) F[i]=F[i-1]*(i%5?i:1)%p2;if(n==m)T2-=calc(2*n,n,5,p2);else for(ll i=m+1;i<n;i++) T2+=calc(n+m,i,5,p2);T2=(T2%p2+p2)%p2;
	printf("%0*lld\n",k,CRT(p1,T1,p2,T2)/2);
}
int main(){
	freopen("1.in","r",stdin);
	while(~scanf("%lld%lld%d",&n,&m,&k)) Solve();
}