FFT+manacher--bzoj3160: 万径人踪灭

传送门

一道隐蔽的 F F T FFT FFT

题目要求回文子序列的个数，不能连续的话就用全部的减去连续的，连续的可以用 m a n a c h e r manacher manacher计算，全部的回文子序列怎么求呢？

因为这个序列中只有 a a a和 b b b，假设 a = 1 , b = 0 a=1,b=0 a=1,b=0，那么两个位置的字符都是 a a a就是相乘为 0 0 0（ b b b的反着来同理）。

这样的话以 i i i为中心的回文子序列就可以看成 s [ i − k ] × s [ i + k ] = 1 s[i-k]\times s[i+k]=1 s[i−k]×s[i+k]=1这样的形式，两边的下标相加为定值，就可以使用 F F T FFT FFT来实现。

也就是传说中的“我卷我自己”。

算答案的时候 i i i位置的值就是 i 2 \frac{i}{2} 2i​位置为中心的对称相等的个数，假设为 x x x，那么这个位置答案就是 2 x + 1 2 − 1 2^{\frac{x+1}{2}}-1 22x+1​−1，所有位置相加再减去 m a n a c h e r manacher manacher求出的回文子串个数就是最后的答案

#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cstring>
#include<cmath>
#define maxn 400005
#define int long long
using namespace std;
const int mod=1e9+7;
const double Pi=acos(-1.0);

struct complex{
	double x,y;
	complex(double xx=0,double yy=0){x=xx,y=yy;}
}a[maxn],b[maxn];
complex operator +(complex a,complex b){return complex(a.x+b.x,a.y+b.y);}
complex operator -(complex a,complex b){return complex(a.x-b.x,a.y-b.y);}
complex operator *(complex a,complex b){return complex(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);}

int n,limit=1,l,rev[maxn],bin[maxn],len[maxn],ans;
char s[maxn],tmp[maxn];

void FFT(complex *F,int type){
	for(int i=0;i<limit;i++)
		if(i<rev[i]) swap(F[i],F[rev[i]]);
	for(int mid=1;mid<limit;mid<<=1){
		complex Wn(cos(Pi/mid),1.0*type*sin(Pi/mid));
		for(int r=mid<<1,j=0;j<limit;j+=r){
			complex w(1,0);
			for(int k=0;k<mid;k++,w=w*Wn){
				complex x=F[j+k],y=w*F[j+mid+k];
				F[j+k]=x+y,F[j+mid+k]=x-y;
			}
		}
	}
}

inline void prework(){
	for(int i=0;i<n;i++)
		if(s[i]=='a') a[i].x=1,b[i].x=0;
		else b[i].x=1,a[i].x=0;
	while(limit<2*n) limit<<=1,l++;
	for(int i=0;i<limit;i++)
		rev[i]=(rev[i>>1]>>1)|((i&1)<<(l-1));
	bin[0]=1;
	for(int i=1;i<=n;i++) bin[i]=(bin[i-1]<<1)%mod;
}

inline int init(char *p){
	tmp[0]='@';
	for(int i=1;i<=2*n;i+=2)
		tmp[i]='#',tmp[i+1]=p[i/2];
	tmp[2*n+1]='#'; tmp[2*n+2]='$';
	return 2*n+1;
}

inline int manacher(char *p,int m){
	int mx=0,pos=0,tot=0;
	for(int i=1;i<=m;i++){
		if(mx>i) len[i]=min(mx-i,len[2*pos-i]);
		else len[i]=1;
		while(p[i-len[i]]==p[i+len[i]]) len[i]++;
		if(len[i]+i>mx) mx=len[i]+i,pos=i;
		tot+=(len[i]/2)%mod;
	} return tot;
}

signed main(){
	scanf("%s",s); n=strlen(s);
	prework();
	FFT(a,1); FFT(b,1);
	for(int i=0;i<=limit;i++) a[i]=a[i]*a[i],a[i]=a[i]+b[i]*b[i];
	FFT(a,-1);
	for(int i=0;i<=2*n-2;i++){
		int xx=(int)((a[i].x+0.5)/limit);
		xx=(xx+1)>>1;
		(ans+=bin[xx]-1)%=mod;
	}
	int m=init(s); ans-=manacher(tmp,m)%mod;
	printf("%d\n",(ans+mod)%mod);
	return 0;
}