3160: 万径人踪灭|FFT|manacher

答案可以转化为所有的回文子序列减去回文子串
回文子串的个数可以用$manacher$来求出
回文子序列的个数可以这样求：
先求出以每个点为中心左右对称的点的个数$x$，那么以这个点为中心的回文子序列的个数就是$2^x-1$，然后现在只需要求出以每个点为中心左右对称的点的个数，就是$a[i+p]==a[i-p]$这种情况的个数，发现和卷积非常类似。。。然后还只有$a,b$两个字母，可以分开求，令每个$a$的位置为$1$,$b$的位置为$0$,然后$FFT$求出卷积，对字母$b$的做法也同理，这样就可以直接求出以每个点为中心的左右对称的点的个数，同时也能求出长度为偶数的情况
手打的$complex$比stl自带的$complex$快了居然有两倍。。。

#include<algorithm>
#include<iostream>
#include<complex>
#include<cstdlib>
#include<cstring>
#include<cstdio>
#include<vector>
#include<queue>
#include<ctime>
#include<cmath>
#include<map>
#include<set>
#define N 266333
#define mod 1000000007
using namespace std;
struct E
{
    double real,i;
    E (double a=0,double b=0){real=a,i=b;}
    E operator+(E a){return (E){real+a.real,i+a.i};}
    E operator-(E a){return (E){real-a.real,i-a.i};}
    E operator*(E a){return (E){real*a.real-i*a.i,real*a.i+i*a.real};}
};
const double pi=acos(-1);
char s[N],str[N];
int n,m,L,p[N],f[N],R[N];
E a[N],b[N];
long long ans;
void manacher()
{
    int mx=0,id=0;
    for(int i=0;i<n;i++)
    {
        if(mx>i)p[i]=min(p[2*id-i],mx-i);else p[i]=1;
        while(i-p[i]>=0&&str[i+p[i]]==str[i-p[i]])p[i]++;
        if(p[i]+i>mx)mx=p[i]+i,id=i;
    }
}
void FFT(E *a,int f)
{
    for(int i=0;i<n;i++)if(i<R[i])swap(a[i],a[R[i]]);
    for(int i=1;i<n;i<<=1)
    {
        E wn(cos(pi/i),f*sin(pi/i));
        for(int k=0;k<n;k+=(i<<1))
        {
            E w(1,0);
            for(int j=0;j<i;j++,w=w*wn)
            {
                E x=a[k+j],y=w*a[k+j+i];
                a[k+j]=x+y,a[k+j+i]=x-y;
            }
        }
    }
    if(f==-1)for(int i=0;i<n;i++)a[i]=(E){a[i].real/n,a[i].i/n};
}
int cal(int y)
{
    long long ans=1,x=2;
    for(;y;x=x*x%mod,y>>=1)
        if(y&1)ans=ans*x%mod;
    return ans;
}
int main()
{
    scanf("%s",s);m=strlen(s);
    str[0]='#';
    for(int i=0;i<m;i++)
        str[2*i+1]=s[i],str[2*i+2]='#';
    n=2*m+1;manacher();

    for(int i=0;i<m;i++)ans-=(p[2*i+1]+1)/2;
    for(int i=0;i<m;i++)ans-=(p[2*i+2]-1)/2;
    m=m*2;
    for(n=2;n<m;n<<=1,L++);
    for(int i=1;i<n;i++)R[i]=(R[i>>1]>>1)|((i&1)<<L);

    for(int i=0;i<n;i++)a[i].real=(s[i]=='a'),a[i].i=0;
    FFT(a,1);
    for(int i=0;i<n;i++)a[i]=a[i]*a[i];
    FFT(a,-1);
    for(int i=0;i<n;i++)f[i]+=((int)(a[i].real+0.1)+(~i&1))/2;


    for(int i=0;i<n;i++)a[i].real=(s[i]=='b'),a[i].i=0;
    FFT(a,1);
    for(int i=0;i<n;i++)a[i]=a[i]*a[i];
    FFT(a,-1);
    for(int i=0;i<n;i++)f[i]+=((int)(a[i].real+0.1)+(~i&1))/2;

    for(int i=0;i<n;i++)
        if(f[i])ans+=cal(f[i])-1;
    cout<<(ans%mod+mod)%mod;
    return 0;
}