两道FFT的应用题

解题报告

---

BZOJ 快速傅里叶之二

题意

计算：

$$C[k]=\sum_{k\le i<n}(a[i]*b[i-k])$$

思路

正好看到《具体数学》上处理和式的Tricks，虽然热身题也不会做，但碾OI题还是很稳的(Orz神犇高教授)...

对于：

$$\sum_{k\le i<n}(a[i]*b[i-k])$$

将b数组倒置，即$b[i] = b'[n-1-i]$，原式变为：

$$c[k] = \sum_{0\le i<n, 0\le n+k-1-i<n}(a[i]*b'[n+k-1-i])$$

化简得到

$$c[k] = \sum_{k\le i\le n+k-1}(a[i]*b'[n+k-1-i])$$

我们想得到形如

$$c'[k] = \sum_{0\le i\le k}a[i]*b'[k-i]$$

的形式。只需要将后式中的$k$用$k+n-1$替换：

$$c'[k+n-1] = \sum_i a[i]*b'[k+n-1-i]$$

条件是

$$0\le k+n-1-i<n\iff k\le i \le n+k-1$$

所以

$$c[k] = c'[k+n-1] = \sum_{k\le i\le n+k-1} a[i]*b'[k+n-1-i]$$

然后用卷积做就好了。

Code

#include <bits/stdc++.h>
using namespace std;
typedef complex<double> Complex;
const int MAXN = 300005;
int rev[MAXN];
Complex A[MAXN];
Complex a[MAXN], b[MAXN], c[MAXN];
int n;
void fft(Complex a[], int n, int flag)
{
    int lgn = int(log2(n)+0.01);
    rev[0] = 0;
    for (int i = 1; i < n; i++)
        rev[i] = (rev[i>>1]>>1)|((i&1)<<(lgn-1));
    for (int i = 0; i < n; i++)
        A[rev[i]] = a[i];
    Complex u, t;
    for (int k = 2; k <= n; k <<= 1) {
        Complex dw = Complex(cos(2*M_PI/k), flag*sin(2*M_PI/k));
        for (int i = 0; i < n; i += k) {
            Complex w = 1;
            for (int j = 0; j < k>>1; j++) {
                u = A[i+j], t = w*A[i+j+(k>>1)];
                A[i+j] = u+t, A[i+j+(k>>1)] = u-t;
                w *= dw;
            }
        }
    }
    for (int i = 0; i < n; i++)
        a[i] = A[i]/Complex(flag==1?1:n,0);
}
int main()
{
    scanf("%d", &n);
    for (int i = 0; i < n; i++) {
        double x, y; scanf("%lf%lf", &x, &y);
        a[i] = Complex(x, 0), b[n-i-1] = Complex(y, 0);
    }
    int nn = 1;
    while (nn < n*2) nn <<= 1;
    fft(a, nn, 1), fft(b, nn, 1);
    for (int i = 0; i < nn; i++) c[i] = a[i]*b[i];
    fft(c, nn, -1);
    for (int i = 0; i < n; i++) printf("%d\n", int(c[i+n-1].real()+0.01));
    return 0;
}

BZOJ3160 万境人踪灭

题意

给定一个ab串，求所有不连续回文子序列的数量和。

思路

由于卷积可以处理关于一个对称轴两侧对称的字符总数，因此将串中的a、b变为0、1，自己和自己做卷积，就可以求出所有对称轴下b所对应数对称的字符总数。然后把a、b变为1、0，求出a的结果，相加即为关于对称轴两边对称的总字符数$A_i$，则$2^{A_i}$就是关于其对称的字符串总数。

那么如何去除连续的呢？跑一遍Manacher就好了。复杂度$O(n\lg n)$

Code

#include <bits/stdc++.h>
using namespace std;
const int MAXN = 300005;
typedef complex<double> Complex;
int rev[MAXN];
Complex A[MAXN];
const long long mod = 1000000007ll;
void fft(Complex a[], int n, int flag)
{
    rev[0] = 0;
    int lgn = int(log2(n)+0.01);
    for (int i = 1; i < n; i++)
        rev[i] = (rev[i>>1]>>1)|((i&1)<<(lgn-1));
    for (int i = 0; i < n; i++)
        A[rev[i]] = a[i];
    Complex u, t;
    for (int k = 2; k <= n; k <<= 1) {
        Complex dw = Complex(cos(2*M_PI/k), sin(flag*2*M_PI/k));
        for (int i = 0; i < n; i += k) {
            Complex w = 1;
            for (int j = 0; j < k>>1; j++) {
                u = A[i+j], t = A[i+j+(k>>1)]*w;
                A[i+j] = u+t, A[i+j+(k>>1)] = u-t;
                w *= dw;
            }
        }
    }
    for (int i = 0; i < n; i++)
        a[i] = A[i]/(flag == 1?1:Complex(n,0));
}
int p[MAXN];
long long manacher(char str[], int n)
{
    str[0] = '^', str[++n] = '#';
    int id = 0, mx = 0;
    long long ans = 0;
    for (int i = 1; i <= n; i++) {
        if (mx >= i) p[i] = min(mx-i, p[id-(i-id)]);
        else p[i] = 0;
        while (str[i+p[i]+1] == str[i-p[i]-1]) p[i]++;
        if (i+p[i] > mx) id = i, mx = i+p[i];
        (ans += p[i]+1) %= mod;
    }
    ans--;
    id = mx = 0;
    for (int i = 1; i < n; i++) {
        if (mx >= i) p[i] = min(mx-i, p[id-(i-id)]);
        else p[i] = 0;
        while (str[i+p[i]+1] == str[i-p[i]]) p[i]++;
        if (i+p[i] > mx) id = i, mx = i+p[i];
        (ans += p[i]) %= mod;
    }
    return ans;
}
long long power(int a, int n)
{
    if (n == 0) return 1;
    long long p = power(a, n>>1);
    (p *= p) %= mod;
    if (n&1) (p *= a) %= mod;
    return p;
}
char str[MAXN];
Complex a[MAXN], c[MAXN];
int cp[MAXN];
int main()
{
    scanf("%s", str+1);
    int len = strlen(str+1), n;
    for (n = 1; n < (len+1)*2; n<<=1);
    long long sub = manacher(str, len), ans = 0;
    for (int i = 1; i <= len; i++) a[i] = str[i] == 'a';
    fft(a, n, 1);
    for (int i = 0; i < n; i++) a[i] *= a[i];
    fft(a, n, -1);
    for (int i = 0; i < n; i++) cp[i] = int(a[i].real()+0.1), a[i] = 0;
    for (int i = 1; i <= len; i++) a[i] = str[i] == 'b';
    fft(a, n, 1);
    for (int i = 0; i < n; i++) a[i] *= a[i];
    fft(a, n, -1);
    for (int i = 0; i < n; i++) cp[i] += int(a[i].real()+0.1);
    for (int i = 0; i < n; i++) (ans += power(2, (cp[i]+1)/2)-1) %= mod;
    cout << ((ans-sub)%mod+mod)%mod << endl;
    return 0;
}