洛谷 P1403 [AHOI2005]约数研究 题解 【 约数个数、约数和（普通方法+筛法）】

时间  2021-08-14
标签 c++ 数组 ide 函数 性能 测试 spa get 数学 栏目 C&C++ 繁體版
原文   原文链接

 

最近本人脑洞大开，发现了一种线性筛约数个数和约数和的神奇方法。
目前网上的方法基本都是利用\(num[i]\)数组记录\(i\)最小的质因子的个数，而后进行转移。
其实能够省去\(num[i]\)数组，直接进行递推。c++

设\(n\)的惟一分解式为：数组

$n={p_1}^{r_1} {p_2}^{r_2} ... {p_k}^{r_k}$

1、\(n\)的约数个数公式：ide

$d(n) = (r_1+1) * (r_2+1) * ... * (r_k+1) $

简写为：函数

$d(n) = \displaystyle\prod_{i=1}^k (r_i+1) $

证实也很简单，以\(p_1\)为例，这个质数因子，能够选择\(0\)个，能够选择\(1\)个，...,最多能够选择\(r_1\)个，就是有\(r_1+1\)种选择的可能性，其它\(p_2,p_3,...,p_k\)都是如此，根据乘法原理，全部的可能性就是\((r_1+1) * (r_2+1) * ... * (r_k+1)\)。性能

举个栗子：
\(180= 2^2 * 3^2 * 5\)测试

约数个数\(=(1+2) * (1+2) * (1+1) =18\)spa

听说小学奥数比赛中常常考查约数个数与约数和的公式，我也没有参加过奥数比赛，不知道是否是真的。ci

在给出具体的实现代码以前，还须要了解一个质因子分解的数学性质：get

一个正整数的质因子，最多只有一个大于其平方根。

证实：用反证法，假设有超过两个质因子大于其平方根，那么两者相乘必定大于该数。得证。数学

下面给出朴素版本的约数个数代码：

#include <bits/stdc++.h>

using namespace std;
typedef long long LL;

/**
 * 根据 X =(p1^a1)* (p2^a2)*.......(pn^an)
   T = (a1+1)*(a2+1)*(a3+1)....*(an+1)(a1,a2,a3...an位 p1,p2,p3....pn的系数）
 */

/**
 * 功能：计算约数个数
 * @param n
 * @return
 */
LL getDivisorCount(LL x) {
    unordered_map<int, int> primes; //key:质数 value:个数
    //求质数因子
    for (int i = 2; i <= x / i; i++)
        while (x % i == 0) x /= i, primes[i]++; //primes[i]表示质数i因子的个数+1

    //若是还有质数，那就加上
    if (x > 1) primes[x]++;
    //公式大法
    LL res = 1;
    for (auto p : primes) res = res * (p.second + 1);
    return res;
}

/**
 感悟：
 这是单个处理约数的办法，能够找出全部质数因子，和每一个质数因子的个数。再经过公式(1+r1)*(1+r2)*...*(1+rk)，就能知道约数的个数了。
 优势：
 直白，公式，暴力

 缺点：
 不能处理在某个区间内筛出约数的场景，挨个计算，性能差，须要再有一个面对区间的约数筛法。
 */
LL res;

int main() {
    int n;
    cin >> n;
    for (int i = 1; i <= n; i++) cout << i << " " << getDivisorCount(i) << endl;
    return 0;
}

2、\(n\)的约数和公式：

$\sigma(n)=({p_1}^{0}+{p_1}^{1}+{p_1}^{2}+...+{p_1}^{r_1})({p_2}^{0}+{p_2}^{1}+{p_2}^{2}+...+{p_2}^{r_2})...({p_k}^{0}+{p_k}^{1}+{p_k}^{2}+...+{p_k}^{r_k}) $

简写为：

$\sigma(n)=\prod\limits_{i=1}^{k}(1+p_{i}+p_{i}^{2}+\cdots+p_{i}^{r_{i}}) $

