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\n",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;
}

 


免責聲明!

本站轉載的文章為個人學習借鑒使用,本站對版權不負任何法律責任。如果侵犯了您的隱私權益,請聯系本站郵箱yoyou2525@163.com刪除。



 
粵ICP備18138465號   © 2018-2025 CODEPRJ.COM