博客
关于我
强烈建议你试试无所不能的chatGPT,快点击我
线性时间 筛素数,求前n个数的欧拉函数值,求前n个数的约数个数
阅读量:7072 次
发布时间:2019-06-28

本文共 3821 字,大约阅读时间需要 12 分钟。

转自 

利用积性函数的优化.

这个文章主要介绍了3算法

1线性时间筛素数

2线性时间求前n个数的欧拉函数值

3线性时间求前n个数的约数个数

 

一、   首先介绍下积性函数。

 

下面是wiki的条目:

 

在非数论的领域,积性函数指有对于任何a,b都有性质f(ab)=f(a)f(b)的函数。

 

在数论中的积性函数。对于正整数n的一个算术函数f(n),当中f(1)=1且当a,b互质,f(ab)=f(a)f(b),在数论上就称它为积性函数。

若某算术函数f(n)符合f(1)=1,且就算a,b不互质,f(ab)=f(a)f(b),称它为完全积性的。

 

例子

φ(n) -欧拉φ函数,计算与n互质的正整数之数目

μ(n) -默比乌斯函数,关于非平方数的质因子数目

gcd(n,k) -最大公因子,当k固定的情况

d(n) -n的正因子数目

σ(n) -n的所有正因子之和

σk(n): 因子函数,n的所有正因子的k次幂之和,当中k可为任何复数。在特例中有:

σ0(n) = d(n) 及

σ1(n) = σ(n)

1(n) -不变的函数,定义为 1(n)=1 (完全积性)

Id(n) -单位函数,定义为 Id(n)=n (完全积性)

Idk(n) -幂函数,对于任何复数、实数k,定义为Idk(n) = nk (完全积性)

Id0(n) = 1(n) 及

Id1(n) = Id(n)

ε(n) -定义为:若n = 1,ε(n)=1;若n > 1,ε(n)=0。有时称为“对于狄利克雷回旋的乘法单位”(完全积性)

(n/p) -勒让德符号,p是固定质数(完全积性)

λ(n) -刘维尔函数,关于能整除n的质因子的数目

γ(n),定义为γ(n)=(-1)ω(n),在此加性函数ω(n)是不同能整除n的质数的数目

所有狄利克雷特性均是完全积性的

 

 

二、再介绍下线性筛素数方法

 

bool notp[mr];//素数判定
__int64 pr[670000],pn,ans;//pr存放素数,pn当前素数个数。
 
void getprime()
{
    pn=0;
    memset(notp,0,sizeof(notp));
    for(int i=2;i<mr;i++)
    {
        if(!notp[i])pr[pn++]=i;
        for(int j=0;j<pn && pr[j]*i<mr;j++)
        {
            notp[pr[j]*i]=1;
            if(i%pr[j]==0)break;
        }
    }
}

 

利用了每个合数必有一个最小素因子

每个合数仅被它的最小素因子筛去正好一次。所以为线性时间。

代码中体现在:

if(i%pr[j]==0)break;

pr数组中的素数是递增的,当i能整除pr[j],那么i*pr[j+1]这个合数肯定被pr[j]乘以某个数筛掉。

因为i中含有pr[j],pr[j]比pr[j+1]小。接下去的素数同理。所以不用筛下去了。

在满足i%pr[j]==0这个条件之前以及第一次满足改条件时,pr[j]必定是pr[j]*i的最小因子。

 

 

三、结合线性筛素数算法的优化算法

基于这个线性筛素数算法,我们可以很容易地得到某个数的最小素因子。

因为当i%pr[j]!=0的时候,最小素因子pr[j]与i互质,满足积性函数的条件,可以直接得到f(i*pr[j])=f(i)*f(pr[j]).

不过当i%pr[j]==0时我们必须根据该积性函数本身的特性进行计算.或者在筛的同时保存并递推些附加信息.总之要O(1)求得f(i*pr[j])及完成递推附加信息.

 

下面的两个例子是欧拉函数phi和约数个数.这两个是最常用和最有优化价值的。

利用上面的性质都可以很容易地把前n个用O(n)时间推出来.

当然,利用这个性质还可以对其他积性函数进行优化,这里仅介绍两个常用和有优化价值的。

 

1)欧拉函数(phi)

传统的算法:

对于某素数p且p|n(n能整除p)

if( (n/p) % i == 0 ) phi(n)=phi(n/p)*i;

else phi(n)=phi(n/p)*(i-1);

 

