约数之和问题

约数之和

0x01 背景题目

0. 定理

算术基本定理（正整数唯一分解定律）：

不考虑排列顺序的情况下，每个正整数都能够以唯一的方式表示成它的质因数的乘积。

$x={p_1}^{k_1}\ast {p_2}^{k_2}\ast {p_3}^{k_3}.....{p_n}^{k_n}$

人话：对于每个大于 1 的自然数，要么本身是质数，要么可以写为两个或以上的质数的积。

哥德巴赫猜想：

任意一个大于 2 的偶数都可以分解为两个质数之和。

弱哥德巴赫猜想：

大于 5 的奇数都可以表示成 3 个素数之和

（实际上就是哥德巴赫猜想的另一种表达方式，因为一个奇数可以表示成一个奇数加上一个偶数。如果我们让这个奇数等于 3，那么就转化为了哥德巴赫猜想）

半质数：

一个可以表示为两个质数的乘积的合数成为半质数。

任意一个合数不一定能分解为两个质数的乘积。

1. 分解质因数

题目描述：

给定 n 个正整数 $a_i$，将每个数分解质因数，并按照质因数从小到大的顺序输出每个质因数的底数和指数。

代码：

#include <iostream>
#include <unordered_map>

using namespace std;

int main()
{
    int n;  cin >> n;
    while (n -- )
    {
        int x;  cin >> x;
        // 看似枚举的是每一个约数，其实枚举的是每一个因子
        // 因子一定是质数，递归应用唯一分解定律即可
        // 例如，4可以被分解为2*2。
        // 所以说，对于16，他永远不会在这里mod4==0
        // 因此4一定被分解了没了（变成两个2）。
        // 另外，最多只有一个大于 sqrt(x) 的因子
        // 否则两个因子的乘积会大于 x
        for(int i = 2; i <= x / i; i ++ )
        {
            if(x % i == 0)
            {
                int s = 0;
                while(x % i == 0)
                {
                    x /= i;
                    s ++ ;
                }
                cout << i << ' ' << s << endl;
            }
        }
        // 可能没有分解完
        // 注意是 if(x>1) 而不是 if(x)
        if(x > 1)   cout << x << ' ' << 1 << endl;
        cout << endl;
    }
    return 0;
}

关于最后 if(x>1) 的解释

最后剩下的数如果大于1， 并不能说明这个数一定是大于sqrt（n）的一个质因数

例如分解96的质因数，96 = 2^5 * 3，当分解完质因数2后，n已经为3，此时i不满足 i <= n / i，会进入后面一个判断，而3并不是大于sqrt（96）的数

因此最后需要特判。

注意，是判断 if(x>1)，因为 1 不是质数。我们不包含它，包含它也没啥用。

参考：

reference here

2. 约数之和

题目描述：

