🔗 UVA-12546 LCM Pair Sum

tags: `數學(Math)`

範例程式碼已於UVA、瘋狂程設(CPE) 上皆測試通過。(此題無 ZeroJudge )

題意

給定 $n=\prod_{i=1}^m p_i^{a_i}$ ，求：

$f(n)=\sum_{\substack{1 \leq x \leq y \leq n \\ \operatorname{lcm}(x, y)=n}}(x+y)$

由於答案可能很大，輸出對 $10^9 + 7$ 取餘數的結果。

思路：乘法原理

雖然直覺上就是用乘法原理，但這裡我並沒有詳細證明從一個質因數推廣到多個質因數的過程。可以參考參考資料2.，也許能為你提供一點想法。

公式推導

首先先忽略 $x \leq y$ $x \leq y$ 的限制，考慮所有的 $(x, y)$ $(x, y)$ 對，由於 $\operatorname{lcm}(x, y)=n$ $lcm (x, y) = n$ ，所以除了 $(n, n)$ $(n, n)$ 以外的所有 $(x, y)$ $(x, y)$ 對都會被計算 $2$ $2$ 次。
- 換句話說， $x = y$ 的情況只會出現在 $(n, n)$ 這一對中，其他的 $(x, y)$ 對中都存在大小關係。
再來先考慮只有一個質因數的情況，即令 $n=p_i^{a_i}$ $n = p_{i}^{a_{i}}$ ，考慮質因數 $p_i$ $p_{i}$ 在 $x$ $x$ 中的次方數 $k$ $k$ ：
1. 如果 $p_i$ 在 $x$ 中的次方數 $k < a_i$ ，那麼 $p_i$ 在 $y$ 中次數必須 $= a_i$ ，此時 $(x, y) = (p_i^k, p_i^{a_i})$ ，故此種情況下 $x$ 對 $f(n)$ 的貢獻為 $\sum_{k=0}^{a_i-1} p_i^k$
2. 如果 $p_i$ 在 $x$ 中的次方數 $k = a_i$ ，那麼 $p_i$ 在 $y$ 中次數可以是 $0$ 到 $a_i$ ，共 $a_i + 1$ 種情況，此時 $(x, y) = (p_i^{a_i}, p_i^{k_y})$ ，故此種情況下 $x$ 對 $f(n)$ 的貢獻為 $(a_i + 1) \cdot p_i^{a_i}$
合併兩種情況，得出 $x$ 對 $f(n)$ 的貢獻為 $\sum_{j=0}^{a_i} p_i^j+a_i \cdot p_i^{a_i}$ 。同理， $y$ 對 $f(n)$ 的貢獻也是一樣的，故需乘 $2$ 。但由於除了 $(n, n)$ 以外的所有 $(x, y)$ 對都會被計算 $2$ 次，所以最後答案為 $(\sum_{j=0}^{a_i} p_i^j+a_i \cdot p_i^{a_i}) + n$ 或是 $(\sum_{j=0}^{a_i-1} p_i^j+(a_i+1) \cdot p_i^{a_i}) + n$
最後利用乘法原理，擴展到有多個質因數的情況，即 $n=\prod_{i=1}^m p_i^{a_i}$ ，答案為 $\prod_{i=1}^m (\sum_{j=0}^{a_i} p_i^j+a_i \cdot p_i^{a_i}) + n$ 。

實作過程

雖然推出了公式，但直接寫肯定是會 TLE 的，我們不能傻傻的一項一項算。
由於要從 $a_i^0$ 算到 $a_i^{a_i}$ ，每一項都會算到，因此我們可以保留前一項的結果，這樣就可以在 $O(a_i)$ 的時間內算出 $\sum_{j=0}^{a_i} p_i^j+a_i \cdot p_i^{a_i}$ 。
此外，在計算過程中也要記得取餘數，避免數字過大。
時間複雜度為： $O(C \times a)$ ，其中 $C$ 為質因數個數， $a$ 為質因數中最大的次方數，對於每個詢問。
空間複雜度為： $O(C)$ ，輸入使用的空間。

MOD = 10**9 + 7
T = int(input())

for tc in range(1, T+1):
    C = int(input())
    PA = [list(map(int, input().split())) for _ in range(C)]
    
    ans, n = 1, 1
    for (p, a) in PA:
        sigma, cur = 1, 1
        for _ in range(1, a+1): # p^1 + p^2 + ... + p^a
            cur = (cur * p) % MOD
            sigma = (sigma + cur) % MOD
        ans *= (sigma + a * cur) % MOD
        n *= cur
    ans = (ans + n) % MOD
    print("Case %d: %d" % (tc, ans))

#include <bits/stdc++.h>
using namespace std;
using LL = long long;
const int MOD = 1e9 + 7;
const int N = 20;

int main() {
    ios::sync_with_stdio(false); cin.tie(nullptr); cout.tie(nullptr);
    int T, C;
    LL n, ans, sigma, cur;
    LL P[N], A[N]; // p_i, a_i
    cin >> T;
    for (int t=1; t<=T; t++) {
        cin >> C;
        for (int i=0; i<C; i++) {
            cin >> P[i] >> A[i];
        }
        ans = 1; n = 1;
        for (int i=0; i<C; i++) {
            sigma = 1; cur = 1;
            for (int j=1; j<=A[i]; j++) {
                cur *= P[i];
                cur %= MOD;
                sigma += cur;
                sigma %= MOD;
            }
            sigma += A[i] * cur % MOD;
            ans *= sigma % MOD;
            ans %= MOD;
            n *= cur;
            n %= MOD;
        }
        ans = (ans + n) % MOD;
        cout << "Case " << t << ": " << ans << endl;
    }
    return 0;
}