【bzoj3160】万径人踪灭 Manacher+FFT

题目描述

 

输入

 

---

题解

Manacher+FFT

显然答案=总数-不合法的方案数。

本题中不合法的方案数即连续的回文子序列，即回文子串。可以使用Manacher来求出该串的回文子串的个数（Manacher都忘得差不多了）

然后求总数，即求对于每个对称轴的回文子序列的个数，我们可以先求出回文子序列的最长长度。

可以想到对于两个字符，相同为1，不同为0；而这又是一个01串，所以可以先把a当作1，使用FFT求卷积；再把b当作1，再求卷积，即可得到回文子序列的最长长度。（FFT也快忘得差不多了）

而对于最长长度为l的回文子序列，它包含了$lceilfrac l2 ceil$个关于该对称轴对称的字符，而每个字符都是可选可不选的，所以对总数的贡献为$2^{lceilfrac l2 ceil}-1$（因为不能一个也不选）。

最后二者相减即为答案。

#include <cmath>
#include <cstdio>
#include <cstring>
#include <algorithm>
#define N 1000010
using namespace std;
const int mod = 1000000007;
const double pi = acos(-1);
struct data
{
	double x , y;
	data() {}
	data(double x0 , double y0) {x = x0 , y = y0;}
	data operator+(const data a)const {return data(x + a.x , y + a.y);}
	data operator-(const data a)const {return data(x - a.x , y - a.y);}
	data operator*(const data a)const {return data(x * a.x - y * a.y , x * a.y + y * a.x);}
}A[N << 2] , B[N << 2];
char str[N] , tmp[N << 1];
int p[N << 1] , num[N << 1] , po[N];
void fft(data *a , int len , int flag)
{
	int i , j , k;
	for(i = k = 0 ; i < len ; i ++ )
	{
		if(i > k) swap(a[i] , a[k]);
		for(j = len >> 1 ; (k ^= j) < j ; j >>= 1);
	}
	for(k = 2 ; k <= len ; k <<= 1)
	{
		data wn(cos(2 * pi * flag / k) , sin(2 * pi * flag / k));
		for(i = 0 ; i < len ; i += k)
		{
			data t , w(1 , 0);
			for(j = i ; j < i + (k >> 1) ; j ++ , w = w * wn)
				t = w * a[j + (k >> 1)] , a[j + (k >> 1)] = a[j] - t , a[j] = a[j] + t;
		}
	}
}
void work(data *a , data *b , int len)
{
	int i;
	fft(a , len , 1) , fft(b , len , 1);
	for(i = 0 ; i < len ; i ++ ) a[i] = a[i] * b[i];
	fft(a , len , -1);
	for(i = 0 ; i < len ; i ++ ) a[i].x /= len;
}
int main()
{
	int n , i , last , mx = 0 , len , ans = 0;
	scanf("%s" , str) , n = strlen(str);
	tmp[0] = '$';
	for(i = 0 ; i < n ; i ++ ) tmp[i * 2 + 2] = str[i];
	for(i = 1 ; i <= n * 2 ; i ++ )
	{
		if(mx >= i) p[i] = min(p[2 * last - i] , mx - i + 1);
		else p[i] = 1;
		while(tmp[i - p[i]] == tmp[i + p[i]]) p[i] ++ ;
		if(i + p[i] - 1 > mx) mx = i + p[i] - 1 , last = i;
	}
	for(i = 1 ; i <= n * 2 ; i ++ ) ans = (ans - p[i] / 2 + mod) % mod;
	for(len = 1 ; len < n * 2 ; len <<= 1);
	for(i = 0 ; i < n ; i ++ )
		if(str[i] == 'a')
			A[i].x = B[i].x = 1;
	work(A , B , len);
	for(i = 0 ; i < len ; i ++ ) num[i] += (int)(A[i].x + 0.1);
	for(i = 0 ; i < len ; i ++ ) A[i] = B[i] = data(0 , 0);
	for(i = 0 ; i < n ; i ++ )
		if(str[i] == 'b')
			A[i].x = B[i].x = 1;
	work(A , B , len);
	for(i = 0 ; i < len ; i ++ ) num[i] += (int)(A[i].x + 0.1);
	for(i = po[0] = 1 ; i < len ; i ++ ) po[i] = (po[i - 1] << 1) % mod;
	for(i = 0 ; i < len ; i ++ ) ans = (ans + po[(num[i] + 1) >> 1] - 1) % mod;
	printf("%d
" , ans);
	return 0;
}