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

读贾志鹏线性筛有感

时间:2015-05-04 21:59:57      阅读:331      评论:0      收藏:0      [点我收藏+]

标签:

先拜大牛。感谢贾志鹏严谨的思维。以及简单清晰的论文描述。

一定要结合论文看。我只是提出我觉得关键的部分。论文在网上随处可见。贾志鹏线性筛。

 

开头两种线性筛的比较。

一种是传统的线性筛。时间复杂度为N*log(log(N))。

另外一种是优化了合数的筛法。文中称作Euler线性筛。

其优化的地方。

举个例子:合数6。 是2的倍数也是3的倍数。 当你用传统的筛法的时候在遍历2的倍数的时候会遍历到6。遍历3的倍数的时候同样也会遍历到6。

而另外一种只会筛出6为2的倍数。3就不会筛6了。

另外个人认为筛法二有一个很重要的思想。当i为合数的时候。其实脑海里不认为是合数。而是素数的乘积。这样就能比较直观地确定这个算法的正确性了。

 

积性函数。

分为完全积性和条件积性。

我们最喜欢的积性。大概就是互素积性了。因为满足互素积性的话。根据算术基本定理。就能够简单做到推广到任意实数。

f(1) = 1 。 这个在我们高中数学题。抽象函数。就已经能简单知道了。

 

欧拉函数。就不再谈了。包括其线性筛的那一步至关重要的证明。也在我其它博文提到过了。

其 欧拉定理和费马小定理的作用。我得再多补充一点。

以及约数和。 n的约数和为 n*φ(n)/2.

 

莫比乌斯函数和容斥定理的关系。

可以发现莫比乌斯函数其实就是容斥定理的映射一般。

莫比乌斯函数 是我们再熟悉不过的了。

首先看 (-1)^r m = p1p2p3p4p5pr 其实就是在模拟容斥定理。

假如一但不是素数。那就为0.

 

两个函数的线性筛。这其实是我们处理问题的基本

问题一:求1~N对质数P的乘法逆元。

关于f(n)为完全积性函数。根据同余定理可以简单获得。要证明的话。减法证同余即可。

P = nt + k

n‘  ≡ n*(t^2)*f(k)^2 (mod P)

这个证明过程很漂亮(很佩服这么顺畅,思维这么清晰)。也是根据同余定理。还有逆元的性质。就能推理的。

这个问题的意义。可以求N!的 mod P 的逆元了。逆元还是很有用的。因为毕竟除法并没有特别好的同余式。(依稀还记得那两个。)

问题二:给T组N,M.依次求出技术分享的值.(N,M<=10^6,T<=10^3)

求解gcd(a,b).把gcd(a,b)当做n.再通过欧拉函数和式。推导过程如下。

技术分享

第二个等式是用d来看待式子的方法来化简和式的。

之后再穷举d即可。论文上有个优化在这里没体现。

 

 

技术分享
#include<stdio.h>
#include<string.h>
#define    N 100

bool mark[N+5];
int prime[N+5];
int num;
int euler[N+5];
int Min(int a,int b){return a>b?a:b;}
void Euler()
{
    int i,j;
    num = 0;
    memset(euler,0,sizeof(euler));
    memset(prime,0,sizeof(prime));
    memset(mark,0,sizeof(mark));
    euler[1] = 1; // multiply function
    for(i=2;i<=N;i++)
    {
        if(!mark[i])
        {    
            prime[num++] = i;
            euler[i] = i-1;
        }
        for(j=0;j<num;j++)
        {
            if(i*prime[j]>N){break;}
            mark[i*prime[j]] = 1;
            if(i%prime[j]==0)
            {
                euler[i*prime[j]] = euler[i]*prime[j];
                break;
            }
            else
            {
                euler[i*prime[j]] = euler[prime[j]]*euler[i];
            }
        }
    }
}

