码迷,mamicode.com
首页 > 其他好文 > 详细

线性筛求约数个数以及约数和

时间:2020-09-24 21:27:32      阅读:46      评论:0      收藏:0      [点我收藏+]

标签:++i   time   筛法   个数   ons   max   线性   long   ngx   

RT 考试的时候口胡出来的,正确性不会证,不过貌似100000内的数都是对的(现在已经会了,在此鸣谢gyz大佬)

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int maxn=5e7+10,maxm=3000000;
int a[maxn],prime[maxm],f[maxn];
bool is_not_prime[maxn];
int n,mod,tot;
ll ans;
void Solve(){
	scanf("%d%d",&n,&mod);
	is_not_prime[0]=is_not_prime[1]=f[1]=1;
	for(register int i=2;i<=n;++i){
		if(!is_not_prime[i]){
			prime[++tot]=i;
			f[i]=2;
		}
		for(register int j=1;j<=tot&&prime[j]*i<=n;++j){
			is_not_prime[i*prime[j]]=1;
			if(i%prime[j]!=0) f[i*prime[j]]=f[i]*f[prime[j]];
			else{
                              f[i*prime[j]]=f[i]*2-f[i/prime[j]];
                              break;
			}
		}
	}
	
}
int main(){
	freopen("in","r",stdin);
//	freopen("out2","w",stdout);
	Solve();
	return 0;
}

关于正确性证明(以下证明来自gyz大佬):
首先根据线性筛原理,当 \(prime[j]\) 可以整除 \(i\) 时,\(prime[j]\)\(i\) 的最小质因子。
为了方便,用\(p\)代替\(prime[j]\)
\(i=p^k\times p_1^{k_1}\times p_2^{k_2}\times\)……\(\times p_n^{k_n}\)
根据唯一分解定理,有\(d[i]=(k+1)\times (k_1+1)\times (k_2+1)\times\)……\(\times (k_n+1)\)
不妨设\(t=(k_1+1)\times (k_2+1)\times\)……\(\times (k_n+1)\)
则有:
\(d[i]=(k+1)\times t\)
\(d[i/p]=k\times t\)
\(d[i\times p]=(k+2)\times t\)
所以
\(d[i\times p]=2\times (k+1)\times t-k\times t=d[i]\times 2-d[i/p]\)
和传统写法相比,少开一个数组记录每个数的最小质因子出现的次数。同时也省去一些对这个数组的维护。
oiwiki上的筛法:

void pre() {
  d[1] = 1;
  for (int i = 2; i <= n; ++i) {
    if (!v[i]) v[i] = 1, p[++tot] = i, d[i] = 2, num[i] = 1;
    for (int j = 1; j <= tot && i <= n / p[j]; ++j) {
      v[p[j] * i] = 1;
      if (i % p[j] == 0) {
        num[i * p[j]] = num[i] + 1;
        d[i * p[j]] = d[i] / num[i * p[j]] * (num[i * p[j]] + 1);
        break;
      } else {
        num[i * p[j]] = 1;
        d[i * p[j]] = d[i] * 2;
      }
    }
  }
}

所以这份代码里面的num数组是可以省去的

线性筛求约数和:

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int maxn=5e7+10,maxm=3000000;
int a[maxn],prime[maxm];
ll num[maxn],f[maxn];
bool is_not_prime[maxn];
int n,mod,tot;
ll ans;
void Solve(){
	scanf("%d%d",&n,&mod);
	is_not_prime[0]=is_not_prime[1]=f[1]=1;
	for(register int i=2;i<=n;++i){
		if(!is_not_prime[i]){
			prime[++tot]=i;
			f[i]=num[i]=i+1;
		}
		for(register int j=1;j<=tot&&prime[j]*i<=n;++j){
			is_not_prime[i*prime[j]]=1;
			if(i%prime[j]!=0){
				f[i*prime[j]]=f[i]*f[prime[j]];
				num[i*prime[j]]=prime[j]+1;
			}
			else{
				f[i*prime[j]]=f[i]/num[i]*(num[i]*prime[j]+1);
               	                num[i*prime[j]]*=prime[j];
               	                num[i*prime[j]]++;
               	                break;
			}
		}
	}
	
}
int main(){
	freopen("in","r",stdin);
//	freopen("out2","w",stdout);
	Solve();
	return 0;
}

(写的貌似没错吧。。。。没调试)
证明见这位神犇的博客
暂时没想出来去掉这份代码里面num数组的方法

线性筛求约数个数以及约数和

标签:++i   time   筛法   个数   ons   max   线性   long   ngx   

原文地址:https://www.cnblogs.com/wwcdcpyscc/p/13662918.html

(0)
(0)
   
举报
评论 一句话评论(0
登录后才能评论!
© 2014 mamicode.com 版权所有  联系我们:gaon5@hotmail.com
迷上了代码!