Sigma（大写Σ，小写σ），是第十八个希腊字母。

仔细观察一下就知道了，\(p_1^0+p_1^1+p_1^2+...+p_1^{r_1}\)是一个等比数列的求和，有公式滴，不用白不用：

$\frac{{p_1}^{a_1+1}-1}{p_1-1} * \frac{{p_2}^{a_2+1}-1}{p_2-1} * ... * \frac{{p_k}^{a_k+1}-1}{p_k-1}$

简写为：

$\prod\limits_{i=1}^{k} \frac{p_i^{a_i+1}-1}{p_i-1}$

举个栗子：
\(180= 2^2 * 3^2 * 5\)

约数和\(=(1+2+4) * (1+3+9 ) * (1+5)=546\)

下面给出朴素版本求约数和的代码：

#include <bits/stdc++.h>

using namespace std;
typedef long long LL;

/**
 * 功能：计算约数之和
 * @param n
 * @return
 */
LL getSumOfDivisors(LL x) {
    //拆出全部质数因子及质数因子个数
    unordered_map<int, int> primes;
    for (int i = 2; i <= x / i; i++)
        while (x % i == 0) {
            x /= i;
            primes[i]++;
        }
    if (x > 1) primes[x]++;

    //计算约数个数
    LL res = 1;
    for (auto p : primes) {
        LL a = p.first, b = p.second;
        LL t = 1;
        while (b--) t = (t * a + 1);
        res = res * t;
    }
    return res;
}

LL res;

int main() {
    //测试用例：180
    //参考答案：546
    int n;
    cin >> n;
    cout<<getSumOfDivisors(n) << endl;
    return 0;
}

小结：
上面给出了求约数个数、约数和的公式，适用于求指定数字的约数个数、约数和。对于从 \(1\)~\(n\)这类问题，效率低下，能够考虑采用相似于筛法的思路解决区间内约数个数与约数和的问题。

3、线性筛法求约数个数

要想线性筛一个积性函数\(f(i)\)，只须要知道4个东西。

扩展知识：
一、积性函数的定义：
积性函数指对于全部互质的整数\(a\)和\(b\)有性质\(f(a⋅b)=f(a)f(b)\)的数论函数。

1: \(f(1)=?\)
　　很显然，咱们这里的两个函数：\(d(1)=1,σ(1)=1;\)
\(d(1)=1\)意味着数字\(1\)的约数个数是\(1\)个。\(σ(1)=1\)意味着数字\(1\) 的约数和是\(1\),只有\(1\)个\(1\)，加在一块儿也是\(1\)。

2: \(f(p)=?\)，其中\(p\)为质数。

　　一样很显然，\(d(p)=2,σ(p)=p+1\)
\(d(p)=2\)意味着数字\(p\)的约数有\(2\)个，一个是\(1\)，另外一个是就是\(p\),不可能再有别的约数了，质数嘛。
\(σ(p)=p+1\) 约数是两个，\(1\)和\(p\)，那么约数和是肯定的，就是\(1+p\)。

三、4:\(f(i⋅p_j)=?\)，\(i\)与\(p_j\)互质或不互质。
先说线性筛约数个数的方法:

(1)若\(i\)与\(p_j\)互质
咱们能够利用积性函数的性质直接推出：\(d(i⋅p_j)=d(i)⋅d(p_j)=2d(i)\)
这里我有两个疑问：
疑问1：为啥\(d(p)\)是积性函数，怎么知道的？
疑问2：为啥\(d(p_j)==2\)?如今理解就是\(p_j\)在本文的筛法中是一个质数，\(i\)是它的倍数，既然\(p_j\)是质数，那么\(d(p_j)=2\),这是上面讨论过的内容。

(2)若\(i\)与\(p_j\)不互质
考虑\(i\)，\(i * p_j\)，\(\frac{i}{p_j}\)的关系(由线性筛过程可知,\(i\)与\(p_j\)不互质即\(p_j|i\))

