LGP5075【JSOI2012】分零食

、

• 题解：

  • 令$F$为欢乐度$f(x) = Ox^2 + Sx + U$的生成函数，常数项为$0$；
  • 令$G(x) = \sum_{i=0}^{A} F^i (x) $
    
  • $ans = [x^M]G;$
  • 模数比较麻烦所以我用的分治求:
    
  • 如果现在要求$0$到$n-1$的$G_{n} = \sum_{i=0}^{n-1}F^{i}$和$F_{n} = F^{n} $，假设n为偶数；
  • 那么分治求出关于$n/2$的答案$G_{\frac{n}{2}}$和$F_{\frac{n}{2}}$
  • $$G_{n} = (F_{\frac{n}{2}}+1)G_{\frac{n}{2}}  , F_{n} = F_{\frac{n}{2}}^2$$
  • 如果$n$是奇数先算用上述操作算$n-1$，再把$F_{n-1}$补加给$G_{n-1}$得到$G_{n}$,最后$F_{n-1}$另外乘一次得到$F_{n}$；
  • 和快速幂的思想差不多；
  • 

     1 #include<bits/stdc++.h>
     2 #define ld double
     3 using namespace std;
     4 const int N=40010;
     5 const ld pi=acos(-1);
     6 int M,P,A,O,S,U,len,L,rev[N];
     7 struct C{
     8     ld x,y;
     9     C(ld _x=0,ld _y=0):x(_x),y(_y){};
    10     C operator +(const C&A)const{return C(x+A.x,y+A.y);}
    11     C operator -(const C&A)const{return C(x-A.x,y-A.y);}
    12     C operator *(const C&A)const{return C(x*A.x-y*A.y,x*A.y+y*A.x);}
    13     C operator /(const ld&A)const{return C(x/A,y/A);}
    14 }f[N],g[N],t[N];
    15 int cal(int x){return (O*x*x+S*x+U)%P;}
    16 void fft(C*a,int f){
    17     for(int i=0;i<len;++i)if(i<rev[i])swap(a[i],a[rev[i]]);
    18     for(int i=1;i<len;i<<=1){
    19         C wn=C(cos(pi/i),f*sin(pi/i));
    20         for(int j=0;j<len;j+=i<<1){
    21             C w=C(1,0);
    22             for(int k=0;k<i;++k,w=w*wn){
    23                 C x=a[j+k],y=w*a[j+k+i];
    24                 a[j+k]=x+y,a[j+k+i]=x-y;
    25             }
    26         }
    27     }
    28     if(!~f)for(int i=0;i<len;++i){
    29         a[i]=a[i]/len;
    30         a[i].x=int(a[i].x+0.1)%P;
    31         a[i].y=0;
    32     }
    33 }
    34 void solve(int A){
    35     if(A==1){g[0].x=1;return;}
    36     solve(A>>1);
    37     fft(f,1);fft(g,1);
    38     for(int j=0;j<len;++j)g[j]=g[j]*(f[j]+C(1,0)),f[j]=f[j]*f[j];
    39     fft(f,-1);fft(g,-1);
    40     for(int j=M+1;j<len;++j)f[j].x=g[j].x=0; 
    41     if(A&1){
    42         for(int j=0;j<=M;++j)g[j].x=(int)(g[j].x+f[j].x+0.1)%P; 
    43         fft(f,1);for(int j=0;j<len;++j)f[j]=f[j]*t[j];
    44         fft(f,-1);for(int j=M+1;j<len;++j)f[j].x=0;
    45     }
    46 }
    47 int main(){
    48 //    freopen("P5075.in","r",stdin);
    49 //    freopen("P5075.out","w",stdout);
    50     scanf("%d%d%d%d%d%d",&M,&P,&A,&O,&S,&U);
    51     for(int i=1;i<=M;++i)t[i].x=f[i].x=cal(i%P);
    52     for(len=1;len<=M<<1;len<<=1,L++);
    53     for(int i=0;i<len;++i)rev[i]=(rev[i>>1]>>1)|((i&1)<<(L-1));
    54     fft(t,1);solve(min(A,M)+1);
    55     printf("%d\n",(int)(g[M].x+0.1)%P);
    56     return 0;
    57 }

    View Code

     

转载于:https://www.cnblogs.com/Paul-Guderian/p/10433414.html