題面

對于正整數 \(n\)，定義 \(f(n)\) 為 \(n\) 所含質因子的最大冪指數。例如 \(f(1960)=f(23×51×72)=3,f(10007)=1,f(1)=0\)。給定正整數 \(a,b\)，求下式的值：

(資料圖片)

\[\sum_{i = 1}^a\sum_{j = 1}^bf(\gcd(a,b))\]

題解

在下文的推導中，為避免歧義，設 \(N = \min(a,b), M = \max(a,b)\)，顯然 \(N \le M\)。

\[\displaystyle \begin{aligned} & \sum_{i = 1}^N \sum_{j = 1}^M f(\gcd(a,b))\\= & \sum_{d = 1}^N f(d) \sum_{i = 1}^N \sum_{j = 1}^M \left[\gcd(i,j) = d\right] \\= & \sum_{d = 1}^N f(d) \sum_{i = 1}^{\left\lfloor\frac{N}plhniruxn\right\rfloor}\sum_{j = 1}^{\left\lfloor\frac{M}plhniruxn\right\rfloor} \left[\gcd(i,j) = 1\right] \\ = & \sum_{d = 1}^N f(d) \sum_{i = 1}^{\left\lfloor\frac{N}plhniruxn\right\rfloor}\sum_{j = 1}^{\left\lfloor\frac{M}plhniruxn\right\rfloor} \varepsilon(\gcd(i,j)) \\= & \sum_{d = 1}^N f(d) \sum_{i = 1}^{\left\lfloor\frac{N}plhniruxn\right\rfloor}\sum_{j = 1}^{\left\lfloor\frac{M}plhniruxn\right\rfloor} \sum _{t \mid \gcd(i, j)} \mu(t)\\= & \sum_{d = 1}^N f(d) \sum_{i = 1}^{\left\lfloor\frac{N}plhniruxn\right\rfloor}\sum_{j = 1}^{\left\lfloor\frac{M}plhniruxn\right\rfloor} \sum _{t \mid i \land t \mid j} \mu(t)\\\end{aligned}\]

設 \(T = dt\)，那么：

\[\displaystyle \begin{aligned} & \sum_{d = 1}^N f(d) \sum_{i = 1}^{\left\lfloor\frac{N}plhniruxn\right\rfloor}\sum_{j = 1}^{\left\lfloor\frac{M}plhniruxn\right\rfloor} \sum _{t \mid i \land t \mid j} \mu(t)\\= & \sum_{T = 1}^N {\left\lfloor\frac{N}{T}\right\rfloor} {\left\lfloor\frac{M}{T}\right\rfloor} \sum_{d \mid T} f(d) \mu(\frac{T}plhniruxn)\end{aligned}\]

設 \(h = f * \mu\)，那么原式可化為：

\[\displaystyle \sum_{T = 1}^N {\left\lfloor\frac{N}{T}\right\rfloor} {\left\lfloor\frac{M}{T}\right\rfloor} h(T)\]

下面考慮求函數 \(h(n)\) 的值，首先對 \(n\) 進行質因數分解：

\[\displaystyle n = \prod_{i = 1}^m p_i^{c_i} \ ( \ p_i \in \mathbb{P}, c_i \ge 1 \ )\]

發現可以把 \(h(n)\) 寫成如下形式：

\[\displaystyle h(n) = \sum_{ab = n} f(a) \mu(b)\]

考慮到莫比烏斯函數 \(\mu(n)\) 的性質：如果 \(n\) 中含有平方質因子那么 \(\mu(n) = 0\)。所以可以得出能產生貢獻的 \(b\) 即滿足 \(\mu(b) \ne 0\) 的 \(b\) 一定滿足：

\[\displaystyle b = \prod_{i = 1}^m p_i^{d_i} \ ( \ p_i \in \mathbb{P}, _i \in \{0, 1\} \ )\]

所以可以得出 \(f(a) = l \lor f(a) - l - 1\)。

設 $l = \max \limits_{i = 1}^m c_i,k = \sum \limits_{i = 1}^m \left[ c_i = l\right] $。接下來按 \(k \ne m\) 和 \(k = m\) 兩種情況分類討論 \(h(n) 的值\)。

當 \(k \ne m\) 時，按 \(f(a) = l\) 和 \(f(a) = l - 1\) 兩種子情況討論。

當 \(f(a) = l\) 時，設在 \(k\) 個滿足 \(c_i = l\) 的質數中選了 \(t\) 個，在另外 \(m - k\) 個質數中選了 \(s\) 個，那么可以得出貢獻為：

\[\displaystyle \begin{aligned} & \sum_{s = 0} ^ {m - k} \sum_{t = 0}^{k - 1} \dbinom{k}{t} \times l \times (-1) ^ {s + t} \times \dbinom{m - k}{s}\\= & \sum_{t = 0}^{k - 1} \dbinom{k}{t} \times l \times (-1) ^ {t} \sum_{s = 0} ^ {m - k} (-1) ^ {s} \times 1^{m - k - s} \dbinom{m - k}{s} \\= & 0\end{aligned}\]

當 \(f(a) = l - 1\) 時，\(k\) 個滿足 \(c_i = l\) 的質數中一定全部選上，設在另外 \(m - k\) 個質數中選了 \(s\) 個，那么可以得出貢獻為：

