POJ - 1845 Sumdiv(分治)

Posted 2021-03-14 zzzzzzy

题意：求$A^{B}$的所有约数之和$mod 9901$

思路：由结论有，一个数$n$进行质因数分解得到$n={p_{1}}^{c_{1}} * {p_{2}}^{c_{2}} *...* {p_{k}}^{c_{k}}$，那么$n$的约数之和为

$$sum=(1+{p_{1}}^{1}+cdots+{p_{1}}^{c_{1}})*(1+{p_{2}}^{1}+cdots +{p_{2}}^{c_{2}})*cdots*(1+{p_{k}}^{1}+cdots+{p_{k}}^{c_{k}})$$

所以对$A$质因数分解后，那么$A^{B}$的约数之和

$$sum=(1+{p_{1}}^{1}+cdots+{p_{1}}^{B*c_{1}})*(1+{p_{2}}^{1}+cdots +{p_{2}}^{B*c_{2}})*cdots*(1+{p_{k}}^{1}+cdots+{p_{k}}^{B*c_{k}})$$

上式中每个括号内都是等比数列，利用分治法对等比数列求和，设$sum(p，c)=1+p+p^2+cdots+p^{c}$

当$c$为奇数时

$$sum(p,c)=(1+p+cdots+p^{frac{c-1}{2}})+(p^{frac{c+1}{2}}+cdots+p^c)=(1+p^{frac{c+1}{2}})*sum(p,frac{c-1}{2})$$

当$c$为偶数时

$$sum(p,c)=(1+p+cdots+p^{frac{c}{2}-1})+(p^{frac{c}{2}}+p^{frac{c}{2}+1}cdots+p^{c-1})+p^c=(1+p^{frac{c}{2}})*sum(p,frac{c}{2}-1)+p^c$$

当$c$等于$0$，结束递归， 返回$1$即可

#include <iostream>
#include <algorithm>
#include <cstdio>
#include <cmath>

using namespace std;

typedef long long ll;

const int N = 110;
const ll mod = 9901;

ll a, b;
ll p[N], c[N], m;

void divide(ll n)
{
    m = 0;
    for (ll i = 2; i <= sqrt(n); i++) {
        if (0 == n % i) {
            p[++m] = i, c[m] = 0;
            while (0 == n % i) n /= i, c[m]++;
        }
    }
    if (n > 1) p[++m] = n, c[m] = 1;
    return;
}

ll power(ll a, ll b, ll p)
{
    ll res = 1;
    while (b) {
        if (b & 1) res = (res * a) % p;
        a = (a * a) % p, b >>= 1;
    }
    return res % p;
}

ll sum(ll p, ll c)
{
    if (0 == c) return 1;
    if (1 == c % 2) {
        ll tp1 = (1 + power(p, (c + 1) / 2, mod)) % mod;
        ll tp2 = sum(p, (c - 1) / 2) % mod;
        return tp1 * tp2 % mod;
    }
    else {
        ll tp1 = (1 + power(p, c / 2, mod)) % mod;
        ll tp2 = sum(p, c / 2 - 1) % mod;
        return (tp1 * tp2 % mod + power(p, c, mod)) % mod;
    }
}

int main()
{
    scanf("%lld%lld", &a, &b);
    divide(a);
    if (0 == a) printf("0
");
    else {
        ll res = 1;
        for (int i = 1; i <= m; i++)
            res = res * sum(p[i], b * c[i]) % mod;
        printf("%lld
", res);
    }
    return 0;
}