P3301 [SDOI2013]方程 记录

· · 题解

你一眼看出来是一个母函数的题目，于是你列出这么一个柿子

[x^m]\prod_{i=1}^{n_1}(x+x^2+\ldots+x^{A_{i}})\prod_{i=n_1+1}^{n_1+n_2}(x^{A_i}+x^{A_i+1}+\ldots)\prod_{i=n_1+n_2+1}^n(1+x+\ldots) =[x^{m-\sum_{i=n_1+1}^{n_1+n_2}A_i-1}](\frac{x}{1-x})^{n-n_1}\prod_{i=1}^{n_1}(x+x^2+\ldots+x^{A_{i}})

然后你发现 m 很大，于是你使用亿些科技

现在是 14:02 ，我看我什么时候写完

现在是 19:33 ，我写完了

#define int long long
int ksm(int a,int b,int p){
    int ans=1;
    for(;b;b>>=1,a=a*a%p)
    if(b&1)ans=ans*a%p;
    return ans;
}
int _exgcd(int a,int b,int &x,int &y){
  if(!b){x=1,y=0;return a;}
  int d=_exgcd(b,a%b,y,x);
  y-=a/b*x;return d;
}
int _exCRT(int b[],int a[],int n){
    int ans=0,g,c,x,y,lcm=1;for(int i=1;i<=n;i++){
        g=_exgcd(lcm,a[i],x,y);c=(b[i]-ans)%a[i];c=(c+a[i])%a[i];
    if(c%g)return -1;x=((x*c/g)%(a[i]/g)+(a[i]/g))%(a[i]/g);
    ans=ans+x*lcm;lcm=lcm*a[i]/g;ans=(ans%lcm+lcm)%lcm;
  }return ans;
} 
int _inv(int x,int p){int a,b;_exgcd(x,p,a,b);return (a+p)%p;}
int f[300000];
int _calc(int n,int x,int p){
    if(n==0)return 1;int ans=ksm(f[p],n/p,p);
    for(int i=n/p*p+1;i<=n;i++)
    if(i%x)ans=i%p*ans%p;return ans*_calc(n/x,x,p)%p;
}
int _mutilucas(int n,int m,int x,int p){
    if(n<m)return 0;
    int cnt=0;for(int i=n;i;i/=x)cnt+=i/x;
    for(int i=m;i;i/=x)cnt-=i/x;for(int i=n-m;i;i/=x)cnt-=i/x;
    f[0]=1;
    for(int i=1;i<=p;i++)if(i%x)f[i]=f[i-1]*i%p;else f[i]=f[i-1];
    return ksm(x,cnt,p)*_calc(n,x,p)%p*_inv(_calc(m,x,p),p)%p*_inv(_calc(n-m,x,p),p)%p;
}
int T,p,n,n1,n2,m;
int A[20],x;
int solve(int n,int m,int x,int p){
    int ans=0;
    for(int i=0;i<1<<n1;i++){
            int cc=n,op=(__builtin_popcount(i)&1)?-1:1;
            for(int j=0;j<n1;j++)
            if(i>>j&1)cc-=A[j];
            ans=(ans+p+op*_mutilucas(cc,m,x,p))%p;
    }
    return ans;
}
int a[20],b[20],c[20],cnt=0;
void shai(int p){
    for(int i=2;i*i<=p;i++)
    if(p%i==0){b[++cnt]=1;c[cnt]=i;while(p%i==0)p/=i,b[cnt]*=i;}
    if(p>1)b[++cnt]=p,c[cnt]=p;
}
int exlucas(int n,int m,int p){
    if(n<m)return 0;
    for(int i=1;i<=cnt;i++)a[i]=solve(n,m,c[i],b[i]);
    return _exCRT(a,b,cnt);
}
signed main(){
#ifndef ONLINE_JUDGE
    freopen("testdata.in","r",stdin);
#endif
    for(cin>>T>>p,shai(p);T;T--){
        cin>>n>>n1>>n2>>m;
        for(int i=0;i<n1;i++)cin>>A[i];
        for(int i=0;i<n2;i++)cin>>x,m-=x-1;
        cout<<exlucas(m-1,n-1,p)<<endl;
    }
#ifndef ONLINE_JUDGE
    cerr<<endl<<(double)clock()/CLOCKS_PER_SEC;
#endif
}

你发现这个多项式很美妙，但是 m\le10^9 ，于是你毅然弃疗。

然后你发现中间的 n_2 个东西很简单，在多项式里你就直接意识到可以用 m 减掉，然后就变成了一般情况，一般情况就直接转化为插板法。

此时你发现上面的 n_1 个玩意儿很难，于是你想起了 \text{devin} 说的考虑容斥，全部转化为一般情况。

几个不满足一式，几个不满足二式，几个不满足一式二式，加加减减。

还有容斥可以直接在 \text{exlucas} 里做就行了。

然后你一交， \text{70pts TLE} 。

于是你开始卡常，你在 \text{calc} 中预处理阶乘，因为要多次用到，然后质因数提前分解，因为全局 p 相等。

不论怎么说你要坚信复杂度跑不满， O(2^{n_1}n^2_1)\times 扩展卢卡斯复杂度。如果想做得更好好像可以用 \text{dfs} 优化掉一个 O(n^2_1)