\[\displaystyle \begin{aligned} & \sum_{s = 0} ^ {m - k} (l - 1) \times (-1) ^ {s} \times \dbinom{m - k}{s}\\= & \sum_{s = 0} ^ {m - k} (l - 1) \times (-1) ^ {s}\times 1^{m - k - s} \dbinom{m - k}{s} \\= & (l - 1) \times \sum_{s = 0} ^ {m - k} (-1) ^ {s}\times 1^{m - k - s} \dbinom{m - k}{s}\\= & 0\end{aligned}\]

當 \(k = m\) 時，有 \(f(a) = l\) 和 \(f(a) = l - 1\) 兩種子情況，可以得出貢獻為：

\[\displaystyle \begin{aligned} & \sum_{s = 0} ^ {m - 1} (l) \times (-1) ^ {s} \times \dbinom{m - k}{s} + (l - 1) \times (-1) ^ {m} \times \dbinom{m - k}{m - k}\\= & \sum_{s = 0} ^ {m} (l) \times (-1) ^ {s} \times \dbinom{m - k}{s} - 1 \times (-1) ^ m\\= & - 1 \times (-1) ^ m \\= & (-1) ^ {m + 1}\end{aligned}\]

綜上可以得出 \(h(n)\) 的計算式：

\[h(n)=\begin{cases} 0 & k \ne m\\(-1)^{m + 1} & k = m\end{cases}\]

下面考慮如何高效的預處理出 \(h(n)\) 的值。

考慮一下歐拉篩在篩出合數 \(\displaystyle n = \prod_{i = 1}^m p_i^{c_i} \ ( \ p_i \in \mathbb{P}, c_i \ge 1 , p_i < p_{i + 1})\) 時的路徑：

\[p_m \rightarrow p_m^2 \rightarrow p_m^3 \rightarrow \cdots \rightarrow p_m^{c_m} \rightarrow \\p_{m - 1} p_m^{c_m} \rightarrow p_{m - 1}^2 p_m^{c_m} \rightarrow \cdots n\]

也就是說一個合數被篩出的路徑是按質因子從大到小的順序篩出的，所以可以開三個數組 \(preCount,nowCount\) 和 \(factorCount\)，分別記錄當前數的最小質因子的冪次，其他所有質因子的冪次和本質不同的質因子的個數。

如果在篩的過程中，設 \(i\) 為當前篩的數， \(j\) 為枚舉的質數，\(t = i \cdot j\)，如果 \(j \nmid i\)，也就是說 \(j\) 不是 \(i\) 的質因子，但是其是 \(t\) 的最小質因子，所以 \(nowCount[t] = 1\)，但是 \(preCount[t]\) 的值要根據 \(nowCount[i]\) 和 \(preCount[i]\) 的值來考慮：如果兩者相等，直接賦值即可；否則直接賦 \(-1\)，也就是說現在這個歐拉篩上的路徑上的數的質因子冪次不可能相等了。如果 \(j \mid i\)，直接累加即可。

Code

#include typedef long long valueType;constexpr valueType maxN = 1e7 + 5;class LineSieve {public:    typedef long long valueType;    typedef std::vector container;private:    valueType size;    container minFactorList;    container primeList;    container preCount, nowCount, factorCount, data, sum;public:    explicit LineSieve(valueType _size_) : size(_size_), minFactorList(_size_ + 1),     preCount(_size_ + 1, 0),nowCount(_size_ + 1, 0),     factorCount(_size_ + 1), data(_size_ + 1),sum(_size_ + 1) {        primeList.reserve((size_t) std::floor(std::log((long double) (_size_ + 1))));        for (valueType i = 2; i <= size; ++i) {            if (minFactorList[i] == 0) {                primeList.push_back(i);                minFactorList[i] = i;                nowCount[i] = 1;                preCount[i] = 0;                factorCount[i] = 1;                data[i] = 1;            }            for (auto const &iter: primeList) {                valueType const t = i * iter;                if (t > size)                    break;                minFactorList[t] = iter;                if (i % iter == 0) {                    nowCount[t] = nowCount[i] + 1;                    preCount[t] = preCount[i];                    factorCount[t] = factorCount[i];                    break;                } else {                    nowCount[t] = 1;                    preCount[t] = (nowCount[i] == preCount[i] || preCount[i] == 0) ? nowCount[i] : -1;                    factorCount[t] = factorCount[i] + 1;                }            }        }        for (int i = 2; i <= size; ++i)            if (nowCount[i] == preCount[i] || preCount[i] == 0)                data[i] = (factorCount[i] & 1) == 1 ? 1 : -1;            else                data[i] = 0;        std::partial_sum(data.begin(), data.end(), sum.begin());    }    valueType ans(valueType x) const {        if (x > size)            throw std::range_error("Larger than Size.");        if (x < 0)            throw std::range_error("Too small.");        return sum[x];    }};int main() {    valueType T;    std::cin >> T;    LineSieve Euler(maxN);    typedef std::function solveFunction;    solveFunction solve = [&Euler](valueType N, valueType M) -> valueType {        if (N > M)            std::swap(N, M);        valueType result = 0;        valueType l = 1, r;        while (l <= N) {            r = std::min(N / (N / l), M / (M / l));            result += (Euler.ans(r) - Euler.ans(l - 1)) * (N / l) * (M / l);            l = r + 1;        }        return result;    };    for (int i = 1; i <= T; ++i) {        valueType N, M;        std::cin >> N >> M;        std::cout << solve(N, M) << "\n";    }    std::cout << std::flush;    return 0;}