Hdu 5730 Shell Necklace(cdq+fft)

Shell Necklace

Time Limit: 16000/8000 MS (Java/Others) Memory Limit: 65536/65536 K (Java/Others)

Problem Description
Perhaps the sea‘s definition of a shell is the pearl. However, in my view, a shell necklace with n beautiful shells contains the most sincere feeling for my best lover Arrietty, but even that is not enough.

Suppose the shell necklace is a sequence of shells (not a chain end to end). Considering i continuous shells in the shell necklace, I know that there exist different schemes to decorate the i shells together with one declaration of love.

I want to decorate all the shells with some declarations of love and decorate each shell just one time. As a problem, I want to know the total number of schemes.

Input
There are multiple test cases(no more than 20 cases and no more than 1 in extreme case), ended by 0.

For each test cases, the first line contains an integer n, meaning the number of shells in this shell necklace, where 1≤n≤105. Following line is a sequence with n non-negative integer a1,a2,…,an, and ai≤107 meaning the number of schemes to decorate i continuous shells together with a declaration of love.

Output
For each test case, print one line containing the total number of schemes module 313(Three hundred and thirteen implies the march 13th, a special and purposeful day).

Sample Input
3
1 3 7
4
2 2 2 2
0

Sample Output
14
54

题意:dp[n] = sigma（a[i]*dp[n-i]）, 给出a1…..an， 求dp[n]。 n为1e5.

思路:这个式子直接计算是n^2的,我们可以考虑一下fft进行优化,
如果直接用n次fft的话,时间复杂度会变为n^2logn,还没有暴力来的快
我们再利用cdq分治进行一下优化,比如现在已经知道了dp[l],dp[l+1]…dp[mid]
我们只要能处理前半部分的答案对后半部分的答案的影响
比如前半部分对dp[mid+1]的响:
dp[l]*a[mid+1l]+dp[l+1]*a[mid-l]+…+dp[mid]*a[1]
dp[x]*a[y]中x+y之和始终为mid+1,这显然是一个裸的fft
每一层fft的总长度都为n,有logn层,所以总的时间复杂度为n*logn*logn

#include<bits/stdc++.h>
using namespace std;
const int MOD=313;
const double PI=acos(-1.0);

struct Complex{
    double x,y;
    Complex(double _x=0,double _y=0){
        x=_x;
        y=_y;
    }
    Complex operator -(const Complex &b)const{
        return Complex(x-b.x,y-b.y);
    }
    Complex operator +(const Complex &b)const{
        return Complex(x+b.x,y+b.y);
    }
    Complex operator *(const Complex &b)const{
        return Complex(x*b.x-y*b.y,x*b.y+y*b.x);
    }
};

void change(Complex y[],int len){
    int i,j,k;
    for(i=1,j=len/2;i<len-1;i++){
        if(i<j)
            swap(y[i],y[j]);
        k=len/2;
        while(j>=k){
            j-=k;
            k/=2;
        }
        if(j<k) j+=k;
    }
}

void fft(Complex y[],int len,int on){
    change(y,len);
    for(int h=2;h<=len;h<<=1){
        Complex wn(cos(-on*2*PI/h),sin(-on*2*PI/h));
        for(int j=0;j<len;j+=h){
            Complex w(1,0);
            for(int k=j;k<j+h/2;k++){
                Complex u=y[k];
                Complex t=w*y[k+h/2];
                y[k]=u+t;
                y[k+h/2]=u-t;
                w=w*wn;
            }
        }
    }
    if(on==-1)
        for(int i=0;i<len;i++)
            y[i].x/=len;
}
const int MAXN=200010;
Complex x1[MAXN],x2[MAXN];
int dp[MAXN],a[MAXN];

void solve(int l,int r){
    if(l==r){
        dp[l]=(dp[l]+a[l])%MOD;
        return ;
    }
    int mid=l+r>>1;
    solve(l,mid);
    int len=1,len1=mid-l+1,len2=r-l+1;
    while(len<len2)     len<<=1;
    for(int i=0;i<len1;i++)
        x1[i]=Complex(dp[i+l],0);
    for(int i=len1;i<len;i++)
        x1[i]=Complex(0,0);
    for(int i=0;i<len2;i++)
        x2[i]=Complex(a[i],0);
    for(int i=len2;i<len;i++)
        x2[i]=Complex(0,0);
    fft(x1,len,1),fft(x2,len,1);
    for(int i=0;i<len;i++)
        x1[i]=x1[i]*x2[i];
    fft(x1,len,-1);
    for(int i=mid+1;i<=r;i++)
        dp[i]=(dp[i]+(int)(x1[i-l].x+0.5))%MOD;
    solve(mid+1,r);
}

int main(){
    int n;
    while(scanf("%d",&n)!=EOF){
        if(n==0)
            break;
        for(int i=1;i<=n;i++)   dp[i]=0;
        for(int i=1;i<=n;i++)
            scanf("%d",&a[i]),a[i]%=MOD;
        solve(1,n);
//        for(int i=1;i<=n;i++)
//            printf("dp[%d] %d\n",i,dp[i]);
        printf("%d\n",dp[n]);
    }
    return 0;
}