这里不妨设\(p_j=p_1,r_j=r_1\)：

$i={p_1}^{r_1} {p_2}^{r_2} .... {p_k}^{r_k}$
$i p_1 = {p_1}^{r_1+1} {p_2}^{r_2} .... {p_k}^{r_k}$
$\frac{i}{p_1}={p_1}^{r_1-1} {p_2}^{r_2} ... {p_k}^{r_k}$
则有：
$d(i)=(r_1+1) (r_2+1) ... (r_k+1)$
$d(i*p_1)=(r_1+2) (r_2+1) ... (r_k+1)$
$d(\frac{i}{p_1})=r_1 (r_2+1) ... (r_k+1)$

设\(d(i)\)中\(r_1+1\)后面的一大坨为\(T\),便可表示为：

即: $T=(r_2+1) (r_3+1) ... (r_k+1) $

则：\(d(i)=(r_1+1)T\) ①

\(d(i⋅p_1)=(r_1+2)T=(r_1+1+1)T=(r_1+1)T+T= d(i)+T\) ②

\(d(\frac{i}{p_1})=r_1T=(r_1+1-1)T=(r_1+1)T-T=d(i)−T\) ③

将\(二、3\)式相加，整理得
$d(i⋅p_1) + d(\frac{i}{p_1}) = d(i)+T +d(i)−T = 2d(i) $

\(d(i⋅p_1)=2d(i)−d(\frac{i}{p_1})\)

由于最初设\(p_j=p_1,r_j=r_1\)：因此，就是

\(d(i⋅p_j)=2d(i)−d(\frac{i}{p_j})\)

说白了，费了这么半天劲，就是为了找出这个递推式。

有了递推式，咱们来看一下线性筛约数个数的代码：

#include <bits/stdc++.h>

using namespace std;
typedef long long LL;

//用途：线性筛约数个数的代码
const int N = 1e6 + 10;
int d[N];      //约数个数数组,记录每一个数字的约数个数
int primes[N]; //质数数组
bool st[N];     //是否是已经被筛掉了
int idx;       //质数数组的下标游标，从1开始，它的最终结果值，就是质数的数量

void getDivisorCount(int n) {
    //1的约数个数是1
    d[1] = 1;
    for (int i = 2; i <= n; i++) { //从2开始到n,进行筛
        if (!st[i])primes[++idx] = i, d[i] = 2; //若是没有被筛掉，那么是质数。记录到质数数组中，同时将d[i]设置为2.由于上文论证了，i为质数时，d[i]=2
        for (int j = 1; i * primes[j] <= n && j <= idx; j++) {//遍历每个已经出现的质数（欧拉筛大法好！），对primes[j]的整数倍进行筛除,i * primes[j] <= n表示是范围内，超过边界就不用了
            st[i * primes[j]] = true;//标识已筛出
            //i与primes[j]不互质
            if (i % primes[j] == 0) {
                d[i * primes[j]] = d[i] * 2 - d[i / primes[j]];//利用费了牛劲才得出的递推公式进行计算，由于是从小到大过来的，依赖项提早算好，因此能够算出来
                break;  //都说这个break值钱，只为第一个质数因子筛掉，欧拉筛的精髓，妙！
            }
            //i与primes[j]互质
            d[i * primes[j]] = d[i] * 2; //积性函数的性质直接推出,在上文中能够找到理由
        }
    }
}

LL res;

int main() {
    int n = 100;
    //线性筛约数个数
    getDivisorCount(n);
    for (int i = 1; i <= n; i++) res += d[i];
    cout << res << endl;
    return 0;
}

4、线性筛法求约数和
一、\(i\) 与\(p_j\)互质时

$\sigma(i\cdot p_{j})=\sigma(i)\cdot \sigma(p_{j})=\sigma(i)(p_{j}+1)$
解读：由于在本文的上下文中，$p_j$是个质数，而质数的$sigma(p)=p+1$，由于它只有两个约数，$1$和$p$,和固然是$p+1$了。

