[数论-欧拉函数的应用]NEFU 1115

description

在高数课上老师给出了一个等式GCD(n-a,n)*GCD(n-b,n)=n^k;

已知,n,k，问同学们有多少组(a,b)满足这个等式。结果对10^9+7取模。

1<=a,b<=n

  

只需考虑k=1的情况。

设g= GCD(n-a,n)。

g*k1=n-a g*k2=n

k1与k2互质，且k1<k2。

  

选择一个身为n的因子g= GCD(n-a,n)，有多少种a可以选择呢？也就是说有多少个k1可以满足上述条件。

这不就是欧拉函数吗，，哈哈!!!，也就是说a共有euler(n/g)种选择。

  

现在对于g1*g2=n，共有euler(g2)* euler (g1)种选择。

  

因为euler(n) = n(1-1/p1)(1-1/p2)……(1-1/pr)，可以想到 euler(g2)* euler(g1)其中包括g1*g2，还有n的质因数，这些质因数有些被乘了两次。所以对于g1，g2，只需要找出他们的最大公约数再巴拉巴拉。。。。，一些细节要考虑啊。

#include <stdio.h> 
#include<cstring> 
using namespace std; 
const int MAXN=10000; 
int prime[MAXN+1]; 
void getPrime() 
{ 
    memset(prime,0,sizeof(prime)); 
    for(int i=2; i<=MAXN; i++) 
    { 
        if(!prime[i]) prime[++prime[0]]=i; 
        for(int j=1; j<=prime[0]&&prime[j]<=MAXN/i; j++) 
        { 
            prime[prime[j]*i]=1; 
            if(i%prime[j]==0) break; 
        } 
    } 
} 
long long factor[100][2]; 
int fatCnt; 
int getFactors(long long x) 
{ 
    fatCnt=0; 
    long long tmp=x; 
    for(int i=1; prime[i]<=tmp/prime[i]; i++) 
    { 
        factor[fatCnt][1]=0; 
        if(tmp%prime[i]==0) 
        { 
            factor[fatCnt][0]=prime[i]; 
            while(tmp%prime[i]==0) 
            { 
                factor[fatCnt][1]++; 
                tmp/=prime[i]; 
            } 
            fatCnt++; 
        } 
    } 
    if(tmp!=1) 
    { 
        factor[fatCnt][0]=tmp; 
        factor[fatCnt++][1]=1; 
    } 
    return fatCnt; 
} 
   
   
   
int n,k; 
long long euler[40001]; 
void getEuler() 
{ 
    memset(euler,0,sizeof(euler)); 
    euler[1] = 1; 
    for(int i = 2; i <= 40000; i++) 
        if(!euler[i]) //i是素数 
            for(int j = i; j <= 40000; j += i) 
            { 
                if(!euler[j]) 
                    euler[j] = j; 
                euler[j] = euler[j]/i*(i-1); 
            } 
} 
int gcd(int a,int b) 
{ 
    if(a==0) return b; 
    return gcd(b%a,a); 
} 
int main() 
{ 
    getPrime(); 
    getEuler(); 
    while(scanf("%d%d",&n,&k)!=EOF) 
    { 
        getFactors(n); 
        long long ret = n; 
        for(int i = 0; i < fatCnt; i++) 
        { 
            ret = ret/factor[i][0]*(factor[i][0]-1); 
        } 
        long long ans=(ret*2)%1000000007; 
        if(k>2) puts("0"); 
        else if(k==2) puts("1"); 
        else if(k==1) 
        { 
            int i; 
            for(i=2; i*i<n; i++) 
            { 
                if(n%i==0) 
                { 
                    //printf("%d %d
",i,n/i); 
                    //printf("%d %d
",n/i,i); 
                    ans+=ret*euler[gcd(i,n/i)]/gcd(i,n/i)*2; 
                    if(ans>1000000007) ans%=1000000007; 
                } 
            } 
            if(i*i==n) 
            { 
                if(n%i==0) 
                { 
                    //printf("%d %d
",i,n/i); 
                    ans+=ret*euler[gcd(i,n/i)]/gcd(i,n/i); 
                    if(ans>1000000007) ans%=1000000007; 
                } 
            } 
            //ans++;ans%=1000000007; 
            printf("%lld
",ans); 
        } 
    } 
    return 0; 
}