FFT_应用和例题

卷积

现有两个定义在 N 上的函数 (f(n),g(n))，定义 (f) 和 (g) 的卷积（convolution）为 (f otimes g)

[(f otimes g)(n) = sum_{i=0}^n f(i)g(n-i) ]

示意图: from http://blog.miskcoo.com/2015/04/polynomial-multiplication-and-fast-fourier-transform#i-17

考虑两多项式 (A, B) 的乘积 (C), (c(x) = sum_{i=0}^{x} a(i) cdot b(x - i))
系数记为卷积形式

于是计算卷积 ((f otimes g)(n)) 就可以把 (f, g) 的值直接作为系数写成两个多项式, 然后 FFT 计算多项式乘积, 得到的系数的前 (n) 项即为所求

BZOJ3527[ZJOI2014]力

题意:
给出 (n) 个数 (q_i) ，给出 (F_j) 的定义如下：

[F_j = sum_{i < j}frac{q_i q_j}{(i - j) ^ 2} - sum_{i > j}frac{q_i q_j}{(i - j) ^ 2} ]

令 (E_j = F_j / q_j) , 求 (E_j)

Sol:

因为知道是卷积的例题了, 所以想着把这个式子往卷积的方向靠

[egin{aligned} E_j = sum_{i < j}frac{q_i}{(i - j) ^ 2} - sum_{i > j}frac{q_i}{(i - j) ^ 2} end{aligned} ]

考虑分母当做系数, 再考虑下标和为 (j) , 写成这样

[egin{aligned} E_j = sum_{i < j}frac{1}{(j - i) ^ 2}q_i - sum_{i > j}frac{1}{(j - i) ^ 2}q_i end{aligned} ]

一开始想把两个一起做, 发现写不出两个函数, 于是考虑分开做
显然 (sum_{i < j}frac{1}{(j - i) ^ 2}q_i) 就是 (f(n) = q_n) 和 (g(n) = frac{1}{n^2}) 的卷积
然后后面一项同理, 把 (f(n)) 翻转一下即可

然后跑 FFT

double ans[MAXN], q[MAXN];

/*
20191212
0859~0922~0939
BZOJ3527 FFT
 */

int main()
{
    scanf("%d", &lena);
    for (int i = 0; i < lena; ++ i) 
    {
	scanf("%lf", &a[i].x); q[lena - 1 - i] = a[i].x;
	b[i].x = (i == 0 ? 0.0 : 1.0 / i / i);
    }
    while ((1 << dgt) < lena * 2) ++ dgt;
    n = 1 << dgt;
    init(n, dgt);
    FFT(b, n ,1);
    FFT(a, n, 1);
    for (int i = 0; i < n; ++ i) a[i] = a[i] * b[i];
    FFT(a, n, -1);
    for (int i = 0; i < lena; ++ i) ans[i] += a[i].x / n;
    for (int i = 0; i < n; ++ i) a[i].x = q[i], a[i].y = 0;
    FFT(a, n, 1);
    for (int i = 0; i < n; ++ i) a[i] = a[i] * b[i];
    FFT(a, n, -1);
    for (int i = 0; i < lena; ++ i) ans[i] -= a[lena - 1 - i].x / n;
    for (int i = 0; i < lena; ++ i) printf("%.3f
", ans[i]);
    return 0;
}
/*
3
1 2 3
 */

字符串匹配

给定两个字符串 (A, B, |A| ge |B|), 用 (B) 取匹配 (A),
那么可以发现对应位置的差恒定, 要转化成卷积形式, 可以将 (B) 翻转, 于是就可以构造卷积了

BZOJ4892[Tjoi2017]DNA

题意:
多测, 给两个字符串 (A, B), 字符集是 ACGT, 匹配的定义是相差不超过 3 个字符, 求 (B) 在 (A) 中匹配的次数
(n le 1e5, T le 10)

Sol:
翻转 (B)
一开始构造了这样的 : (ans_i = sum_{j=0}^i (a_j - b_{i - j}) ^ 2 cdot a_j cdot b_{i-j})
然后单独计算 (B) 中四个字符的贡献, 36 个 FFT

其实不需要这么套路, 反正都单独计算了, 可以更加钦点一点
(ans_i = sum_{j=0}^i a_j cdot b_{i-j})
(B) 中是枚举的字符就 1 否则 0 , (A) 中相反
这样得到的某一位卷积就是不一样的个数

int ans[MAXN];
void read(int * a, int & len)
{
    scanf("%s", tmp);
    len = strlen(tmp);
    for (int i = 0; i < len; ++ i) 
	if (tmp[i] == 'A') a[i] = 0;
	else if (tmp[i] == 'G') a[i] = 1;
	else if (tmp[i] == 'C') a[i] = 2;
	else a[i] = 3;
}
int fpow(int a, int b)
{
    int ret = 1;
    for (int i = 1; i <= b; ++ i) ret *= a;
    return ret;
}

void solve(int x)
{
    for (int i = 0; i < n; ++ i) ca[i] = cb[i] = 0;
    for (int i = 0; i < lena; ++ i)
    {
	if (sa[i] == x) ca[i] = 0;
	else ca[i] = 1;
    } 
    for (int i = 0; i < lenb; ++ i)
    {
	if (sb[i] == x) cb[i] = 1;
	else cb[i] = 0;
    } 
    for (int i = 0; i < n; ++ i) 
	a[i] = (Cpx) { ca[i], 0 }, b[i] = (Cpx) { cb[i], 0 };
    FFT(a, n, 1); FFT(b, n, 1);
    for (int i = 0; i < n; ++ i) a[i] = a[i] * b[i];
    FFT(a, n, -1);
    for (int i = 0; i < n; ++ i) 
	ans[i] += (int)(a[i].x / n + 0.5);
}

int main()
{
    int T;
    scanf("%d", &T);
    while (T --)
    {
	read(sa, lena); 
	read(sb, lenb);
	reverse(sb + 0, sb + lenb);
	n = 1, dgt = 0;
	while (n < lena + lenb) n <<= 1, ++ dgt;
	init(n, dgt);
	for (int i = 0; i < n; ++ i) ans[i] = 0;
	for (int i = 0; i < 4; ++ i) solve(i);
	int cnt = 0;
	for (int i = lenb - 1; i < lena; ++ i) 
	    if (ans[i] <= 3) ++ cnt;
	printf("%d
", cnt);
    }
    return 0;
}

/*
2
ATCGCCCTA
CTTCA
ACGT
GTCA
 */