这个传统算法的性质正好用在筛素数算法中.

p为n的最小素因子,当n/p包含该因子p,则phi(n)=phi(n/p)*i;否则phi(n)=phi(n/p)*(i-1);

p即pr[j], n/p即i, n即i*pr[j]了.

 

 

2)约数个数(divnum)

约数不能像phi那么自然,但还是有不错的方法.

约数个数有个性质

divnum(n)=(e1+1)*(e2+1)...(ei表示n的第i个质因数的个数.)

传统方法就是对每个数分解质因数,获得各因数个数再用上式.

 

开一个空间e[i]表示最小素因子的次数

这次说直接点:

筛到i 第j个素数

 

对于divnum

如果i|pr[j] 那么 divnum[i*pr[j]]=divsum[i]/(e[i]+1)*(e[i]+2) //最小素因子次数加1

否则 divnum[i*pr[j]]=divnum[i]*divnum[pr[j]] //满足积性函数条件

 

对于e

如果i|pr[j]  e[i*pr[j]]=e[i]+1; //最小素因子次数加1

否则 e[i*pr[j]]=1; //pr[j]为1次

 

Eular函数代码如下:

#include 
#include
#include
#define MAXN 1000000using namespace std;int p[MAXN+5], pri[1000000], idx = -1;int phi[MAXN+5];void GetPrime(){ for (int i = 2; i <= MAXN; ++i) { if (!p[i]) { // 说明i是一个素数 pri[++idx] = i; } for (int j = 0; j <= idx && pri[j]*i <= MAXN; ++j) {
// 遍历所有的素因子 p[pri[j]*i] = 1; if (i % pri[j] == 0) { // 如果i能够整除pri[j]那么i*pri[j+1]就一定被pri[j]数整除 break; } } } printf("idx = %d\n", idx); for (int i = 0; i <= 100; ++i) { printf("%d ", pri[i]); } puts("");}// 现在计算phi[N=p^e]的取值,由于p是一个素数,我们可以知道在 1 ~ p^e 中的数有 1*p, 2*p, 3*p ... [p^(e-1)]*p // 一共有p^(e-1)个不与p^e互质,所以phi[N] = p^e - p^(e-1) = [p^(e-1)]*(p-1) = << phi[N/p]*p >>// N = p1^e1 * p2^e2 *.... 中 如果 e1 = 1 的话,那么我们由欧拉函数积性直接等于 phi[i/p1] * phi[p1]// 如果 e1 != 1,那么我们就直接利用上面尖括号包括的式子得到 phi[N] = phi[i/p1] * p1void Eular(){ for (int i = 2; i <= MAXN; ++i) { if (!p[i]) { phi[i] = i - 1; // 如果是一个素数,欧拉函数就是 i-1 continue; } for (int j = 0; pri[j] * pri[j] <= i; ++j) { // 进入的 i 一定是一个合数 if (i % pri[j] == 0) { // 如果pri[j]是i的一个素因子 if (i / pri[j] % pri[j] == 0) { // 并且这个素因子还有两个以上的指数 phi[i] = pri[j] * phi[i/pri[j]]; // 那么就直接乘以这个素因子 } else { phi[i] = phi[pri[j]] * phi[i/pri[j]]; } break; } } printf("phi[%d] = %d\n", i, phi[i]); getchar(); } return;}int main(){ GetPrime(); Eular();// system("pause"); return 0; }

 

 

转载地址:http://ktzml.baihongyu.com/

你可能感兴趣的文章
AIR文件操作(二):使用文件对象操作文件和目录
查看>>
mongodb的删除语句
查看>>
Linux 实用技巧
查看>>
hadoop-
查看>>
团队博客
查看>>
4.Struts2中Action的三种访问方式
查看>>
delphi 10.3 IOS中英文错位
查看>>
初识RabbitMQ
查看>>
“Linux内核分析”实验二
查看>>
WCF使用MSMQ通信
查看>>
前后端分离djangorestframework——权限组件
查看>>
redis去重方案
查看>>
内连接和外连接
查看>>
RT-Thread--内核基础
查看>>
PC键盘在Mac下Command/Option键切换
查看>>
数字签名和验签的详细过程
查看>>
漫谈《信号与系统》
查看>>
POJ 1742 Coins(多重背包,优化)
查看>>
内容不随模态框滚动
查看>>
Flume+Kafka+SparkStreaming+Hbase+可视化(二)
查看>>