Processing math: 100%

//好吧标题是扯蛋的

本文全文 9624 个字(不包含这句声明),全都是 XZY 一个字一个字码出来的 QAQ
如果需要转载非常欢迎,但请注明出处:集合卷积——从 FMT 到 FWT —— MiNa!
谢谢!= ̄ω ̄=

1. 啥是集合卷积

先看看多项式乘法卷积:

Ck=i+j=kAi×Bj

这玩意可以用 FFT 做

集合卷积,我们先看或卷积:

Ck=i|j=kAi×Bj

就像这个一样,我们把 FFT 做的那个卷积里的 “乘” 改成 “与”、“或”、“异或” 等等,这种关于集合的逻辑运算的卷积就是集合卷积。

2. 或卷积

先把式子列出来:

Ck=i|j=kAi×Bj

暴力做法就是枚举 i,j,判断 i|j是否等于 k,总复杂度 O(n3)

那怎么办捏.jpg

当然就是用 FMT 啦(快速莫比乌斯变换)(这货和莫比乌斯变换貌似没有任何关系)

我们类比 FFT 的那种转换成点值表达,然后 O(n)每项乘到一起的思想,构造一波~

在 FFT 里,对于一个 n项多项式,我们会把它转成 n个点来表示它

在 FMT 里,对于一个 n项多项式 F,我们这样转换:

Fi=jiFj

这里的 ji表示 j的二进制表示中的每个 1都在 i的二进制表示中出现。

j&i=j

为啥这样构造呢?

Ck=i|j=kAi×Bj

Ai=jiAj

Bi=jiBj

Ci=jiCj

Ai×Bi=miAm×niBn=(m|n)iAm×Bn=Ci

嘿嘿嘿,真开心,也就是 Ai×Bi=Ci,这个过程是 O(n)哒!

那么问题来了:如何计算 Fi=jiFj

我们可以枚举子集嘿嘿嘿。。。

QwQ

好吧我们还是有 O(nlog2n)的做法的,依然是——分治

(TMD 刚刚不小心按到刷新键了刚写的全没了沃日。。。)

对于长度为 N的多项式 FN2的整数次幂,和 FFT 一个套路,不够次数的项补零),设其前 N2项为多项式 L,后 N2项为多项式 RLR的次数都是 N1

假设我们已经求得了 LR,怎么求得 F呢?

Li=jiFj(i<N2)

Ri=jiFj+N2(i<N2)

可以发现 RiF中的下标等于 LiF中的下标加 N2

也就是 LiF中的下标的二进制的第 log2(N)1位一定为 0Ri的一定为 1

Li所包含的集合,Ri也应当包含,因此合并的时候,Ri应当加上 Li

所以就有:

F={L},{R+L}

