【清华集训2017模拟】Catalan
2017-08-20 22:04
316 查看
Description
求Cnmod3814697265625(518)其中Cn为卡特兰数第n项n<=10^18,T<=10
Solution
这么大的组合数取模啊。。。。以前真没见过首先我们知道Ans=Cn2nn+1
根据套路我们只需要把n!写成5^e*f的形式,然后就可以用逆元直接算了
但是如何计算?
首先我们把所有5的倍数提取出来,那么有⌊n5⌋个5,并且剩下的数变成了⌊n5⌋!
这样就可以递归下去计算了
那么如何处理剩余的部分?
我们在处理乘法的时候忽略5的倍数,设L=m−−√o=59,那么我们相当于要求
∏i=0⌊nL⌋−1∏j=0L−1(iL+j)
假如我们把i,l当做常数,那么我们可以发现后面的东西可以写成(a*i*L+b)的形式
a和b只和L有关,可以预处理
那么这样我们处理出a,b,再处理出一个前缀积,这个东西就可以O(1)计算了
预处理的复杂度是O(L)的
多出来的那一部分直接处理就好了
但是这样只能处理n<=mo的情况,n>mo的情况处理不了,怎么办呢?
我们发现n>mo的时候我们可以把1~n的数按mo分类,其中每一组在模意义下都是同余的
这样直接快速幂计算即可
总复杂度O(L+T log^2N)
Code
#include <cstdio> #include <cstring> #include <algorithm> #define fo(i,a,b) for(int i=a;i<=b;i++) using namespace std; typedef long long ll; typedef double db; const ll mo=3814697265625; const int L=1953125; struct note{ int e;ll f; note(int _e=0,ll _f=0) {e=_e;f=_f;} }; ll f[L][2],sum[L],n; ll mult(ll a,ll b) { ll a1=a/L,b1=a%L,a2=b/L,b2=b%L; ll tmp=a1*a2%mo*L%mo*L%mo; (tmp+=a1*L%mo*b2%mo)%=mo; (tmp+=a2*L%mo*b1%mo)%=mo; (tmp+=b1*b2%mo)%=mo; return tmp; } ll mi(ll x,ll y) { ll z=1; for(;y;y/=2,x=mult(x,x)) if (y&1) z=mult(z,x); return z; } ll F(ll n) { ll a=n/L,b=n%L,ans=1; if (a) ans=mult(ans,sum[a-1]); ans=mult(ans,(mult(mult(f[b][1],a),L)+f[b][0])%mo); return ans; } note solve(ll n) { if (!n) return note(0,1); note x=solve(n/(ll)5),ans; if (n>mo) { ll a=n/mo,b=n%mo; ans.f=mult(x.f,mult(mi(F(mo),a),F(b))); } else ans.f=mult(x.f,F(n)); ans.e=(ll)(x.e+n/(ll)5)%18; return ans; } void exgcd(ll a,ll b,ll &x,ll &y) { if (!b) { x=1;y=0; return; } ll xx,yy; exgcd(b,a%b,xx,yy); x=yy;y=xx-a/b*yy; } ll get_inv(ll x) { ll a,b; exgcd(x,mo,a,b); return ((a%=mo)+=mo)%=mo; } int main() { freopen("catalan.in","r",stdin); freopen("catalan.out","w",stdout); f[0][0]=1; fo(i,1,L-1) { if (!(i%5)) { f[i][0]=f[i-1][0]; f[i][1]=f[i-1][1]; continue; } f[i][0]=mult(f[i-1][0],(ll)i); f[i][1]=(f[i-1][0]+mult(f[i-1][1],(ll)i))%mo; } sum[0]=f[L-1][0]; fo(i,1,L-1) sum[i]=mult(sum[i-1],(mult(mult(f[L-1][1],(ll)i),(ll)L)+f[L-1][0])%mo); int ty; for(scanf("%d",&ty);ty;ty--) { scanf("%lld",&n); note a=solve(2*n); note b=solve(n); (a.e+=18-2*b.e%18)%=18; ll ni=get_inv(b.f); a.f=mult(a.f,mult(ni,ni)); ll nn=n+1; for(;!(nn%5);nn/=5) (a.e+=17)%=18; a.f=mult(a.f,get_inv(nn)); ll c=mult(a.f,mi(5,a.e)); printf("%lld\n",c); } }
相关文章推荐
- JZOJ 5485. 【清华集训2017模拟11.26】字符串
- [JZOJ5500]【清华集训2017模拟12.10】营养餐
- 【JZOJ5295】【清华集训2017模拟】Create
- 【清华集训2017模拟】Sequence
- JZOJ5489. 【清华集训2017模拟11.28】海明距离
- jzoj5498 【清华集训2017模拟12.10】大佬的难题 巧妙容斥
- 【清华集训2017模拟】Sequence
- JZOJ 5483. 【清华集训2017模拟11.26】简单路径
- 【清华集训2017模拟11.26】字符串
- 5296. 【清华集训2017模拟】Sequence 树套树
- JZOJ5483. 【清华集训2017模拟11.26】简单路径
- JZOJ 5500. 【清华集训2017模拟12.10】营养餐
- 【清华集训2017模拟11.29】K小数查询
- 【JZOJ5316】【清华集训2017模拟8.19】merge
- [jzoj]5483. 【清华集训2017模拟11.26】简单路径
- JZOJ 5489. 【清华集训2017模拟11.28】海明距离
- 【清华集训2017模拟12.10】回文串
- JZOJ5484. 【清华集训2017模拟11.26】快乐树
- 【清华集训2017模拟12.09】塔
- jzoj5484 【清华集训2017模拟11.26】快乐树 (发现性质的dp)