Define ∩ as bitwise AND, ∪ as bitwise OR, and ⊕ as bitwise XOR. The notation i ⊆ x means the set of bits in binary representation i is a subset of x, i.e., i ∩ x = i.

SOS DP (Sum Over Subsets)

Given an array a, define F(x) = Σ a_i for all i ⊆ x. The goal is to compute F(x) for all x.

Define dp[x][k] as the sum of a_i for all i that differ from x only in the lower k bits and satisfy i ⊆ x. Bits are zero-indexed.

Transition is as follows:

If the k-th bit of x is 0: dp[x][k] = dp[x][k-1]
If the k-th bit of x is 1: dp[x][k] = dp[x][k-1] + dp[x ^ (1 << k)][k-1]

This can be implemented with the second dimension optimized out:

int n = total_bits;
for (int bit = 0; bit < n; bit++) {
    for (int mask = (1 << bit); mask < (1 << n); mask++) {
        if (mask & (1 << bit)) {
            f[mask] += f[mask ^ (1 << bit)];
        }
    }
}

Complexity is O(n * 2^n), where n is the number of bits.

To compute sums over supersets instead of subsets, modify the condition:

for (int bit = 0; bit < n; bit++) {
    for (int mask = (1 << bit); mask < (1 << n); mask++) {
        if (!(mask & (1 << bit))) {
            f[mask] += f[mask ^ (1 << bit)];
        }
    }
}

This SOS DP concept can be applied to optimize multidimensional prefix sums.

Fast Walsh-Hadamard Transform (FWT)

FWT computes convolution for bitwise operations. Given sequences A and B, compute C where: C_i = Σ A_j * B_k for j ⊙ k = i, with ⊙ being AND, OR, or XOR.

Bitwise OR Convolution

We seek C_i = Σ A_j * B_k for j ∪ k = i. Define F(X)_i = Σ X_j for all j ⊆ i. Then: FA = F(A), FB = F(B), and FC_i = FA_i * FB_i. Now, FC_i = Σ C_j for all j ⊆ i. To recover C, apply the inverse SOS DP:

for (int bit = n - 1; bit >= 0; bit--) {
    for (int mask = (1 << bit); mask < (1 << n); mask++) {
        if (mask & (1 << bit)) {
            fc[mask] -= fc[mask ^ (1 << bit)];
        }
    }
}

Bitwise AND convolution follows similarly using superset sums.

Bitwise XOR Convolution

Define popcnt(x) as the number of 1s in x. Let x ∘ y = popcnt(x ∩ y) % 2. The property (x ∘ y) ⊕ (x ∘ z) = x ∘ (y ⊕ z) holds. Define F(X)_i = Σ_{i∘j=0} X_j - Σ_{i∘j=1} X_j. Then FC_i = FA_i * FB_i and FC_i = Σ_{i∘j=0} C_j - Σ_{i∘j=1} C_j.

To compute F(X):

for (int bit = 0; bit < n; bit++) {
    for (int mask = 0; mask < (1 << n); mask++) {
        if (mask & (1 << bit)) {
            int t = f[mask];
            f[mask] = (f[mask] + f[mask ^ (1 << bit)]) % MOD;
            f[mask ^ (1 << bit)] = (f[mask ^ (1 << bit)] - t + MOD) % MOD;
        }
    }
}

The inverse transform:

int inv2 = (MOD + 1) / 2; // modular inverse of 2
for (int bit = n - 1; bit >= 0; bit--) {
    for (int mask = 0; mask < (1 << n); mask++) {
        if (mask & (1 << bit)) {
            int t = fc[mask], u = fc[mask ^ (1 << bit)];
            fc[mask] = 1LL * (t - u + MOD) * inv2 % MOD;
            fc[mask ^ (1 << bit)] = 1LL * (t + u) * inv2 % MOD;
        }
    }
}

Example implementation for all three operations:

#include<cstdio>
const int MAXN = 1 << 17, MOD = 998244353, INV2 = 499122177;
int n, size;
int a[MAXN], b[MAXN], c[MAXN];
int fa[MAXN], fb[MAXN], fc[MAXN];

int main() {
    scanf("%d", &n);
    size = 1 << n;
    for (int i = 0; i < size; i++) scanf("%d", &a[i]);
    for (int i = 0; i < size; i++) scanf("%d", &b[i]);

    // OR convolution
    for (int i = 0; i < size; i++) fa[i] = a[i], fb[i] = b[i];
    for (int bit = 0; bit < n; bit++)
        for (int mask = 0; mask < size; mask++)
            if (mask & (1 << bit)) {
                fa[mask] = (fa[mask] + fa[mask ^ (1 << bit)]) % MOD;
                fb[mask] = (fb[mask] + fb[mask ^ (1 << bit)]) % MOD;
            }
    for (int i = 0; i < size; i++) fc[i] = 1LL * fa[i] * fb[i] % MOD;
    for (int bit = n - 1; bit >= 0; bit--)
        for (int mask = 0; mask < size; mask++)
            if (mask & (1 << bit))
                fc[mask] = (fc[mask] - fc[mask ^ (1 << bit)] + MOD) % MOD;
    for (int i = 0; i < size; i++) printf("%d ", fc[i]);
    puts("");

    // AND convolution
    for (int i = 0; i < size; i++) fa[i] = a[i], fb[i] = b[i];
    for (int bit = 0; bit < n; bit++)
        for (int mask = 0; mask < size; mask++)
            if (!(mask & (1 << bit))) {
                fa[mask] = (fa[mask] + fa[mask ^ (1 << bit)]) % MOD;
                fb[mask] = (fb[mask] + fb[mask ^ (1 << bit)]) % MOD;
            }
    for (int i = 0; i < size; i++) fc[i] = 1LL * fa[i] * fb[i] % MOD;
    for (int bit = n - 1; bit >= 0; bit--)
        for (int mask = 0; mask < size; mask++)
            if (!(mask & (1 << bit)))
                fc[mask] = (fc[mask] - fc[mask ^ (1 << bit)] + MOD) % MOD;
    for (int i = 0; i < size; i++) printf("%d ", fc[i]);
    puts("");

    // XOR convolution
    for (int i = 0; i < size; i++) fa[i] = a[i], fb[i] = b[i];
    for (int bit = 0; bit < n; bit++)
        for (int mask = 0; mask < size; mask++)
            if (mask & (1 << bit)) {
                int t = fa[mask];
                fa[mask] = (fa[mask] + fa[mask ^ (1 << bit)]) % MOD;
                fa[mask ^ (1 << bit)] = (fa[mask ^ (1 << bit)] - t + MOD) % MOD;
                t = fb[mask];
                fb[mask] = (fb[mask] + fb[mask ^ (1 << bit)]) % MOD;
                fb[mask ^ (1 << bit)] = (fb[mask ^ (1 << bit)] - t + MOD) % MOD;
            }
    for (int i = 0; i < size; i++) fc[i] = 1LL * fa[i] * fb[i] % MOD;
    for (int bit = n - 1; bit >= 0; bit--)
        for (int mask = 0; mask < size; mask++)
            if (mask & (1 << bit)) {
                int t = fc[mask], u = fc[mask ^ (1 << bit)];
                fc[mask] = 1LL * (t - u + MOD) * INV2 % MOD;
                fc[mask ^ (1 << bit)] = 1LL * (t + u) * INV2 % MOD;
            }
    for (int i = 0; i < size; i++) printf("%d ", fc[i]);
    puts("");
    return 0;
}