HDU-5730 Shell Necklace (CDQ+FFT/NTT)

HDU-5730(CDQ+FFT/NTT)

题意：将长度为\(n\)的序列分成若干段，每段\([l,r]\)的权值为\(a_{r-l+1}\)，一种分法的权值为所有段的乘积，求所有可能的分法的权值和

根据题意可以得到简单\(dp\)

\(dp_0=1,dp_i=\sum_0^{i-1}dp_j \cdot a_{i-j}\)

可以看到是一个\(i-j\)形式的作差卷积

但是直接卷积我们无法保证先求出了\(dp_j\)，所有可以用\(CDQ\)分治优化，复杂度\(n\log^2 n\)

具体实现见代码

#include<bits/stdc++.h>
using namespace std;

#define reg register
typedef long long ll;
#define rep(i,a,b) for(reg int i=a,i##end=b;i<=i##end;++i)
#define drep(i,a,b) for(reg int i=a,i##end=b;i>=i##end;--i)

template <class T> inline void cmin(T &a,T b){ ((a>b)&&(a=b)); } 
template <class T> inline void cmax(T &a,T b){ ((a<b)&&(a=b)); } 

char IO;
int rd(){
	int s=0,f=0;
	while(!isdigit(IO=getchar())) if(IO=='-') f=1;
	do s=(s<<1)+(s<<3)+(IO^'0');
	while(isdigit(IO=getchar()));
	return f?-s:s;
}

const double PI=acos(-1);

const int N=(1<<18)+4,P=313;

int n;
int a[N],dp[N];
struct Cp{
	double x,y;
	Cp operator + (const Cp t) const { return (Cp){x+t.x,y+t.y}; }
	Cp operator - (const Cp t) const { return (Cp){x-t.x,y-t.y}; }
	Cp operator * (const Cp t) const { return (Cp){x*t.x-y*t.y,x*t.y+y*t.x}; }
} A[N],B[N];

Cp w[N];
int rev[N];
void FFT(int n,Cp *a,Cp *w){
	rep(i,1,n-1) if(rev[i]>i) swap(a[i],a[rev[i]]);
	for(reg int i=1;i<n;i<<=1) {
		int len=i*2;
		for(reg int l=0;l<n;l+=len) {
			for(reg int j=l;j<l+i;++j) {
				Cp t=a[j+i]*w[n/len*(j-l)];
				a[j+i]=a[j]-t;
				a[j]=a[j]+t;
			}
		}
	}
}

void Solve(int l,int r){ 
	if(l==r) return;
	int mid=(l+r)>>1;
	Solve(l,mid); //先求出左边所有的
    
	int R=1,c=-1;
	while(R<=r-l+1) R<<=1,c++;
	rep(i,1,R) rev[i]=(rev[i>>1]>>1)|((i&1)<<c);
	w[0]=(Cp){1,0};
	w[1]=(Cp){cos(2*PI/R),sin(2*PI/R)};
	rep(i,2,R) w[i]=w[i-1]*w[1];
	rep(i,l,mid) A[i-l]=(Cp){1.0*dp[i],0};
	rep(i,1,r-l+1) B[i]=(Cp){1.0*a[i],0};
	FFT(R,A,w),FFT(R,B,w);
	rep(i,0,R) w[i].y=-w[i].y,A[i]=A[i]*B[i];
	FFT(R,A,w);
	rep(i,mid+1,r) dp[i]=(dp[i]+ll(A[i-l].x/R+0.5))%P;
	rep(i,0,R) A[i].x=A[i].y=B[i].x=B[i].y=0;
    //让左边已经求出来的向右边转移
    
	Solve(mid+1,r);
    //让右边之间转移
}

int main(){
	while(~scanf("%d",&n) && n) {
		rep(i,1,n) a[i]=rd()%P,dp[i]=0;
		dp[0]=1;
		Solve(0,n);
		printf("%d\n",dp[n]);
	}
}