二、\(i\) 与\(p_j\)不互质时
这里对于上下文中，若是不互质，那么确定是\(p_j|i\),不存在其它有公约数的状况。
这时考虑\(\sigma(i)\),\(\sigma(\frac{i}{p_{1}})\),\(\sigma(i\cdot p_{1})\)三者之间的关联性。

疑问：
你要说考虑\(\sigma(i)\),\(\sigma(i\cdot p_{1})\)的关系，我能理解，由于递推嘛，确定要找二者之间的关系，怎么想到与\sigma(\frac{i}{p_{1}})也有关联呢？这是大神的思路啊？仍是有其它鲜为人知的推导过程，我没看过，因此孤陋寡闻了~

设：

$i=p_{1}^{r_{1}}p_{2}^{r_{2}}\cdots p_{k}^{r_{k}}$
$i\cdot p_{1}=p_{1}^{r_{1}+1}p_{2}^{r_{2}}\cdots p_{k}^{r_{k}}$
$\frac{i}{p_{1}}=p_{1}^{r_{1}-1}p_{2}^{r_{2}}\cdots p_{k}^{r_{k}}$

则有：
\(\sigma(i)=(1+p_{1}+\cdots+p_{1}^{r_{1}})(1+p_{2}+\cdots+p_{2}^{r_{2}})\cdots(1+p_{k}+\cdots+p_{k}^{r_{k}})\)

\(\sigma(i\cdot p_{1})=(1+p_{1}+\cdots+p_{1}^{r_{1}+1})(1+p_{2}+\cdots+p_{2}^{r_{2}})\cdots(1+p_{k}+\cdots+p_{k}^{r_{k}})\)

\(\sigma(\frac{i}{p_{1}})=(1+p_{1}+\cdots+p_{1}^{r_{1}-1})(1+p_{2}+\cdots+p_{2}^{r_{2}})\cdots(1+p_{k}+\cdots+p_{k}^{r_{k}})\)

同理：设后面的那一大串为\(T\).
即：\(T=(1+p_{2}+\cdots+p_{2}^{r_{2}})\cdots(1+p_{k}+\cdots+p_{k}^{r_{k}})\)
则：
\(\sigma(i)=(1+p_{1}+\cdots+p_{1}^{r_{1}})T\)
\(\sigma(i\cdot p_{1})=(1+p_{1}+\cdots+p_{1}^{r_{1}+1})T=(1+p_{1}+\cdots+p_{1}^{r_{1}})T + p_{1}^{r_{1}+1}T = \sigma(i)+p_{1}^{r_{1}+1}T\)
\(\sigma(\frac{i}{p_{1}})=(1+p_{1}+\cdots+p_{1}^{r_{1}-1})T=(1+p_{1}+\cdots+p_{1}^{r_{1}-1} +p_{1}^{r_{1}} -p_{1}^{r_{1}})T =(1+p_{1}+\cdots+p_{1}^{r_{1}-1} +p_{1}^{r_{1}})T-p_{1}^{r_{1}}T =\sigma(i)-p_{1}^{r_{1}}T\)

整理一下：
\(\sigma(i)=(1+p_{1}+\cdots+p_{1}^{r_{1}})T\) ①
\(\sigma(i\cdot p_{1})=(1+p_{1}+\cdots+p_{1}^{r_{1}+1})T=\sigma(i)+p_{1}^{r_{1}+1}T\) ②
\(\sigma(\frac{i}{p_{1}})=(1+p_{1}+\cdots+p_{1}^{r_{1}-1})T=\sigma(i)-p_{1}^{r_{1}}T\) ③

为了消元，去掉\(T\),因此③*\(p_1\)+②
获得：
\(\sigma(i\cdot p_{1})+p_{1}\sigma(\frac{i}{p_{1}})=(p_{1}+1)\sigma(i)\)
整理，即：
\(\large{\sigma(i\cdot p_{1})=(p_{1}+1)\sigma(i)-p_{1}\sigma(\frac{i}{p_{1}})}\)

