zoukankan      html  css  js  c++  java
  • [P5488] 差分与前缀和

    Description

    给定一个序列,求它的 (k) 阶前缀和与 (k) 阶差分。

    Solution

    前缀和可以看做与 (frac 1 {1-x}) 卷积的形式,多项式快速幂加多项式求逆即可。

    差分可以看做与 (1-x) 卷积的形式(注意这里方向要反一下),多项式快速幂即可。

    #include <bits/stdc++.h>
    using namespace std;
    
    #define int long long
    const int N = 262150;
    const int mod = 1004535809;
    
    int qpow(int p,int q) {return (q&1?p:1)*(q?qpow(p*p%mod,q/2):1)%mod;}
    
    namespace cipolla {
    inline int le(int x) {return qpow(x,(mod-1)/2);}
    int w;
    struct comp {
        int x,y;
        comp(int a=0,int b=0) {x=a;y=b;}
    };
    comp operator + (comp a,comp b) {return comp((a.x+b.x)%mod,(a.y+b.y)%mod);}
    comp operator - (comp a,comp b) {return comp((a.x-b.x+mod)%mod,(a.y-b.y+mod)%mod);}
    comp operator * (comp a,comp b) {return comp((a.x*b.x+a.y*b.y%mod*w)%mod,(a.x*b.y+a.y*b.x)%mod);}
    comp operator ^ (comp a,int b) {comp o(1,0); for(;b;a=a*a,b>>=1) if(b&1) o=o*a; return o;}
    int calc(int x) {
        x%=mod;
        int a;
        while(true) {
            a=rand();
            w=(a*a-x+mod)%mod;
            if(le(w)==mod-1) break;
        }
        comp s=comp(a,1)^((mod+1)/2);
        return min(s.x,mod-s.x);
    }
    }
    
    namespace po {
    int rev[N],inv[N],w[N],sz;
    void presolve(int l) {
        int len=1;
        sz=0;
        while(len<l) len<<=1, ++sz;
        for(int i=1;i<len;i++) {
            inv[i]=(i==1?1:inv[mod%i]*(mod-mod/i)%mod);
            rev[i]=(rev[i>>1]>>1)|((i&1)<<(sz-1));
        }
        int wn=qpow(3,(mod-1)/len);
        w[len/2]=1;
        for(int i=len/2+1;i<len;i++) w[i]=w[i-1]*wn%mod;
        for(int i=len/2-1;i;i--) w[i]=w[i<<1];
    }
    int pre(int l) {int g; for(g=1;g<l;g<<=1); return g;}
    void ntt(int *a,int o,int n) {
        static unsigned long long s[N];
        int t=sz-__builtin_ctz(n),x;
        for(int i=0;i<n;i++) s[rev[i]>>t]=a[i];
        for(int l=1;l<n;l<<=1) for(int i=0;i<n;i+=l<<1) for(int j=0;j<l;j++) {
            x=s[i+j+l]*w[j+l]%mod;
            s[i+j+l]=s[i+j]+mod-x;
            s[i+j]+=x;
        }
        for(int i=0;i<n;i++) a[i]=s[i]%mod;
        if(o) {
            x=qpow(n,mod-2);
            for(int i=0;i<n;i++) a[i]=a[i]*x%mod;
            reverse(a+1,a+n);
        }
    }
    void mult(int n,int *x,int *y,int *z) {
        static int a[N],b[N];
        int l=pre(n<<1);
        for(int i=0;i<l;i++) {
            a[i]=(i<n?x[i]:0);
            b[i]=(i<n?y[i]:0);
        }
        ntt(a,0,l); ntt(b,0,l);
        for(int i=0;i<l;i++) z[i]=a[i]*b[i]%mod;
        ntt(z,1,l);
        for(int i=n;i<l;i++) z[i]=0;
    }
    void inve(int len,int *a,int *b) {
        if(len==1) *b=qpow(*a,mod-2);
        else {
            inve((len+1)/2,a,b);
            static int c[N];
            int n=pre(len<<1);
            for(int i=0;i<n;i++) i<len?c[i]=a[i]:b[i]=c[i]=0;
            ntt(b,0,n);
            ntt(c,0,n);
            for(int i=0;i<n;i++) b[i]=((b[i]+b[i]-b[i]*b[i]%mod*c[i])%mod+mod)%mod;
            ntt(b,1,n);
            for(int i=len;i<n;i++) b[i]=0;
        }
    }
    void sqrt(int n,int *a,int *b) {
        if(n==1) *b=cipolla::calc(*a);
        else {
            sqrt((n+1)/2,a,b);
            static int c[N];
            inve(n,b,c);
            mult(n,a,c,c);
            for(int i=0;i<n;i++) b[i]=(b[i]+c[i])*inv[2]%mod;
        }
    }
    void deri(int n,int *a,int *b) {
        for(int i=0;i<n-1;i++) b[i]=a[i+1]*(i+1)%mod;
        b[n-1]=0;
    }
    void inte(int n,int *a,int *b) {
        for(int i=n-1;i>0;--i) b[i]=a[i-1]*inv[i]%mod;
        b[0]=0;
    }
    void loge(int n,int *a,int *b) {
        static int c[N];
        inve(n,a,b);
        deri(n,a,c);
        mult(n,b,c,b);
        inte(n,b,b);
    }
    void expr(int n,int *a,int *b) {
        if(n==1) *b=1;
        else {
            expr((n+1)/2,a,b);
            static int c[N];
            loge(n,b,c);
            for(int i=0;i<n;i++) c[i]=(a[i]-c[i]+mod)%mod;
            c[0]=(c[0]+1)%mod;
            mult(n,b,c,b);
        }
    }
    }
    
    int n,k,t,a[N],b[N],c[N],d[N];
    
    signed main()
    {
        ios::sync_with_stdio(false);
    
        string ks;
        cin>>n>>ks>>t;
    
        for(int i=0;i<ks.length();i++)
        {
            k=(k*10+ks[i]-'0')%mod;
        }
    
        po::presolve((n+1)<<1);
        for(int i=0;i<n;i++) cin>>a[i];
    
        if(t==0)
        {
            for(int i=0;i<n;i++) b[i]=1;
            po::loge(n,b,c);
            for(int i=0;i<n;i++) c[i]=(c[i]*k)%mod;
            po::expr(n,c,d);
            for(int i=0;i<n;i++) c[i]=0;
            po::mult(n,a,d,c);
            for(int i=0;i<n;i++) cout<<c[i]<<" ";
        }
        else
        {
            b[0]=1;
            b[1]=mod-1;
            po::loge(n,b,c);
            for(int i=0;i<n;i++) c[i]=(c[i]*k)%mod;
            po::expr(n,c,d);
            for(int i=0;i<n;i++) c[i]=0;
            po::mult(n,a,d,c);
            for(int i=0;i<n;i++) cout<<c[i]<<" ";
        }
        cout<<endl;
        return 0;
    }
    
    
  • 相关阅读:
    个人觉得在前台比较通用的校验的提示方法
    平时面试总结
    sturts2批量导入excel中的信息入库大致代码
    在B/S系统中得到spring的ApplicationContext对象的一点小技巧
    Oracle同义词创建及其作用
    Oracle中创建dblink的方法
    真实项目中struts2 导出excel文件
    十二款jQuery插件(编辑器,图片,验证等)
    华丽到暴表的网站 cnn ecosphere
    Mysql 临时表
  • 原文地址:https://www.cnblogs.com/mollnn/p/13758238.html
Copyright © 2011-2022 走看看