ZOJ 3868 GCD 期望：容斥原理與數論優化完全解析

問題重述

給定 n 個整數 a₁..aₙ（1≤n,k≤10⁶，1≤aᵢ≤10⁶），等概率隨機選取非空子集，求子集 GCD 的 k 次方期望 E。輸出 E × (2ⁿ - 1) mod 998244353（即期望分子）。

數學建模

設 f(d) 為所有元素都是 d 的倍數的非空子集數量。則 GCD 恰為 d 的子集數 g(d) = f(d) - Σ g(kd) (k≥2)（容斥/莫比烏斯反演）。

目標答案：Ans = Σ g(d) × dᵏ。

算法流程

1. 統計每個值出現次數 cnt[v]
2. 對每個 d (1..maxA)，計算 f(d) = 2^{Σ cnt[kd]} - 1
3. 從大到小 DP：g(d) = f(d) - Σ_{k≥2} g(kd)
4. 累加 Ans = Σ g(d) × dᵏ mod MOD

關鍵優化

• 預計算 pow2[i] = 2ⁱ mod MOD 和 pow_k[d] = dᵏ mod MOD
• 倍數枚舉複雜度 Σ maxA/d = O(maxA log maxA)
• 總 maxA ≤ 2×10⁶，單測試用例 ~2×10⁷ 操作，極快

完整代碼

#include 
using namespace std;
using ll = long long;
const int MOD = 998244353;
const int MAXA = 2000000;

ll mod_pow(ll a, ll e) {
    ll r = 1;
    while (e) {
        if (e & 1) r = r * a % MOD;
        a = a * a % MOD;
        e >>= 1;
    }
    return r;
}

int cnt[MAXA + 5];
ll f[MAXA + 5], g[MAXA + 5], pw2[MAXA + 5], pwk[MAXA + 5];

int main() {
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
    
    pw2[0] = 1;
    for (int i = 1; i <= MAXA; ++i) pw2[i] = (pw2[i-1] * 2) % MOD;
    
    int T; cin >> T;
    while (T--) {
        int n, k; cin >> n >> k;
        fill(cnt, cnt + MAXA + 1, 0);
        int maxv = 0;
        for (int i = 0, x; i < n; ++i) {
            cin >> x;
            ++cnt[x];
            maxv = max(maxv, x);
        }
        
        for (int d = 1; d <= maxv; ++d) {
            int c = 0;
            for (int m = d; m <= maxv; m += d) c += cnt[m];
            f[d] = (pw2[c] - 1 + MOD) % MOD;
            pwk[d] = mod_pow(d, k);
        }
        
        ll ans = 0;
        for (int d = maxv; d >= 1; --d) {
            g[d] = f[d];
            for (int m = 2*d; m <= maxv; m += d) {
                g[d] = (g[d] - g[m] + MOD) % MOD;
            }
            ans = (ans + g[d] * pwk[d]) % MOD;
        }
        cout << ans << '\n';
    }
    return 0;
}

樣例驗證

輸入：1 5 1 1 2 3 4 5

• 非空子集共 31 個
• GCD=1: 26 個，GCD=2: 2 個，GCD=3,4,5: 各 1 個
• 期望分子 = 26×1 + 2×2 + 3 + 4 + 5 = 42 → 輸出 42 ✓

核心技巧總結