Bzoj5019: [Snoi2017]遗失的答案

将 (L) 唯一分解为 (p_1^{a_1}p_2^{a_2}...p_k^{a_k})
对 G 也分解为 (p_1^{b_1}p_2^{b_2}...p_k^{b_k})。
称 (a_i , b_i) 分别为 (p_i) 这个质因子幂次的上下界。
显然为了满足 (gcd) 为 (G) 且 (lcm) 为 (L),对于每一个 (p_i) ,就至少要
有一个数触其上界,有一个数触其下界。
那么就可以拿一个长为 (2^k) 的二进制状态表示每个质因子是否已
触其上界,是否已触其下界。
然后前后缀对于 (L) 的约数且是 (G) 的倍数的数作一遍 (DP)
直接 (FWT) 合并成为答案
还有更加优秀的容斥做法
首先问题已经转化为给定一些集合，求出或为全集的方案数
设 (f_s) 表示或为 (s) 的方案数
直接求不方便
设 (g_s=sum_{isubset s}f_i) 表示或为 (s) 的子集的方案数
那么 (g_s=2^{cnt_s}-1)((cnt_s) 表示 (s) 的子集个数)
那么容斥得到 (f_s=sum_{isubset s}(-1)^{|s|-|i|}g_i)
对于钦定了一定要选集合 (x)，只需要把 (x) 的超集 (-1) 然后再次做上面的容斥
记得 (g_s=2^{cnt_s-1})

# include <bits/stdc++.h>
using namespace std;
typedef long long ll;
 
namespace IO {
    const int maxn(1 << 21 | 1);
 
    char ibuf[maxn], *iS, *iT, c;
    int f;
 
    inline char Getc() {
        return iS == iT ? (iT = (iS = ibuf) + fread(ibuf, 1, maxn, stdin), (iS == iT ? EOF : *iS++)) : *iS++;
    }
 
    template <class Int> inline void In(Int &x) {
        for (f = 1, c = Getc(); c < '0' || c > '9'; c = Getc()) f = c == '-' ? -1 : 1;
        for (x = 0; c >= '0' && c <= '9'; c = Getc()) x = (x << 1) + (x << 3) + (c ^ 48);
        x *= f;
    }
}
 
using IO :: In;
 
const int maxn(805);
const int maxm(16);
const int mod(1e9 + 7);
 
int n, g, l, bl, q, p[maxm], mx[maxm], mn[maxm], tot, st, cnt;
int size[1 << maxm], ans[maxn], pw[1 << maxm], bit[1 << maxm];
 
struct Fac {
    int x, y;
 
    inline bool operator <(Fac b) const {
        return x < b.x;
    }
} a[maxn];
 
inline void Inc(int &x, int y) {
    x += y;
    if (x >= mod) x -= mod;
}
 
inline void Pre_Work(int x) {
    if (x > n || x % g) return;
    register int i, v;
    a[++cnt].x = x;
    for (i = 0; i < tot; ++i) {
        v = 0;
        while (x % p[i] == 0) x /= p[i], ++v;
        if (v == mn[i]) a[cnt].y |= 1 << i;
        if (v == mx[i]) a[cnt].y |= 1 << (i + tot);
    }
}
 
int main() {
    In(n), In(g), In(l), In(q);
    register int qq, i, j, x, y, k, tp;
    if (l % g) {
        for (qq = 1; qq <= q; ++qq) puts("0");
        return 0;
    }
    x = l;
    for (i = 2; i * i <= x; ++i)
        if (x % i == 0) {
            while (x % i == 0) x /= i;
            p[tot++] = i;
        }
    if (x > 1) p[tot++] = x;
    st = 1 << (tot << 1), x = g, y = l;
    for (i = 0; i < tot; ++i) {
        while (x % p[i] == 0) x /= p[i], ++mn[i];
        while (y % p[i] == 0) y /= p[i], ++mx[i];
    }
    for (i = 1; i * i <= l; ++i)
        if (l % i == 0) {
            Pre_Work(i);
            if (i * i != l) Pre_Work(l / i);
        }
    sort(a + 1, a + cnt + 1);
    for (i = 1; i <= cnt; ++i) ++size[a[i].y];
    for (i = 1; i < st; i <<= 1)
        for (tp = i << 1, j = 0; j < st; j += tp)
            for (k = 0; k < i; ++k) size[j + k + i] += size[j + k];
    for (pw[0] = 1, j = size[st - 1], i = 1; i <= j; ++i) pw[i] = (pw[i - 1] + pw[i - 1]) % mod;
    for (i = 1; i < st; ++i) bit[i] = bit[i >> 1] + (i & 1);
    for (i = 1; i <= cnt; ++i)
        for (j = 0; j < st; ++j)
            if ((a[i].y & j) == a[i].y)
                Inc(ans[i], (bit[(st - 1) ^ j] & 1) ? mod - pw[size[j] - 1] : pw[size[j] - 1]);
    for (qq = 1; qq <= q; ++qq) {
        In(x);
        if (l % x || x % g) puts("0");
        else {
            y = lower_bound(a + 1, a + cnt + 1, (Fac){x, 0}) - a;
            a[y].x != x ? puts("0") : printf("%d
", ans[y]);
        }
    }
    return 0;
}