int main()
{
    int i;
    int M1,M2;
    Euler();
    for(i=0;i<num;i++)printf("%d ",prime[i]);
    printf("\n");
    for(i=1;i<=N;i++)printf("%d ",euler[i]);
    printf("\n");
    M1 = 2;
    M2 = 3;
    int sum = 0;
    int min = Min(M1,M2);
    for(i=1;i<=min;i++)
    {
        sum += euler[i]*(M1/i)*(M2/i);
    }
    printf("%d\n",sum);
}
the second problem test

 

 

 

问题三:给T组N,M.依次求出技术分享的值.(N,M<=10^6,T<=10^3)

在证明之前,先证明以下式子。

 技术分享

 

问题的解决推导。

技术分享

第一个等式。lcm(a,b) = ab/gcd(a,b).

第二个等式。令d=gcd(a,b)。

第三个等式。转化为d的视角。(这个手法经常有)。

第四个等式。转化为莫比乌斯函数。

第五个等式。利用上述的等式来转化。注意d和d‘

第六个等式。论文中提到的有趣的化简性质。

第七个等式。至今未明白。

之后就是如论文中介绍的。g(d) 为积性函数。线性筛之。

总体上算法还是N的。

问题四:给T组N,M.依次求出技术分享的值.(N,M<=10^6,T<=10^3)

实质上就是求 技术分享其中x和y互素。的对数。

我们是时候需要有和式化成的思想了。[gcd(a,b)=1]真是漂亮的莫比乌斯函数的和式的结果。

技术分享

 

 

技术分享
#include<stdio.h>
#include<string.h>
#define    N 100

bool mark[N+5];
int prime[N+5];
int num;
int mobi[N+5];
int Min(int a,int b){return a>b?a:b;}
void Mobi()
{
    int i,j;
    num = 0;
    memset(mobi,0,sizeof(mobi));
    memset(prime,0,sizeof(prime));
    memset(mark,0,sizeof(mark));
    mobi[1] = 1; // multiply function
    for(i=2;i<=N;i++)
    {
        if(!mark[i])
        {    
            prime[num++] = i;
            mobi[i] = -1;
        }
        for(j=0;j<num;j++)
        {
            if(i*prime[j]>N){break;}
            mark[i*prime[j]] = 1;
            if(i%prime[j]==0)
            {
                mobi[i*prime[j]] = 0;
                break;
            }
            else
            {
                mobi[i*prime[j]] = mobi[prime[j]]*mobi[i];
            }
        }
    }
}

int main()
{
    int i;
    int M1,M2;
    Mobi();
    for(i=0;i<num;i++)printf("%d ",prime[i]);
    printf("\n");
    for(i=1;i<=N;i++)printf("%d ",mobi[i]);
    printf("\n");
    M1 = 2;
    M2 = 3;
    int sum = 0;
    int min = Min(M1,M2);
    for(i=1;i<=min;i++)
    {
        sum += mobi[i]*(M1/i)*(M2/i);
    }
    printf("%d\n",sum);
}
Test problem forth

 

 

 

 

问题五:给T组N.依次求出技术分享的值.(N<=10^6,T<=10^3)

 

其实根据问题三.可以直接获得该化简出来的式子的。

技术分享

然后解法和问题三一样。

但是论文上寻找积性f(n)直接筛出答案。

首先佩服一下其思维的紧密。一个变量啊。就寻找积性函数。这个转化真是清晰而又巧。

技术分享

画个图就能知道 -n 是用来去重复的统计的。

技术分享.

f(n)是积性的。具体证明如论文上解释。

技术分享

第一个等式:d = gcd(n,i)

第二个等式:k = i/d.且全部都除以d.gcd(a,b)=d转化成求互素(gcd(a,b)=1)的问题。

第三个等式:令d=n/d。是对应的。  其实在第二个等式就能看出是欧拉函数求约数和问题了。

第四个等式:不解释了吧。

第五个等式:手算一下容易得。

读贾志鹏线性筛有感

标签:

原文地址:http://www.cnblogs.com/Milkor/p/4474835.html

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