POJ2480《Longge's problem》题解

Description

Longge is good at mathematics and he likes to think about hard mathematical problems which will be solved by some graceful algorithms. Now a problem comes: Given an integer N(1 < N < 2^31),you are to calculate ∑gcd(i, N) 1<=i <=N. 

"Oh, I know, I know!" Longge shouts! But do you know? Please solve it. 

Input

Input contain several test case. 
A number N per line. 

Output

For each N, output ,∑gcd(i, N) 1<=i <=N, a line

Sample Input

2
6

Sample Output

3
15

大概只有我这种人不会做。。。

看了lyd大神（和我的名字好像。。。）的题解，有些懂了： 大神链接。

主要考察一个叫积性函数的东西：f(x*y)=f(x)*f(y)，中学数学经常见到。

积性函数还有个性质就是积性函数的和也是积性函数。

可以得到一个式子f(n)=f(p1^a1)*f(p2^a2)*f(p3^a3)*...*f(pk^ak)。不过目前我们不知道这个式子有什么用，6666666666666666。

再看本题，如果n和m互质，那么gcd(i,n*m)=gcd(i,n)*gcd(i,m)。∴gcd(i,n)是积性函数==>Σgcd(i,n)是积性函数。

设f(n)=Σgcd(i,n)，问题转换为求出所有f(pi^ai)。

下面来求f(pi^ai)：

首先明确，如果p是n的约数，那么满足gcd(i,n)==p的i的个数是Φ(n/p)。

证明：gcd(i,n)==p。设i=k*p,n=m*p。

则gcd(k,m)=1，也就是k和m互质，要求出满足条件的i的个数，就是求出i所对应的k的个数，即求m的欧拉函数，m=n/p，所以满足条件的i的个数就是Φ(m)=Φ(n/p)。

好了，用这种原则来求f(pi^ai)，就是枚举pi^ai的每个约数(其实就是pi^t,t<=ai)，然后求使gcd(i,pi^ai)==pi^t满足的i的个数，所以 个数*pi^t就是要求的答案的一部分。

公式： f(pi^ai)=Φ(pi^ai)+pi*Φ(pi^(ai-1))+pi^2*Φ(pi^(ai-2))+...+pi^(ai-1)* Φ(pi)+pi^ai*Φ(1)；

这里可以把求欧拉函数的部分化简，因为f(pi^ai)中只有一个约数pi，所以：

Φ(pi^ai)=pi^ai-pi^(ai-1)。

证明：小于pi^ai的正整数个数为p^ai - 1个；

其中，和pi^ai不互质的正整数有(pi*1,pi*2,...,pi*(pi^(ai-1)-1) )共计 pi^(ai-1)-1个。

∴ Φ(pi^ai)=pi^ai -1 -(pi^(ai-1)-1)=pi^ai-pi^(ai-1)。

然后整理得：f(pi^ai)= pi^ai*(1+ai*(1-1/pi))。

所以， f(n)=n*(1+a1*(1-1/p1))*(1+a2*(1-1/p2))*.…… =n* π(ai*pi+pi-ai)/pi；

本题解决！！！

写了两个小时，写完了，发现自己好弱，大神的一句话我都得反复琢磨还要问老师各种证明才能弄懂，差距很大，好大的压力。。。

代码如下：

#include<iostream>  
#include<cmath>  
#include<cstdio>  
#include<cstring>  
#include<algorithm>  
using namespace std;  
typedef long long ll;  
ll N,ans;  
void solve()  
{  
    ans=N;
    for(ll i=2;i*i<=N;i++)  
        if(N%i==0)  
        {  
            ll j=0;  
            while(N%i==0) N/=i,j++;  
            ans/=i;  
            ans*=i*j+i-j;
			//ans*=1+j-j/i;//用上面的求法而不用这个是因为ans/i能整除而j/i不能整除，精度有问题。
        }
    if(N>1)//此时N的意义就是原来的N的一个大于sqrt(N)的质因子，再算一遍就好了。
    {  
        ans/=N;  
        ans*=N+N-1;  
    }  
    printf("%lld\n",ans);  
}  
int main()  
{  
    //freopen("a.in","r",stdin);  
    while(scanf("%lld",&N)!=EOF) solve();  
    return 0;  
}