poj2480(利用欧拉函数的积性求解)

题目链接:  http://poj.org/problem?id=2480

题意:∑gcd(i, N) 1<=i <=N,就这个公式,给你一个n,让你求sum=gcd(1,n)+gcd(2,n)+gcd(3,n)+…………gcd(n-1,n)+gcd(n,n),(1<=n<2^31)是多少?

放心吧!!!暴力肯定是做不出来的,如果你数论只会gcd(和我一样),那还是学点东西再来挑战这个题吧!

 

 这个题需要用到欧拉函数的知识……

欧拉函数的定义:对正整数n,欧拉函数是小于n的正整数中与n互质的数的数目(我们定义φ(1)=1

欧拉函数的的通式:φ(n)=n*(1-1/p1)*(1-1/p2)*(1-1/p3)……*(1-1/ps)(p1,p2,p3,……ps均是n的质因子)

欧拉函数有这么几个比较重要的性质:

性质1:如果n是质数p的k次幂,那么φ(n)=p^k-1*(p-1)

性质2:欧拉函数是积性函数——若m,n互质,φ(mn)=φ(n)*φ(m),积性函数和完全积性函数有区别,有兴趣可以自己百度一下

性质3:当n为奇数的时候,φ(2n)=φ(n),这一点是由性质2推出来的,因为2必定和所有的奇数都是互质的,所以而φ(2)=1。所以得出这个结果

性质4:n的所有质因子之和等于φ(n)*n/2(这不算性质,只能算延伸)。

好了,大体介绍完了欧拉函数,我们可以开始来看看这个题怎么做了。

首先要知道gcd(i,n)是积性函数(当n固定时),也就是说gcd(i*j,n)=gcd(i,n)*gcd(j,n)(这里还有一个限制条件,就是i,j互质,所以gcd并非完全积性函数)

一开始我不是利用积性函数的性质做的,但是也用到了欧拉函数,因为一眼就和欧拉函数有关一下就手撕了,相当于半暴力吧,344ms,好丢人,这里也说下是怎么做的。

我们枚举gcd(i,n)的所有情况即n的所有因子都有可能是他和其他数的最大公因数。我们假设M是n与i的最大公因数,那么所有与i互质且小于i的数与M的乘积我们设这个数为j,与n的最大公因数都为M,即gcd(j,n)=M,.这里所有与i互质且小于i的数也就是i的欧拉函数φ(i)而i=n/M。

但是我们直接1-n去枚举n的所有因子设为M,来枚举最大公约数的所有情况答案需要n的复杂度,在2^31次方的情况下是会超时的,所以我们采用折半枚举,具体看代码吧

//Author: xiaowuga
#include<iostream>
#include<cmath>
#define maxx INT_MAX
#define minn INT_MIN
#define inf 0x3f3f3f3f
const long long N=1000;
using namespace std;
typedef long long LL;
LL euler(LL n){
    LL res=n;
    for(LL i=2;i*i<=n;i++){
        if(n%i==0){
            res=res/i*(i-1);
            while(n%i==0) n/=i;
        }
    }
    if(n>1)  res=res/n*(n-1);
    return res;
}
int main(){
    LL n;
    while(scanf("%lld",&n)!=EOF){
        LL ans=0;
        for(LL i=1;i*i<=n;i++){
            if(n%i==0){
                ans+=euler(n/i)*i;
                if(i*i<n) ans+=euler(i)*(n/i);
            }
        }
        printf("%lld
",ans);
    }
    return 0;
}

好了现在我们需要来学习真正的姿势了,我也是刚学的,利用gcd是积性函数的性质,根据前文说的,我们有这样的结论:n>1时 n=p1^a1*p2^a2*...*ps^as,p为n的质因子,那么f(n)是积性函数的充要条件是f(1)=1,及f(n) = f(p1^a1)*f(p2^a2)*...f(pr^ar),所以只要求f(pi^ai)就好。现在来看具体做法。

f(pi^ai) =  Φ(pi^ai)+pi*Φ(pi^(ai-1))+pi^2*Φ(pi^(ai-2))+...+pi^(ai-1)* Φ(pi)+ pi^ai *Φ(1)

根据性质1,我们可以做出如下化简

f(pi^ai)=[pi^(ai-1)*(pi-1) ] +  [pi*pi^(ai-2)*(pi-1)]  +  [pi^2*pi^(ai-3)*(pi-1)]  +  [pi^3*pi^(ai-4)*(pi-1)]....[pi^(ai-1)*pi^(ai-ai)*(pi-1)]+pi^ai   ①

然后对①提取公因子(pi-1)

f(pi^ai)=(pi-1){[pi^(ai-1) ] +  [pi*pi^(ai-2)]  +  [pi^2*pi^(ai-3)]  +  [pi^3*pi^(ai-4)]....[pi^(ai-1)*pi^(ai-ai)]+[pi^ai/(pi-1)]}  

紧接着我们发现出了最后一项每个[]每个方括号内乘积都等于pi^(ai-1),所以对②提取公因子pi^(ai-1)

f(pi^ai)=(pi-1)*pi^(ai-1)*{ai+[pi/(pi-1)]} 

然后把(pi-1)/pi放进括号里得

f(pi^ai)=pi^(ai)*{1+ai*(pi-1)/pi

这个 是单个f(pi^ai)的公式,我们提取所有的pi^(ai)相乘实际上就是n了,所以我们可以得到f(n)的公式:f(n)=n*∏(1+ai*(pi-1)/pi

然后我们看代码吧!

//Author: xiaowuga
#include<iostream>
#include<cmath>
#define maxx INT_MAX
#define minn INT_MIN
#define inf 0x3f3f3f3f
const long long N=1000; 
using namespace std;
typedef long long LL;
int main() {
    ios::sync_with_stdio(false);cin.tie(0);
    int n;
    while(cin>>n){
        LL i,sqr,p,a,ans;
        ans=n;
        sqr=floor(sqrt(1.0*n));
        for(int i=2;i<=sqr;i++){
            if(n%i==0){
                a=0;p=i;
                while(n%p==0){
                    a++;n/=p;
                }
                ans=ans+ans*a*(p-1)/p;
            }
        }
        if(n!=1) ans=ans+ans*(n-1)/n;
        cout<<ans<<endl;
    }        
    return 0;
}
原文地址:https://www.cnblogs.com/xiaowuga/p/7161513.html