题面

牧羊人 A 和牧羊人 B 总是很无聊，所以他们要玩一个游戏。A 有 \(a\) 只羊，B 有 \(b\) 只羊。他们想要知道 \(a^b\) 的因子和是多少。这就很为难两个牧羊人了，由于答案太大，你能不能告诉我答案取模 \(9901\) 的数。

Example In #1

2 3

Example Out #1

15

对于 \(100%\) 的数据，\(0 \leq a, b \leq 50000000\)。

题解

举例：

\[\begin{align*} S &= (1 + 2 + 4)(1 + 3 + 9)\\ &= (2 + 4 + 3 + 9) + (2\times3) + (2\times9) + (4\times3) + (4\times9)\\ &= 7\times13 = 91 \end{align*} \]

可以看出 \(S\) 即为 \(2^3\times3^2 = 72\) 的因数和。
推广到一般：

\[\begin{align} S = (1 + p_1^{1b} + p_1^{2b} + \dots + p_1^{k_1b}) \times \dots \times(1 + p_m^{1b} + p_m^{2b} + \dots + p_m^{k_mb}) \end{align} \]

等比数列求和公式：

\[ S_n = \frac{a_1 \times (1 - q^n)}{1 - q} = \frac{a_1 - a_nq}{1 - q} = \frac{a_nq - a_1}{1 - q} \]

\(n\) 为项数，\(S_n\) 为和，\(q\) 为公比，\(a_1\) 为首项，\(a_n\) 为末项。
此题中：

\[\begin{align*} &n = m (m 为质因数个数)\\ &q = p_i\\ &a_1 = 1\\ &a_n = p_i^{k_ib} \end{align*} \]

代入式中，可见：

\[\begin{align} S &= \prod_{i = 1}^{m} \frac{p_i ^ {k_ib + 1} - 1}{p_i - 1} \end{align} \]

其中，除法取模用费马小定理求乘法逆元。

\[\begin{align} S &= \prod_{i = 1}^{m} (p_i ^ {k_ib + 1} - 1) \times (p_i - 1) ^ {(mod - 2)} \end{align} \]

#include <iostream>
#include <utility>
#include <vector>
#include <cmath>

const int MOD = 9901;
using i64 = long long;

i64 q_pow(i64 a, int b)
{
    i64 res = 1;
    while (b > 0)
    {
        if (b % 2 == 1)
            res = res * a % MOD;
        a = a * a % MOD;
        b /= 2;
    }

    return res;
}

bool is_prime(int n)
{
    for (int i = 2; i <= std::sqrt(n); i++)
        if (n % i == 0)
            return false;
    return true;
}

int main()
{
    int a, b;
    std::cin >> a >> b;

    std::vector<std::pair<int, int>> fac;
    for (int i = 2; i <= std::sqrt(a); i++)
    {
        if (is_prime(i) && a % i == 0)
        {
            fac.emplace_back(i, 0);
            while (a % i == 0)
            {
                a /= i;
                fac.back().second++;
            }
            fac.back().second *= b;
        }
    }
    if (a != 1)
        fac.emplace_back(a, b);

    i64 ans = 1;
    for (auto f : fac)
    {
        i64 res = (q_pow(f.first, f.second + 1) - 1) % MOD;
        res = res * q_pow(f.first - 1, MOD - 2) % MOD;
        ans = ans * res % MOD;
    }

    std::cout << ans;
    return 0;
}