给定 $n$个正整数 $a_i$，请你输出这些数的乘积的约数之和，答案对 $10^9+7$ 取模。(a[i] <= $2\ast 10^9$, n<=100）

思路：

通过数据范围可以看出 a[i] 非常大，因此暴力相乘再 $O(\sqrt{n})$ 也是不可以的，因为即便是 longlong 也无法存下 $1e{9}^{100}$ 这种规模的数据。

我们可以转化思路：根据基本算术定理将一个数 x 表示为 ${\sum }_{i=1}^n{p}_i^{k_i}$，其中，p是质数，k是指数。

那么，$p_i^0,p_i^1,...,p_i^{k_i}$ 都是 x 的约数。

另外，如果 $p_1^0,p_1^1$ 和 $p_2^0,p_2^1$ 都是 x 的约数的话，他们的组合（乘积）也是 x 的约数。但是他们的和不一定是 x 的约数。

A%x==0 && A%y==0, --> A%(x*y)==0

A%x==0 && A%y==0, !–> A%(x+y)==0

那么 x 的约数之和就是 ${\prod }_{i=1}^n{\sum }_{j=0}^n{p}_i^j$

现在，如何求 $s(n)={\sum }_{j=0}^n{p}_i^j$ ?

秦九韶算法！

$s(n)=p^0+s(n-1)\ast p$，$s(0)=p^0=1$

代码：

我们上面讨论的是一个数的情况下求约数之和，对于多个数的乘积是相同的，因为多个数的乘积其实就是质因子的质数相加。

$p^a\ast p^b=p^{a+b}$

但是多个数之和就不能这么求了。

#include <iostream>
#include <unordered_map>

using namespace std;

const int mod = 1e9 + 7;
typedef long long ll;

unordered_map<int,int> primes;

int main()
{
    int n;  cin >> n;
    while (n -- )
    {
        int x;  cin >> x;
        for(int i = 2; i <= x / i; i ++ )
        {
            while(x % i == 0)
            {
                primes[i] ++ ;
                x /= i;
            }
        }
        if(x > 1)   primes[x] ++ ;
    }
    ll res = 1;
    for(auto &x : primes)
    {
        ll p = x.first, k = x.second;
        ll s = 1;
        while(k -- )  s = (s * p + 1) % mod;
        res = res * s % mod;
    }
    cout << res << endl;
    return 0;
}

参考：

详细参考：here

0x02 题目描述

算法提高课

0x03 思路

其实和背景题目的第二题差不多，只不过需要一点优化。

先考虑下 $A^B$，求约数我们是肯定不能暴力枚举约数的（实际上求约数的方法就是固定的，就是我们前面介绍的方法）。

我们首先需要分解质因数，那么，如何分解 $A^B$ 的因数呢？很简单，先分解 $A$ 的因数，由 $(ABC{)}^k=A^k{B}^k{C}^k$ 即可求的 $A^B$ 的因数，实际上质因数个数不变，修改一下指数就行了。

然后就是求 $\sum {p_i}^{k_i}$，在前面，我们介绍了秦九韶算法，它的时间复杂度是 $O(N)$ 的（N 是指数的大小）。但在这里，我们的指数可能非常大，因为 B 的范围是 $O(5\ast 10^7)$，因此我们需要优化这一步。详细内容参考 $0x05$。

另外，你可能会问为啥不用等比数列求和公式？因为等比数列求和涉及到除法，需要用到逆元。

0x04 代码

#include <iostream>
#include <unordered_map>

using namespace std;

typedef long long ll;
const int mod = 9901;

int A, B;
unordered_map<int,int> primes;

// 分解质因数
void divide(int n)
{
    for(int i = 2; i <= n / i; i ++ )
    {
        while(n % i == 0)
        {
            primes[i] ++ ;
            n /= i;
        }
    }       
    // 注意是 if(n>1) 不是 if(n!=0)
    if(n > 1)   primes[n] ++ ;
}

// 快速幂
int qsm(int a, int b)
{
    int res = 1;
    while(b)
    {
        if(b & 1)   res = (ll)res * a % mod;
        b >>= 1;
        a = (ll)a * a % mod;
    }
    return res;
}

// 求等比数列和
int sum(int p, int k)
{
    if(k == 1)  return 1; // 边界
    if(k % 2 == 0)
        return (ll)(qsm(p, k / 2) + 1) * sum(p, k / 2) % mod;
    // 奇数，去掉一项，并且只能去掉最后一项
    // 因为我们要保证连续性
    return ((ll)qsm(p, k - 1) + sum(p, k - 1)) % mod;
    return 0;
}

int main()
{
    cin >> A >> B;
    // special case
    if(A == 0) 
    {
        cout << 0 << endl;
        return 0;
    }
    divide(A);
    int res = 1;
    for(auto &it : primes) 
    {
        int p = it.first, k = it.second * B;
        res = (ll)res * sum(p, k + 1) % mod;
    }
    cout << res << endl;
    return 0;
}

0x05 参考

reference