传送门

---

数字最小公倍数为\(L\)的充分条件是所有数都是\(L\)的约数，而\(10^8\)内最多约数的数的约数也只有\(768\)个。所以我们先暴力找到所有满足是\(L\)的约数、\(G\)的倍数的数。

接下来注意到题目的\(\gcd\)和\(lcm\)的限制等价于对于每一个质数限制所有数在该质数指数上的\(\min\)和\(\max\)。在\(10^8\)内质数数量最多的数只有\(8\)个质数，所以我们对于第一步中求出的数用一个二进制数记录下它每一个质数的指数是否等于限制的\(\min\)和\(\max\)。那么问题就变成了选择子集满足按位或为全集的方案数。

对于这个问题考虑容斥，即强制某些位置在或的过程当中不被取到，那么答案就是\(\sum\limits_{S \subseteq U} (-1)^{|U-S|} 2^{\sum\limits_{j \subseteq S}1}\)。

后面的\(\sum\limits_{j \subseteq S}1\)可以通过高维前缀和求得。

最后我们需要求出所有数字的答案，那么我们只需要先算出所有数都可选择的答案，然后把当前数从高维前缀和中删掉计算其他数可以选择的答案相减即可。

复杂度\(O(2^{2\omega(n)} d(n)+Q)\)，其中\(\omega(n),d(n)\)分别是质数个数、约数个数。

#include<bits/stdc++.h>
using namespace std;

const int MOD = 1e9 + 7; map < int , int > ans;
int prm , Pow2[800] , cnt[1 << 16] , N , G , L , Q; vector < int > num , val;

int poww(long long a , int b){
    int times = 1;
    while(b){
        if(b & 1) times = times * a % MOD;
        a = a * a % MOD; b >>= 1;
    }
    return times;
}

void calc(int id , int x){
    int p = G , q = L , cnt1 = 0 , cnt2 = 0;
    while(p % x == 0){++cnt1; p /= x;}
    while(q % x == 0){++cnt2; q /= x;}
    for(int j = 0 ; j < num.size() ; ++j){
        int cur = num[j] , cnt = 0;
        while(cur % x == 0){++cnt; cur /= x;}
        if(cnt == cnt1) val[j] |= 1 << id;
        if(cnt == cnt2) val[j] |= 1 << (id + prm);
    }
}

int main(){
    cin >> N >> G >> L >> Q;
    if(L % G){for(int i = 1 ; i <= Q ; ++i) cout << "0\n"; return 0;}

    for(int i = 1 ; i <= N && i * i <= L ; ++i)
        if(L % i == 0){
            if(i % G == 0) num.push_back(i);
            if(L / i <= N && L / i != i && (L / i) % G == 0) num.push_back(L / i);
        }

    int tp = L; val.resize(num.size());
    for(int i = 2 ; i * i <= tp ; ++i) if(tp % i == 0){++prm; while(tp % i == 0) tp /= i;}
    if(tp - 1) ++prm; tp = L; int tc = 0;
    for(int i = 2 ; i * i <= tp ; ++i)
        if(tp % i == 0){calc(tc++ , i); while(tp % i == 0) tp /= i;}
    if(tp - 1) calc(tc++ , tp);
 
    for(int i = 0 ; i < val.size() ; ++i) ++cnt[val[i]];
    for(int i = 0 ; i < 2 * prm ; ++i)
        for(int j = 0 ; j < 1 << (2 * prm) ; ++j)
            if(j >> i & 1) cnt[j] += cnt[j ^ (1 << i)];
    Pow2[0] = 1; for(int i = 1 ; i <= num.size() ; ++i) Pow2[i] = (Pow2[i - 1] << 1) % MOD;
    int all = 0;
    for(int j = 0 ; j < 1 << (2 * prm) ; ++j)
        all = (all + (__builtin_popcount(j) & 1 ? MOD - 1ll : 1ll) * Pow2[cnt[j]]) % MOD;
    
    for(int i = 0 ; i < num.size() ; ++i){
        int sum = 0;
        for(int j = 0 ; j < 1 << (2 * prm) ; ++j) if((val[i] & j) == val[i]) --cnt[j];
        for(int j = 0 ; j < 1 << (2 * prm) ; ++j)
            sum = (sum + (__builtin_popcount(j) & 1 ? MOD - 1ll : 1ll) * Pow2[cnt[j]]) % MOD;
        ans[num[i]] = (all + MOD - sum) % MOD;
        for(int j = 0 ; j < 1 << (2 * prm) ; ++j) if((val[i] & j) == val[i]) ++cnt[j];
    }

    while(Q--){int x; cin >> x; cout << ans[x] << endl;}
    return 0;
}