BZOJ[3992][SDOI2015]序列统计 生成函数+NTT

33 篇文章 0 订阅
8 篇文章 0 订阅

传送门ber~

首先了解一下指标
看我瞎bb也可以
因为原根 g g 满足gi,gj(i,j(1,MOD1),ij)互不相同
则可以给每个数 i i 定义一个指标indi表示模意义下的 loggi log g ⁡ i ,并且在区间 [1,φ(MOD)] [ 1 , φ ( M O D ) ] 中是互不相同的
log l o g 类似,指标也满足 indijindi+indj i n d i ∗ j ≡ i n d i + i n d j 就可以把乘法弄成加法了
题目要求的 (a×b)modM=x ( a × b ) m o d M = x 等价于 (ind[a]+ind[b])modφ(M)=ind[x] ( i n d [ a ] + i n d [ b ] ) m o d φ ( M ) = i n d [ x ]
所以弄个生成函数再NTT,因为有取模,所以要在每次NTT后将 φ(m) φ ( m ) 项后面的系数挪到前面,剩下的快速幂解决就好了

代码如下:

#include<algorithm>
#include<ctype.h>
#include<cstdio>
#include<cmath>
#define MOD 1004535809
#define N 50020
#define int long long
using namespace std;
const bool DFT=false,IDFT=true;
const int G=3;
inline int read(){
    int x=0,f=1;char c;
    do c=getchar(),f=c=='-'?-1:f; while(!isdigit(c));
    do x=(x<<1)+(x<<3)+c-'0',c=getchar(); while(isdigit(c));
    return x*f;
}
typedef long long LL;
void exgcd(LL a,LL b,LL &x,LL &y){
    if(!b){
        x=1;y=0;
        return;
    }
    exgcd(b,a%b,y,x);
    y=y-a/b*x;
    return;
}
inline int GetNi(int k){
    LL x,y;
    exgcd(k,MOD,x,y);
    return (x+MOD)%MOD;
}
inline int qpow(LL x,LL k,LL p){
    LL sum;
    for(sum=1;k;k>>=1,x=x*x%p)
        if(k&1) sum=sum*x%p;
    return sum;
}
bool b[N];
int ind[N],pos[N],sum[N],p[N],q[N],a[N];
int n,m,x,s,t,g,len,inv_len,inv_g;
inline void NTT(int a[],bool mode){
    for(int i=0;i<len;i++)
        if(i<pos[i])
            swap(a[i],a[pos[i]]);
    int g=mode?inv_g:G;
    for(int i=2,mid=1;i<=len;i<<=1,mid<<=1){
        int wm=qpow(g,(MOD-1)/i,MOD);
        for(int j=0;j<len;j+=i){
            for(int k=j,w=1;k<mid+j;k++,w=(LL)w*wm%MOD){
                int l=a[k],r=(LL)w*a[k+mid]%MOD;
                a[k]=(l+r)%MOD;a[k+mid]=(l-r+MOD)%MOD;
            }
        }
    }
    if(mode==IDFT)
        for(int i=0;i<len;i++)
            a[i]=a[i]*inv_len%MOD;
    return;
}
inline void Mul(int a[],int b[],int c[]){
    for(int i=0;i<len;i++)
        p[i]=a[i],q[i]=b[i];
    NTT(p,DFT);NTT(q,DFT);
    for(int i=0;i<len;i++)
        c[i]=(LL)p[i]*q[i]%MOD;
    NTT(c,IDFT);
    for(int i=m-1;i<=m*2-4;i++)
        c[i-m+1]=(c[i-m+1]+c[i])%MOD,c[i]=0;
    return;
}
bool check(int x,int k) {
    for(int i=1;i<k;i++) b[i]=false;
    for(int i=1;i<k;i++){
        int t=qpow(x,i,k);
        if(b[t]) return false;
        b[t]=true;
    }
    return true;
}
inline int GetG(int k){
    int x=1;
    while(!check(x,k)) x++;
    return x;
}
main(){
    inv_g=GetNi(G);
    n=read();m=read();x=read();s=read();
    g=GetG(m);
    for(int i=0,now=1;i<m-1;i++,(now*=g)%=m) ind[now]=i;
    for(int i=1;i<=s;i++){
        t=read()%m;
        if(t) a[ind[t]]=1;
    }
    for(len=1;len<m<<1;len<<=1);
    for(int i=0;i<len;i++){
        pos[i]=pos[i>>1]>>1;
        if(i&1) pos[i]|=len>>1;
    }
    inv_len=GetNi(len);
    sum[0]=1;
    while(n){
        if(n&1) Mul(sum,a,sum);
        Mul(a,a,a);
        n>>=1;
    }
    printf("%d",sum[ind[x]]);
return 0;
}
  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值