再次费了牛劲，获得了第二个递推关系式！

下面给出代码：

#include <bits/stdc++.h>

using namespace std;
typedef long long LL;

//用途：线性筛约数和的代码
const int N = 1e6 + 10;
int primes[N];  //质数数组
int idx;        //质数数组下标游标
bool st[N];     //是否被筛出
int sigma[N];   //约数和数组

void getDivisorSum(int n) {
    sigma[1] = 1; //σ(1)=1,由于1的约数只有1，约数和就是1，这是递推的起点
    for (int i = 2; i <= n; i++) { //倍数
        if (!st[i])primes[++idx] = i, sigma[i] = i + 1; //若是是质数，那么放入到质数数组中，而且质数的约数和是i+1
        for (int j = 1; i * primes[j] <= n && j <= idx; j++) {//遍历每个已经出现的质数（欧拉筛大法好！），对primes[j]的整数倍进行筛除,i * primes[j] <= n表示是范围内，超过边界就不用了
            st[i * primes[j]] = true;   //标识为质数
            if (i % primes[j] == 0) {    //i与primes[j]不互质
                //利用费了牛劲才得出的递推公式进行计算，由于是从小到大过来的，依赖项提早算好，因此能够算出来
                sigma[i * primes[j]] = sigma[i] * (primes[j] + 1) - primes[j] * sigma[i / primes[j]];
                break;//都说这个break值钱，只为第一个质数因子筛掉，欧拉筛的精髓，妙！
            }
            //i与primes[j]互质
            sigma[i * primes[j]] = sigma[i] * (primes[j] + 1);//约数和是当前质数+1，由于一共两个约数，一个是1一个是本身，和固然是p+1
        }
    }
}

LL res;

int main() {
    int n = 100;
    //线性筛约数和
    getDivisorSum(n);
    for (int i = 1; i <= n; i++) res += sigma[i];
    cout << res << endl;
    return 0;
}

事实上它还能够再短一点(附上约数个数和约数和放在一块儿的版本)：

#include <bits/stdc++.h>

using namespace std;
typedef long long LL;

/**
 * 功能：线性筛出约数个数与约数和
 * Tag:模板，约数个数，约数
 */
const int N = 1e6 + 10;
int n;
int primes[N];      //质数数组
int idx;            //质数数组下标游标
bool st[N];          //是否已被筛出
int d[N];           //约数个数数组
int sigma[N];       //约数和数组

void get_divisor(int n) {
    //积性函数的出发值
    d[1] = sigma[1] = 1;

    for (int i = 2; i <= n; i++) {  //倍数
        if (!st[i])primes[++idx] = i, d[i] = 2, sigma[i] = i + 1;
        for (int j = 1; i * primes[j] <= n & j <= idx; j++) {
            st[i * primes[j]] = true;
            d[i * primes[j]] = d[i] << 1;
            sigma[i * primes[j]] = sigma[i] * (primes[j] + 1);
            if (i % primes[j] == 0) {
                d[i * primes[j]] -= d[i / primes[j]];
                sigma[i * primes[j]] -= primes[j] * sigma[i / primes[j]];
                break;
            }
        }
    }
}

LL res;

int main() {
    cin >> n;
    //开始筛约数个数，约数和
    get_divisor(n);
    //输出约数个数和
    for (int i = 1; i <= n; i++) res += d[i];
    cout << res << endl;

    //输出约数和
    cout << sigma[n] << endl;

    return 0;
}
相关文章
相关标签/搜索
每日一句
    每一个你不满意的现在，都有一个你没有努力的曾经。
本站公众号
   欢迎关注本站公众号,获取更多信息
联系我们 最近搜索 最新文章 沪ICP备13005482号-10 MyBatis教程 SQL 教程 MySQL教程 Java 教程 Thymeleaf 教程 Hibernate教程 Spring教程 Redis教程