1. 程式人生 > >Luogu P4199/BZOJ3160 [2013湖北互測week1]萬徑人蹤滅

Luogu P4199/BZOJ3160 [2013湖北互測week1]萬徑人蹤滅

題目大意

給定一個長度為 n n 的字串 S S ,求有多少個非子串的子序列滿足所有字元和在原串中的位置都關於某條對稱軸對稱。
資料範圍  1

n 1 0 5 1\leqslant n\leqslant 10^5

題解

可以很容易發現,答案等於迴文子序列(位置也迴文)的個數減去迴文子串的個數。迴文子串的個數很容易用 m

a n a c h e r manacher 演算法或者 PAM
\text{PAM}
搞定,關鍵是前面。
為了相容 m a n a c h e r manacher 演算法,並且敘述簡便,這裡定義一個長度為 l = 2 n 2 l=2n-2 的字串 P P 表示中間用特殊字元 # \# 隔開的字串,則令 t i = j = 0 m i n ( j i , i ) [ P i j = P i + j ] t_i=\sum\limits_{j=0}^{min(j-i,i)} [P_{i-j}=P_{i+j}]
i i 為偶數則答案就是 2 t i + 1 1 2^{t_i+1}-1 (包括中心在內每對位置選或不選)。
i i 為奇數則答案就是 2 t i 1 2^{t_i}-1 (每對位置選或不選)。
t i t_i 的計算可以用 FFT 加速(同基站)。

程式碼

// luogu-judger-enable-o2
#include<bits/stdc++.h>
using namespace std;
#define maxn 1000100
#define maft 1310730
const double pi=acos(-1.0);
long long p[maxn*2];
void manacher(const char *_s){
    static char tmp[maxn*2];
    long long max_right=0,mid=0,cnt=-1;
    for(int i=0;_s[i];i++)
        tmp[++cnt]=_s[i],tmp[++cnt]='#';
    tmp[cnt]='@';
    for(int i=0;i<cnt;i++){
        if(i<max_right)
            p[i]=min(max_right-i,p[mid+mid-i]);
        else p[i]=1;
        for(;tmp[i-p[i]]==tmp[i+p[i]];p[i]++);
        if(p[i]+i>max_right){
            max_right=p[i]+i;
            mid=i;
        }
    } for(int i=1;i<=cnt;i++)
        p[i]=(p[i]+(i%2==0))>>1;
}
struct comp{
    double x,y;
    comp(double a=0,double b=0):x(a),y(b){}
    comp operator+(const comp &t){return comp(x+t.x,y+t.y);}
    comp operator-(const comp &t){return comp(x-t.x,y-t.y);}
    comp operator*(const comp &t){return comp(x*t.x-y*t.y,y*t.x+x*t.y);}
};
long long limit,l,r[maft];
void reset(long long n,long long m){
    for(limit=1,l=0;limit<=n+m;limit<<=1,l++);
    for(int i=1;i<limit;i++)
        r[i]=((r[i>>1]>>1)|((i&1)<<(l-1)));
}
void fft(comp *t,long long ty){
    for(int i=1;i<limit;i++)
        if(i<r[i]) swap(t[i],t[r[i]]);
    for(long long mid=1;mid<limit;mid<<=1){
        comp wn(cos(pi/mid),ty*sin(pi/mid));
        for(int j=0,R=(mid<<1);j<limit;j+=R){
            comp w(1,0);
            for(int k=0;k<mid;k++,w=w*wn){
                comp x=t[j+k],y=w*t[j+k+mid];
                t[j+k]=x+y,t[j+k+mid]=x-y;
            }
        }
    }
}
char s[maxn];
comp A[maft],B[maft];
long long sa[maft],sb[maft],S[maft];
long long ans,pows[maxn];
const long long mod=1e9+7;
int main(){
    scanf("%s",s);
    pows[0]=1;
    for(int i=1;i<maxn;i++)
        pows[i]=pows[i-1]*2%mod;
    manacher(s);
    const long long l=strlen(s);
    long long n=2*l-2;
    reset(l,l);
    for(int i=0;i<l;i++)
        if(s[i]=='a') A[i].x=1;
        else if(s[i]=='b') B[i].x=1;
    fft(A,1);fft(B,1);
    for(int i=0;i<limit;i++)
        A[i]=A[i]*A[i],B[i]=B[i]*B[i];
    fft(A,-1);fft(B,-1);
    for(int i=0;i<=n;i++){
        sa[i]=(long long)floor(A[i].x/limit+0.1);
        sb[i]=(long long)floor(B[i].x/limit+0.1);
        if((i&1)==0){
            if(s[i>>1]=='a')
                --sa[i];
            else if(s[i>>1]=='b')
                --sb[i];
        } sa[i]>>=1,sb[i]>>=1;
        ans+=pows[sa[i]+sb[i]+(!(i&1))]-p[i]-1+mod;
        ans%=mod;
    } printf("%lld",ans);
    return 0;
}