这里的 “,” 表示相接,比如 {1,2,3},{4,5,6}=1,2,3,4,5,6,“+” 表示按位相加(Ri+Li

有了正变换,我们考虑怎么从 F推到 F,也就是逆变换

(网上很多逆变换说得云里雾里看不懂什么鬼。。。)

通过容斥我们能找到这个:

Fi=jiFj×(1)bitcount(ij)

这里的 bitcount(a)表示 a在二进制下的 “1” 的个数

显然这个式子和正向 FMT 是一模一样的,只不过多了一个 1这个常数项

假设对于当前的 F,已经知道了 LR

我们前面已经说过了,RiF中的下标在二进制下只比 Li的多了一个 log2(N)1位的 1

也就是说:LR异号。。。

所以根据上面的:F={L},{R+L}

可以得出:F={L},{RL}

(就是原本的 Ri+Li变成了 Ri+(Li)

这样我们就能写出 FMT 的或卷积代码啦!

void FMT_OR(int (&a)[NS], int N, int f)
{
    for (int l = 1; l < N; l <<= 1)
        for (int i = 0; i < N; i += (l << 1))
            for (int j = i; j < i + l; j += 1)
                if (f == 1) a[j + l] = pls(a[j + l], a[j]);
                else a[j + l] = mns(a[j + l], a[j]);
}
C++

正向 FMT 就是 FMT_OR(a, N, 1),逆向的话(根据习惯)就是 FMT_OR(a, N, -1)

啊,多美妙啊!

3. 与卷积

与和或实质是一样的 —— 沃-兹基朔德

其实吧。。。

a & b = !((!a) | (!b))

因此与卷积只要把或卷积的下标 01 翻转就行了(0 变 1,1 变 0)。

所以构造出来的 FMT 式子长这样:

Fi=jiFj

也就是把 ji的子集变成了 ji的超集 OvO

然后捏

然后就是在合并的时候,LR的二进制少一个 1,因此 L包含 R。。。式子变成了:

F={L+R},{L}

emmmm… 逆向变换就是:

F={LR},{L}

好背吧

代码长这样:

void FMT_AND(int (&a)[NS], int N, int f)
{
    for (int l = 1; l < N; l <<= 1)
        for (int i = 0; i < N; i += (l << 1))
            for (int j = i; j < i + l; j += 1)
                if (f == 1) a[j] = pls(a[j], a[j + l]);
                else a[j] = mns(a[j], a[j + l]);
}
C++

4. 异或卷积

好吧我也不会证。。。

听 boshi 大佬说这玩意只能构造出来证明他是正确的,没法推导

(或者说是他也不会证)

异或就是对称差卷积。。。用的是 FWT

(然而实际上和 FMT 没啥区别,只不过发明的人不同而已。。。)

(而且长得都一样,只不过把中间的字母倒过来了)

把式子写一下吧:

正向:

F={L+R},{LR}

逆向:

F={L+R2},{LR2}

实际上根据毕姥爷的说法,我们对一个多项式 F做两遍 FWT,得到的是 F×N,也就是每个系数变大了多项式项数倍。。。

这个东西对模数比较奇怪的情况下还是有用的。。。

代码长这样:

void FWT_XOR(int (&a)[NS], int N, int f)
{
    for (int l = 1; l < N; l <<= 1)
        for (int i = 0; i < N; i += (l << 1))
            for (int j = i, t1, t2; j < i + l; j += 1)
            {
                t1 = a[j], t2 = a[j + l];
                a[j] = pls(t1, t2), a[j + l] = mns(t1, t2);
                if (f == -1) a[j] = mul(a[j], IV2), a[j + l] = mul(a[j + l], IV2);
            }
}
C++

5. 子集卷积

啥是子集卷积?

就是:Ck=ij=k&&ij=Ai×Bj

看起来就是或卷积的增强版,条件增加了卷积的俩货不能有交集

那怎么办捏.jpg×2

对于 Fi,我们增加一维使得它表示为 Fbitcount(i),i

bitcount(a)依然是上文提到的数二进制下 1的个数)

这样有什么好处呢?

我们增加了这一维,那么卷积形式就成了:

Cbitcount(k),k=ij=k&&bitcount(i)+bitcount(j)=bitcount(k)Abitcount(i),i×Bbitcount(j),j

如果我们设 Fx,i=0(bitcount(i)x)

那么上面那个式子也可以写成:

Cbitcount(k),k=ij=kAbitcount(i),i×B(bitcount(k)bitcount(i)),j

那很好,上面那个式子就是最普通的或卷积了

我们只需要枚举 bitcount(k)bitcount(i),再做或卷积即可。

但是这样卷积以后我们并不能保证卷出来的多项式 Cbitcount(k),x中的 x符合 bitcount(x)=bitcount(k)。。。

因此枚举一遍 x,如果 bitcount(x)bitcount(k),则将 Cbitcount(k),x置为 0

总复杂度是 O(n22n)的 QAQ

代码见例题 4 QwQ

6. 例题

1. Hard Nim BZOJ – 4589

毕姥爷讲课题

O(n)筛一波素数,然后对于素数搞个生成函数 G,若 p为素数,则 Gp=1,否则 Gp=0

接着快速幂求 Gn即可

答案为 Gn0

代码:

#include <bits/stdc++.h>

#define P (1000000007)
#define MS (150005)
#define IV2 (500000004)

using namespace std;

int pls(int a, int b) {return a + b >= P ? a + b - P : a + b;}
int mns(int a, int b) {return a - b < 0 ? a - b + P : a - b;}
int mul(int a, int b) {return 1ll * a * b % P;}

int n, m, prm[MS], pnt, p[MS], res[MS];

bool ntp[MS];

void FWT(int (&a)[MS], int N, int f)
{
    for (int l = 1; l < N; l <<= 1)
        for (int i = 0; i < N; i += (l << 1))
            for (int j = i, t1, t2; j < i + l; j += 1)
            {
                t1 = a[j], t2 = a[j + l];
                a[j] = pls(t1, t2), a[j + l] = mns(t1, t2);
                if (f == -1) a[j] = mul(a[j], IV2), a[j + l] = mul(a[j + l], IV2);
            }
}

int main(int argc, char const* argv[])
{
    for (int i = 2; i <= MS - 5; i += 1)
    {
        if (!ntp[i]) prm[++pnt] = i;
        for (int j = 1, k; j <= pnt; j += 1)
        {
            k = i * prm[j];
            if (k > MS - 5) break;
            ntp[k] = 1;
            if (i % prm[j] == 0) break;
        }
    }
    while (~scanf("%d%d", &n, &m))
    {
        memset(p, 0, sizeof(p));
        for (int i = 2; i <= m; i += 1) p[i] = !ntp[i];
        int N = 1; while (N <= m) N <<= 1;
        FWT(p, N, 1);
        for (int i = 0; i < N; i += 1) res[i] = 1;
        while (n)
        {
            if (n & 1) for (int i = 0; i < N; i += 1) res[i] = mul(res[i], p[i]);
            for (int i = 0; i < N; i += 1) p[i] = mul(p[i], p[i]);
            n >>= 1;
        }
        FWT(res, N, -1), printf("%d\n", res[0]);
    }
    return 0;
}
C++

2.【模板】快速沃尔什变换 LUOGU – 4717

与或异或卷积模板题 OvO

代码:

#include <bits/stdc++.h>

#define P (998244353)
#define NS (200000)
#define IV2 (499122177)

using namespace std;

template<typename _Tp> inline void IN(_Tp& dig)
{
    char c; bool flag = 0; dig = 0;
    while (c = getchar(), !isdigit(c)) if (c == '-') flag = 1;
    while (isdigit(c)) dig = dig * 10 + c - '0', c = getchar();
    if (flag) dig = -dig;
}

int pls(int a, int b) {return a + b >= P ? a + b - P : a + b;}
int mns(int a, int b) {return a - b < 0 ? a - b + P : a - b;}
int mul(int a, int b) {return 1ll * a * b % P;}

void FMT_OR(int (&a)[NS], int N, int f)
{
    for (int l = 1; l < N; l <<= 1)
        for (int i = 0; i < N; i += (l << 1))
            for (int j = i; j < i + l; j += 1)
                if (f == 1) a[j + l] = pls(a[j + l], a[j]);
                else a[j + l] = mns(a[j + l], a[j]);
}

void FMT_AND(int (&a)[NS], int N, int f)
{
    for (int l = 1; l < N; l <<= 1)
        for (int i = 0; i < N; i += (l << 1))
            for (int j = i; j < i + l; j += 1)
                if (f == 1) a[j] = pls(a[j], a[j + l]);
                else a[j] = mns(a[j], a[j + l]);
}

void FWT_XOR(int (&a)[NS], int N, int f)
{
    for (int l = 1; l < N; l <<= 1)
        for (int i = 0; i < N; i += (l << 1))
            for (int j = i, t1, t2; j < i + l; j += 1)
            {
                t1 = a[j], t2 = a[j + l];
                a[j] = pls(t1, t2), a[j + l] = mns(t1, t2);
                if (f == -1) a[j] = mul(a[j], IV2), a[j + l] = mul(a[j + l], IV2);
            }
}

int N, bs, A[NS], B[NS], C[NS];

int main(int argc, char const* argv[])
{
    IN(bs), N = (1 << bs);
    for (int i = 0; i < N; i += 1) IN(A[i]);
    for (int i = 0; i < N; i += 1) IN(B[i]);
    // or
    FMT_OR(A, N, 1), FMT_OR(B, N, 1);
    for (int i = 0; i < N; i += 1) C[i] = mul(A[i], B[i]);
    FMT_OR(C, N, -1);
    for (int i = 0; i < N; i += 1) printf("%d ", C[i]);
    putchar(10), FMT_OR(A, N, -1), FMT_OR(B, N, -1);
    // and
    FMT_AND(A, N, 1), FMT_AND(B, N, 1);
    for (int i = 0; i < N; i += 1) C[i] = mul(A[i], B[i]);
    FMT_AND(C, N, -1);
    for (int i = 0; i < N; i += 1) printf("%d ", C[i]);
    putchar(10), FMT_AND(A, N, -1), FMT_AND(B, N, -1);
    // xor
    FWT_XOR(A, N, 1), FWT_XOR(B, N, 1);
    for (int i = 0; i < N; i += 1) C[i] = mul(A[i], B[i]);
    FWT_XOR(C, N, -1);
    for (int i = 0; i < N; i += 1) printf("%d ", C[i]);
    return 0;
}
C++

3. [HAOI2015] 按位或 BZOJ – 4036

实际上我并不是很知道这题我是怎么 AC 的

我反正是这样想的:

Ans=i=1i×(PiPi1)

(注意这里只对询问取到 2n1有效,因为如果在第 i次之前就得到了 2n1,那么后面无论选任何数,得到的都依然是 2n1,也就是如果在第 i次前得到了 2n1,那么到第 i次,这种情况下得到 2n1的概率是 1。而对于别的数字则不一定。)

即:

Ans=P+2P22P+3P33P2+PP1

=PP2P3P1+P

因为 P这玩意在 2n1处是向 1无限接近的(除非原问题无解)

所以我们可以从 P里提取出 1出来,P就变成 1

那么原式就是:

Ans=(1P)+(1P2)+(1P3)++(1P1)+1

然后 1P这玩意就向 0 收敛了

真是太棒了,我们就能写出式子了!搞个级数求和:

Ans=1PP1

恩,然后发现确实能 AC。。。真开心。。。

(最后说一下其实这个貌似用 min-max 容斥会好想很多。。。)

代码:

#include <bits/stdc++.h>

#define NS (1 << 20)
#define lowbit(a) (a & -a)

using namespace std;

template<typename _Tp> inline void IN(_Tp& dig)
{
    char c; bool flag = 0; dig = 0;
    while (c = getchar(), !isdigit(c)) if (c == '-') flag = 1;
    while (isdigit(c)) dig = dig * 10 + c - '0', c = getchar();
    if (flag) dig = -dig;
}

void FMT(double (&a)[NS], int N, int f)
{
    for (int l = 1; l < N; l <<= 1)
        for (int i = 0; i < N; i += (l << 1))
            for (int j = i; j < i + l; j += 1)
                a[j + l] = a[j + l] + f * a[j];
}

int N, bs, rem;

double p[NS];

int main(int argc, char const* argv[])
{
    IN(bs), N = (1 << bs);
    for (int i = 0; i < N; i += 1)
    {
        scanf("%lf", &p[i]);
        if (p[i] > 0)
        {
            int j = i;
            while (j) rem |= lowbit(j), j -= lowbit(j);
        }
    }
    if (rem != N - 1) puts("INF"), exit(0);
    FMT(p, N, 1);
    for (int i = 0; i < N - 1; i += 1) p[i] = p[i] / (p[i] - 1);
    p[N - 1] = 0, FMT(p, N, -1);
    printf("%.8f", p[N - 1] + 1);
    return 0;
}
C++

4. [WC2018] 州区划分 BZOJ – 5153

子集卷积模板题(雾

先状压一波,Ps表示点集 s内的点权之和

然后若 s内无欧拉回路,则 Gs=Ps

否则 Gs=0

Fs为集合 s的答案

则有递推式 Fs=TSFST×(GTPS)p

(这里的小写 p是题目给出的 p次方)

然后把 Ps提出来:

Fs=1PpsTSFST×GpT

这很显然就是个子集卷积,搞搞就 A 了。

代码:

#pragma GCC optimize("Ofast,no-stack-protector")

#include <bits/stdc++.h>

#define P (998244353)
#define NS (22)
#define PS (1 << 21)

using namespace std;

template<typename _Tp> inline void IN(_Tp& dig)
{
    char c; bool flag = 0; dig = 0;
    while (c = getchar(), !isdigit(c)) if (c == '-') flag = 1;
    while (isdigit(c)) dig = dig * 10 + c - '0', c = getchar();
    if (flag) dig = -dig;
}

int pls(int a, int b) {return a + b >= P ? a + b - P : a + b;}
int mns(int a, int b) {return a - b < 0 ? a - b + P : a - b;}
int mul(int a, int b) {return 1ll * a * b % P;}

int qpow(int a, int b)
{
    int res = 1; a %= P;
    while (b)
    {
        if (b & 1) res = mul(res, a);
        a = mul(a, a), b >>= 1;
    }
    return res;
}

void FMT(int (&a)[PS], int N, int f)
{
    for (int l = 1; l < N; l <<= 1)
        for (int i = 0; i < N; i += (l << 1))
            for (int j = i; j < i + l; j += 1)
                if (f == 1) a[j + l] = pls(a[j + l], a[j]);
                else a[j + l] = mns(a[j + l], a[j]);
}

int n, m, p, w[NS], bc[PS], fa[NS], sum[PS], iv[PS], fum[NS][PS], f[NS][PS];

bool g[NS][NS];

int findf(int a) {return fa[a] == a ? a : fa[a] = findf(fa[a]);}

int main(int argc, char const* argv[])
{
    IN(n), IN(m), IN(p);
    for (int i = 1, a, b; i <= m; i += 1)
        IN(a), IN(b), a--, b--, g[a][b] = g[b][a] = 1;
    for (int i = 0; i < n; i += 1) IN(w[i]);
    for (int i = 1; i < (1 << n); i += 1)
    {
        bc[i] = bc[i >> 1] + (i & 1);
        for (int j = 0; j < n; j += 1)
            if ((i >> j) & 1) sum[i] = pls(sum[i], w[j]);
        sum[i] = qpow(sum[i], p), iv[i] = qpow(sum[i], P - 2);
        if (bc[i] == 1) continue;
        for (int j = 0, deg; j < n; j += 1) if ((i >> j) & 1) fa[j] = j;
        for (int j = 0, deg; j < n; j += 1)
            if ((i >> j) & 1)
            {
                deg = 0;
                for (int k = 0; k < n; k += 1)
                    if (g[j][k] && ((i >> k) & 1)) deg++, fa[findf(k)] = findf(j);
                if (deg & 1) {fum[bc[i]][i] = sum[i]; goto end;}
            }
        for (int j = 0, cnt = 0; j < n; j += 1)
            if (((i >> j) & 1) && fa[j] == j)
                if (++cnt == 2) {fum[bc[i]][i] = sum[i]; goto end;}
        end : continue;
    }
    f[0][0] = 1, FMT(f[0], 1 << n, 1);
    for (int i = 1; i <= n; i += 1) FMT(fum[i], 1 << n, 1);
    for (int i = 1; i <= n; i += 1)
    {
        for (int j = 0; j < i; j += 1)
            for (int k = 0; k < (1 << n); k += 1)
                f[i][k] = pls(f[i][k], mul(f[j][k], fum[i - j][k]));
        FMT(f[i], 1 << n, -1);
        for (int j = 0; j < (1 << n); j += 1)
            if (bc[j] == i) f[i][j] = mul(f[i][j], iv[j]);
            else f[i][j] = 0;
        if (i != n) FMT(f[i], 1 << n, 1);
    }
    printf("%d\n", f[n][(1 << n) - 1]);
    return 0;
}
C++
分类: 文章

XZYQvQ

炒鸡辣鸡的制杖蒟蒻一枚QvQ

2 条评论

EternalTron · 2018年8月21日 8:19 上午

资瓷,大佬真厉害啊 OwO。

    XZYQvQ · 2018年8月21日 8:43 上午

    突然开% 猝不及防

    其实之前想 fAke 您的都忍住了 OvO

发表回复

Avatar placeholder

您的邮箱地址不会被公开。 必填项已用 * 标注