线性筛求 约数个数 与 约数和

线性筛求 约数个数约数和

线性筛,顾名思义,就是欧拉筛,在线性时间内可以求出答案,也就是O(N)的时间,非常牛X的效率。

一、约数个数

根据数字唯一分解定理,设

n=p1r1p2r2p3r3...pkrk

对于每个n的约数,一定是由以上质因数p1pk相乘而来,根据乘法原理,每个质因数都可以选择0riri+1个选择。

n的约数个数即为

d(n)=(r1+1)×(r2+1)×...×(rk+1)=i=1k(ri+1)

筛的过程中需要保存n的最小质因子的出现个数即r1,原因后面会说明:

  • d[i]记录i的约数个数
  • num[i]记录i的最小质因数个数

分情况讨论:

(一)、i是质数

{d[i]=2num[i]=1

  • d[i]:i是质数,所以约数只有1和自身,约数个数为2
  • num[i]:i是质数,最小质因子就是自己,个数是1

(二)、i取模枚举的第j个质数为0primes[j]整除i

Q1:我们在研究哪些数字之间的递推关系?
A:既然是递推,一般都是已知小的推出大的,这里应该是已知d[i]推出d[iprimes[j]]的值。

Q2:d[iprimes[j]]是啥?
A:
primes[j]i的因子,iprimes[j]唯一分解不会产生新的质因子,依然是p1,p2,...pk
primes[j]从小到大枚举,所以一定是i的最小质因子!所以p1r1+1,按唯一分解定理展开:

d[iprimes[j]]=(1+r1+1)......(1+rk)

Q3:d[iprimes[j]]怎么从d[i]转移?

这个时候就可以用到我们之前维护的num[i]了。

转移也就非常简单了:

d[iprimes[j]]=d[i]/(num[i]+1)(num[i]+2)

解释:

{d(i)=(1+r1)......(1+rk)d(iprimes[j])=(1+r1+1)......(1+rk)

使用瞪眼大法(观察法),找出两者之间的递推关系,很明显,就是 ①式除以1+r1,再乘上1+r1+1即可转化为 ②式,这下明白为什么要记录num[i]最小质因子的个数了吧~

Q4:那么num[iprimes[j]]是不是也需要从num[i]进行转移呢?
A:num[iprimes[j]]也要转移,加上最小质因子primes[j]的贡献也就是:

num[iprimes[j]]=num[i]+1

(三)、当前数取模枚举的第j个质数不为0,即 primes[j] 不能整除 i

首先我们知道iprimes[j]这个数中之前一定不包含primes[j]这个质因数,

那么约数个数要加上prime[j]的,也就是:

d[iprimes[j]]=(1+r1)...(1+rk)(1+1)=d[i]2=d[i]d[primes[j]]

然后对于最小质因子,因为j是从小到大枚举的,所以iprimes[j]这个数的最小质因子也就是primes[j]

所以就可以得到:

num[iprimes[j]]=1

综上,就可以写出筛质因数个数的代码了:

#include <iostream>
#include <cstdio>
#include <cstring>
using namespace std;
const int N = 1010;

int primes[N], cnt; // primes[]存储所有素数
bool st[N];         // st[x]存储x是否被筛掉
int d[N];           // d[x]表示x的约数个数
int num[N];         // num[x]表示x的最小质因数的个数
int n;

//欧拉筛法+求约数个数
void get_primes(int n) {
    d[1] = 1; // 1的约数只有1个,这个比较特殊

    for (int i = 2; i <= n; i++) {
        if (!st[i]) {
            primes[cnt++] = i;
            // i是质数
            d[i] = 2;   //约数个数是2个,一个是1,另一个是i
            num[i] = 1; //最小质因子个数是1,最小质因子就是自己i
        }

        for (int j = 0; i * primes[j] <= n; j++) {
            st[i * primes[j]] = true;
            if (i % primes[j] == 0) {
                d[i * primes[j]] = d[i] / (num[i] + 1) * (num[i] + 2);
                num[i * primes[j]] = num[i] + 1;
                break;
            } else {
                // d[i * primes[j]] = d[i] * d[primes[j]]; 等价于下面的代码 
                d[i * primes[j]] = d[i] * 2;
                num[i * primes[j]] = 1;
            }
        }
    }
}

int main() {
    scanf("%d", &n);
    get_primes(n);
    //输出1~n之间所有数字的约数个数
    for (int i = 1; i <= n; i++) printf("%d %d\n", i, d[i]);
    return 0;
}

二、约数和

sd[i]表示i的约数和

根据算数基本定理得:约数和参考链接

sd[n]=(p10+p11+p12++p1r1)(p20+p21+p22++p2r2)(pk0+pk1+pk2++pkrk)

那么根据这个式子就可以开始干了。。。


num[i]表示最小质因子的那一项:

num[i]=p10+p11+p12++p1r1

分情况讨论:

(一)、i是质数

sd[i]=num[i]=i+1

解释:

  • i是质数,它的约数只有自己和1,约数和是i+1
  • i是质数,质因子分解也只能分解出自己,sd[i]=i0+i1=1+i

(二)、i取模枚举的质数primes[j]不为0,即primes[j]不能整除i

注:下面的讨论中使用pj代替primes[j]

ipj里原先没有pj这一项,加上这一项之后:

sd[i]=(p10+p11+p12++p1r1)(p20+p21+p22++p2r2)(pk0+pk1+pk2++pkrk)

sd[ipj]=(p10+p11+p12++p1r1)(p20+p21+p22++p2r2)(pk0+pk1+pk2++pkrk)(pj0+pj1)

pj 是质数 sd[pj]=pj0+pj1

sd[ipj]=sd[i]sd[pj] 注:积性函数

同时更新一下num[ipj]

num[ipj]=num[pj]=pj0+pj1=1+pj

解释:
因为质因子从小到大枚举,那么ipj的最小质因子就应该是pj,那么num[ipj]也就应该等于num[pj]

(三)、i取模枚举的质数j为0,即 primes[j]整除i

那么sd[ipj]中的第一项也就是num[ipj]也一定是pj的第一项

解释:(等比数列变形小技巧)
1+pi+pi2++piri这个时候要变成1+pi+pi2++piri+piri+1,那么只需要所有的都乘以一个pi,然后再加一个1就好了。

(1+pi+pi2++piri)pi+1=1+pi+pi2++piri+piri+1

结论:

num[ipj]=num[i]pj+1

这样,我们又可以开始愉快的写代码啦...

#include <bits/stdc++.h>
using namespace std;
typedef long long LL;

const int N = 1e3 + 10;

int primes[N], cnt;
bool st[N];

int num[N]; //最小质因子p1组成的等比序列 p1^0+p1^1+...+p1^r1
int sd[N]; //约数和

int n;

void get_primes(int n) {
    sd[1] = 1; // 1的约数只有自己,约数和是1

    for (int i = 2; i <= n; i++) {
        if (!st[i]) {
            primes[cnt++] = i;
            //质数
            sd[i] = num[i] = i + 1; // 1和自身是约数,约数和为i+1; 同时因为i是质数,所以只有一个分拆项,并且分拆项内容=pj^0+pj^1=1+pj=1+i
        }
        for (int j = 0; primes[j] * i <= n; j++) {
            st[i * primes[j]] = true;

            if (i % primes[j]) {
                sd[i * primes[j]] = sd[i] * sd[primes[j]]; //积性函数
                num[i * primes[j]] = primes[j] + 1;
            } else {
                sd[i * primes[j]] = sd[i] / num[i] * (num[i] * primes[j] + 1);
                num[i * primes[j]] = num[i] * primes[j] + 1;
                break;
            }
        }
    }
}

int main() {
    scanf("%d", &n);
    get_primes(n);
    //约数和要不要自己本身,这里是不想要自己本身
    for (int i = 1; i <= n; i++) printf("%d ", sd[i] - i);
    puts("");
    return 0;
}
posted @   糖豆爸爸  阅读(933)  评论(0编辑  收藏  举报
相关博文:
阅读排行:
· 被坑几百块钱后,我竟然真的恢复了删除的微信聊天记录!
· 没有Manus邀请码?试试免邀请码的MGX或者开源的OpenManus吧
· 【自荐】一款简洁、开源的在线白板工具 Drawnix
· 园子的第一款AI主题卫衣上架——"HELLO! HOW CAN I ASSIST YOU TODAY
· Docker 太简单,K8s 太复杂?w7panel 让容器管理更轻松!
历史上的今天:
2021-09-02 等比数列求和公式
2017-09-02 访问Github慢的解决办法
2017-09-02 CentOS6.9下安装 Pika 2.2.5(新增了拷贝安装版本的办法+对于PID的位置及数据库位置的理解)
2016-09-02 鞋业管理系统定期执行任务
2016-09-02 针对各地项目icomet停止服务的临时处理办法
Live2D
点击右上角即可分享
微信分享提示