zoukankan      html  css  js  c++  java
  • [原创题] 雾

    题目:http://162.105.80.126/contest/%E3%80%90%E5%BC%B1%E7%9C%81%E8%83%A1%E7%AD%96%E3%80%91Round%20%233/%E9%9B%BE

    感谢绍兴一中三位神犇前来捧(nue)场,感谢山东神犇。。。。。。

    让我们来看一下题目,

    $$Ans = sum_{i=1}^{n}sum_{j=1}^{i-1}prod_{k=1}^{min(i-j,j)}[sum_{l=1}^k a(l)]^{varphi(j-k)+varphi(i-j-k)}prod_{i=min(i-j,j)+1}^{max(i-j,j)}f(i,j,k)$$
    若$i-j>j$ , $f(i,j,k)=[sum_{l=1}^ka(l)]^{varphi(i-j-k)}$

    否则 $f(i,j,k)=[sum_{l=1}^ka(l)]^{varphi(j-k)}$。

    注意到$sum_{l=1}^k a(l)$和其他的式子并没有关系,所以可以预先前缀和,下面我们所说的$a(i)$都是前缀和后的$a(i)$。
    $$Ans = sum_{i=1}^{n}sum_{j=1}^{i-1}prod_{k=1}^{min(i-j,j)}a(k)^{varphi(j-k)+varphi(i-j-k)}prod_{i=min(i-j,j)+1}^{max(i-j,j)}f(i,j,k)$$

    $10 points$: 按照定义$O(n^3)$

    显然这个式子是难以优化的,我们考虑加以变形,注意到$f(i,j,k)$函数的定义恰好和前面的连乘的内容相似,我们可以将$f(i,j,k)$代入化简,然后将左侧的连乘拆开为
    $$prod_{k=1}^{min(i-j,j)}a(k)^{varphi(j-k)}cdot a(k)^{varphi(i-j-k)}$$
    然后做变形
    $$prod_{k=1}^{min(i-j,j)}a(k)^{varphi(j-k)} prod_{k=1}^{min(i-j,j)}a(k)^{varphi(i-j-k)}$$
    然后发现 $f(i,j,k)$ 其实可以接到左侧式子或者右侧式子(取决于$i-j$和$j$的大小关系),整理得到:
    $$Ans = sum_{i=1}^{n}sum_{j=1}^{i-1}prod_{k=1}^{j}a(k)^{varphi(j-k)}prod_{i=1}^{i-j}a(k)^{varphi(i-j-k)}$$

    然后,仔细观察后面那一大坨连乘和前面的求和无关。
    设$C(n) = prod_{i=1}^{n}a(i)^{varphi(n-i)}$
    然后发现原式变成了:
    $$Ans = sum_{i=1}^{n}sum_{j=1}^{i-1}{C(i)cdot C(i-j)}$$
    (注意这里的卷积形式,要在首项补$0$,保证$j$是从$i-1$结束)

    $30 points$:
    假如你不会$FFT$也不会$NTT$那么你可以$O(n^2)$预处理,计算时直接用。
    总复杂度:$O(n^2)$

    标准的卷积形式,可以用$NTT$ $O(nlogn)$ 解决。
    然后回来看$C(n)$,左右用原根代换:
    $$ gn^{c'(n)} = prod_{i=1}^{n}gn^{a'(i) cdot {varphi(n-i)}}$$
    求出 $C'(n)$ 之后我们就可以 $O(nlogn)$ 快速幂求出 $C(n)$ 了。
    而注意到其指数恰为卷积的形式,但是根据费马小定理,应该取$mod-1$的余,模数不太对,不能用 $NTT%$,而 $FFT$ 又炸精度。

    考虑别的方法,由于数据中的指数是很小的,最大的情况下数字可以是$4 cdot 115000000000$,取一个超级大的质数,我是$1899956092796929$,原根为$7$,然后就可以了吗,朴素乘法会溢出,快速乘或者用$long double$的大数相乘取余函数。(然而快速乘多$logn$;$long double$乘有概率挂,所以不用是上策)

    此外,确定一个数字是原根的多少次方时可以用$BSGS$ $O(nsqrt{n})$注意这里的n小于等于$5cdot 10^4$.

    总效率$O(nlog^2n)$还是不能通过所有的数据呀,预计$40$~$50$分。

    $50points$:
    最后,我们发现我们可以选两个$10^9$级的满足条件的质数,最后$CRT$合并。
    然后就可以$50$了(此处应$orz$ $yjp$)

    $70points$:
    继续改进我们的算法,可以发现更本不用$BSGS$,完全可以预处理出所有的值,用一个$map$或$hash$存下来,取时$O(logn)$查询,总体复杂度:$O(nlogn)$。

    $100points$:
    可以发现,之所以$70$分,并不是因为复杂度,而是$NTT$太慢了,注意到指数很小,所以可以$FFT$,或者卡卡常数即可得到全分。

    注意精度。

    #include <cstdio>
    #include <cmath>
    #include <ctime>
    #include <algorithm>
    #include <map>
    
    #define N 400010
    #define LL unsigned long long
    #define mod 1004535809LL
    #define gn 3
    #define mod2 998244353LL
    
    using namespace std;
    
    const double pi=acos(-1.0);
    
    struct Ex{
        double r, i;
        Ex(){}
        Ex(double r, double i):r(r), i(i){}
        friend Ex operator+(const Ex &a,const Ex &b){return Ex(a.r+b.r,a.i+b.i);}
        friend Ex operator-(const Ex &a,const Ex &b){return Ex(a.r-b.r,a.i-b.i);}
        friend Ex operator*(const Ex &a,const Ex &b){return Ex(a.r*b.r-a.i*b.i,a.r*b.i+a.i*b.r);}
    }a0[N],b0[N],c0[N];
    
    map<LL,int> MT;
    
    inline LL qmult(LL a,LL b,LL P){
        LL ans=0;
        for(;b;b>>=1,a=(a<<1)%P)
            if(b&1) ans=(ans+a)%P;
        return ans;
    }
    
    inline void swap(LL &a,LL &b){
        if(a^b) a^=b,b^=a,a^=b;
    }
    
    inline LL qpow(LL x,LL n,LL P){
        LL ans=1;
        for(;n;n>>=1,x=x*x%P)
            if(n&1) ans=ans*x%P;
        return ans;
    }
    
    int R[N];
    
    inline LL Mod(LL x,LL P){
        if(x>=P) x-=P;
        if(x>=P) x-=P;
        if(x>=P) x%=P;
        return x;
    }
    
    LL wt[N];
    
    inline void fnt(LL *x,int n,int t,LL P){
        for(int i=0;i<n;i++) if(i<R[i]) swap(x[i],x[R[i]]);
        for(int m=1;m<n;m<<=1){
            LL wn=qpow(gn,(P-1)/(m<<1),P);
            if(t==-1) wn=qpow(wn,P-2,P);
            wt[0]=1;
            for(int i=1;i<m;i++) wt[i]=Mod(wt[i-1]*wn,P);
            for(int k=0;k<n;k+=(m<<1))
                for(int i=0;i<m;i++){
                    LL &A=x[i+m+k];
                    LL &B=x[i+k],t=Mod(A*wt[i],P);
                    A=Mod(B-t+P,P); B=Mod(B+t,P);
                }
        }
        if(t==-1){
            LL tmp=qpow(n,P-2,P);
            for(int i=0;i<n;i++)
                x[i]=x[i]*tmp%P;
        }
    }
    
    inline void fft(Ex x[],int n,int t){
        for(int i=0;i<n;i++) if(i<R[i]) swap(x[i],x[R[i]]);
        for(int m=1;m<n;m<<=1){
            Ex wn(cos(pi/m*t),sin(pi/m*t));
            for(int k=0;k<n;k+=(m<<1)){
                Ex wt(1,0);
                for(int i=0;i<m;i++,wt=wt*wn){
                    Ex &A=x[k+m+i],&B=x[k+i],t=wt*A;
                    A=B-t; B=B+t;
                }
            }
        }
        if(t==-1) for(int i=0;i<=n;i++) x[i].r/=n;
    }
    
    LL a1[N],b1[N],c1[N];
    
    inline void multpoly2(LL *A,LL *B,LL *C,int n){
        int nt,L=0;
        for(nt=1;nt<=n+n;nt<<=1) L++;
        for(int i=0;i<nt;i++) R[i]=(R[i>>1]>>1)|((i&1)<<(L-1));
        for(int i=0;i<nt;i++) a0[i].r=b0[i].r=c0[i].r=0;
        for(int i=0;i<nt;i++) a0[i].i=b0[i].i=c0[i].i=0;
        for(int i=0;i<n;i++) a0[i].r=A[i],b0[i].r=B[i];
        fft(a0,nt,1);
        fft(b0,nt,1);
        for(int i=0;i<nt;i++) c0[i]=a0[i]*b0[i];
        fft(c0,nt,-1);
        for(int i=0;i<n;i++) C[i]=(c0[i].r+0.5);
    }
    
    inline void multpoly1(LL *A,LL *B,LL *C,int n,LL P){
        int nt,L=0;
        for(nt=1;nt<=n+n;nt<<=1) L++;
        for(int i=0;i<nt;i++) R[i]=(R[i>>1]>>1)|((i&1)<<(L-1));
        for(int i=0;i<nt;i++) a1[i]=b1[i]=c1[i]=0;
        for(int i=0;i<n;i++) a1[i]=A[i],b1[i]=B[i];
        fnt(a1,nt,1,P);
        fnt(b1,nt,1,P);
        for(int i=0;i<nt;i++) c1[i]=a1[i]*b1[i]%P;
        fnt(c1,nt,-1,P);
        for(int i=0;i<n;i++) C[i]=c1[i];
    }
    
    inline LL CRT(LL a1,LL a2){
        LL ans=0,M=mod*(LL)mod2;
        ans+=qmult(qmult(a1,mod2,M),qpow(mod2,mod-2,mod),M);
        ans+=qmult(qmult(a2,mod,M),qpow(mod,mod2-2,mod2),M);
        return (ans%M+M)%M;
    }
    
    int n,phi[N];
    LL ans=0,a[N],b[N];
    LL c[N],ap[N];
    double ti;
    
    int main(){
        scanf("%d",&n);
        phi[1]=1;
        for(int i=2;i<=n;i++){
            if(phi[i]) continue;
            for(int j=i;j<=n;j+=i){
                if(!phi[j]) phi[j]=j;
                phi[j]=phi[j]/i*(i-1);
            }
        }
        for(int i=0;i<=100000;i++)
            MT[qpow(gn,i,mod)]=i;
        for(int i=0;i<n;i++) scanf("%lld",&a[i]);
        for(int i=1;i<n;i++) a[i]=(a[i]+a[i-1])%mod;
        for(int i=0;i<n;i++) ap[i]=MT[a[i]];
        for(int i=0;i<n;i++) b[i]=phi[i]%23;
        multpoly2(ap,b,c,n);
        for(int i=0;i<n;i++) c[i]=c[i]%(mod-1);
        for(int i=0;i<n;i++) c[i]=qpow(gn,c[i],mod);
        for(int i=n-1;~i;i--) c[i+1]=c[i],c[i]=0;
        multpoly1(c,c,c,n+1,mod);
        for(int i=0;i<=n;i++) ans=(ans+c[i])%mod;
        printf("%lld
    ",ans);
        return 0;
    }
    View Code
  • 相关阅读:
    SPOJ 694 (后缀数组) Distinct Substrings
    POJ 2774 (后缀数组 最长公共字串) Long Long Message
    POJ 3693 (后缀数组) Maximum repetition substring
    POJ 3261 (后缀数组 二分) Milk Patterns
    UVa 1149 (贪心) Bin Packing
    UVa 12206 (字符串哈希) Stammering Aliens
    UVa 11210 (DFS) Chinese Mahjong
    UVa (BFS) The Monocycle
    UVa 11624 (BFS) Fire!
    HDU 3032 (Nim博弈变形) Nim or not Nim?
  • 原文地址:https://www.cnblogs.com/lawyer/p/4564231.html
Copyright © 2011-2022 走看看