【BZOJ3160】 万径人踪灭

万径人踪灭

题目大意

给一个全都是a和b的序列，求以某点为中心的对称的子序列的总数量（子序列不能完全连续），长度在FFT范围内

Solution

首先我们答案就是总共的对称子序列数量减去回文串的个数

回文串个数用马拉车或者PAM都可以解决

考虑怎么解决总共的对称子序列数量

我们发现，如果可以求出以某个点（可能是夹缝）为中心的对称的点对的数量，我们就可以求出满足某点为中心的对称子序列的数量，这个数量是(2^{组数} -1)，-1意思是不能都不取

然后这个过程可以多项式乘法模拟

我们将原序列对应成两个多项式A和B，A中的每一项满足：若序列这一位为a，那么这一项就是1，否则为0

B则为A取逆

A、B各自平方，我们发现它们就是以某点为中心，且字符为a或b的对称的组数（此时若次数为奇数，那么就是夹缝为中心，否则就是以某字符为中心）

所以把A、B系数加起来，像上面一样统计答案，就愉快地解决了这么个入门题

code:

#include<bits/stdc++.h>
using namespace std;
struct comp{
	double x,y;
	comp(double xx=0,double yy=0){
		x=xx,y=yy;
	}
}a[3000010],b[3000010];
comp operator +(comp a,comp b){return comp(a.x+b.x,a.y+b.y);}
comp operator -(comp a,comp b){return comp(a.x-b.x,a.y-b.y);}
comp operator *(comp a,comp b){return comp(a.x*b.x-a.y*b.y,a.x*b.y+b.x*a.y);}
int rev[3000010];
int lim=1;
const double pi=acos(-1);
void fft(comp *A,int type){
	for(int i=0;i<lim;++i){
		if(i<rev[i])swap(A[i],A[rev[i]]);
	}
	for(int mid=1;mid<lim;mid<<=1){
		comp wn=comp(cos(pi/mid),type*sin(pi/mid));
		for(int j=0,R=mid<<1;j<lim;j+=R){
			comp w=comp(1,0);
			for(int k=0;k<mid;++k,w=w*wn){
				comp x=A[j+k],y=w*A[j+k+mid];
				A[j+k]=x+y;
				A[j+k+mid]=x-y;
			}
		}
	}
	if(type==-1){
		for(int i=0;i<lim;++i){
			A[i].x/=lim;
		}
	}
}
const int mod=(1e9+7);
char s[3000010];
char as[3000010];
int p[3000010];
int len;
void manacher(){
	int mid=0,r=0;
	for(int i=1;i<=len;++i){
		if(i<=r){
			p[i]=min(p[mid*2-i],r-i+1);
		}
		else p[i]=1;
		while(as[i-p[i]]==as[i+p[i]]&&i-p[i]>=0&&i+p[i]<=len)p[i]++;
		if(i+p[i]>r){
			r=i+p[i]-1;
			mid=i;
		}
	}
}
int num[3000010];
int ans=0;
int poww[3000010];
int main(){
	scanf("%s",s+1);
	len=strlen(s+1);
	int tmp=0;
	while(lim<=len*2){
		lim<<=1;
		tmp++; 
	}
	for(int i=0;i<=lim;++i){
		rev[i]=(rev[i>>1]>>1)|((1&i)<<(tmp-1));
	}
	as[0]=2;
	for(int i=1;i<=len;++i){
		as[(i-1)*2+1]=s[i]=='a';
		as[i*2]=2;
	}
	len<<=1;
	manacher();
	for(int i=1;i<=len;++i){
		ans=(ans-(p[i]/2)+mod)%mod;
	}
	for(int i=0;i<len/2;++i){
		a[i].x=b[i].x=(s[i+1]=='a');
	}
	fft(a,1);
	fft(b,1);
	for(int i=0;i<lim;++i){
		a[i]=a[i]*b[i];
	}
	fft(a,-1);
	for(int i=0;i<lim;++i){
		num[i]=num[i]+int(a[i].x+0.5);
	}
	for(int i=0;i<lim;++i){
		a[i]=comp(0,0);
		b[i]=comp(0,0);
	}
	for(int i=0;i<len/2;++i){
		a[i].x=b[i].x=(s[i+1]=='b');
	}
	fft(a,1);
	fft(b,1);
	for(int i=0;i<lim;++i){
		a[i]=a[i]*b[i];
	}
	fft(a,-1);
	for(int i=0;i<lim;++i){
		num[i]+=int(a[i].x+0.5);
	}
	poww[0]=1;
	for(int i=1;i<lim;++i){
		poww[i]=(poww[i-1]<<1)%mod;
	}
	for(int i=0;i<lim;++i){
		ans=(ans+poww[(num[i]+1)/2]-1)%mod;
	}
	printf("%d
",ans);
}