【BZOJ 2818】 gcd

Description

给定整数N，求1<=x,y<=N且Gcd(x,y)为素数的
数对(x,y)有多少对.

Input

一个整数N

Output

如题

Sample Input

4

Sample Output

4

HINT

hint

对于样例(2,2),(2,4),(3,3),(4,2)

1<=N<=10^7

【分析】

        设有x,y满足条件，有(x,y)=p，p为质数。那么(x/p,y/p)=1。我们定义f[k]表示(x,y)=1的方案数，1<=x<=k,1<=y<=k。

那么答案很显然为sigma(f[N/p])，p为不超过N的质数。

        即枚举(x,y)=p的每一个p，去求得方案数，然后累加。

        那么问题就在于如何求得f[k]了，很简单，f[k]=1+2*sigma(phi[j])，2<=j<=k。意为由j和与j互质的数来组成(x,y)=1，然后排除(1,1)=1被重复计算的情况。

        最后用筛法在线性时间内求得phi[]即可。

【代码】

#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cmath>
#include <iostream>
#include <algorithm>
#include <vector>

using namespace std;

#define rep(i,a,b) for(int i=(a);i<=(b);++i)
#define rrep(i,b,a) for(itn i=(b);i>=(a);--i)
#define ll long long 
#define sf scanf
#define pf printf
#define sz(x) x.size()
#define p_b push_back
#define MAXN 10000010

int N;
int phi[MAXN];
ll sum[MAXN];
vector <int> pr;

void Phi_table(int n){
	memset(phi,0,sizeof phi);
	phi[1]=1;
	rep(i,2,n){
		if(!phi[i]){
			pr.p_b(i);
			for(int j=i;j<=n;j+=i){
				if(!phi[j])
					phi[j]=j;
				phi[j]=phi[j]/i*(i-1);
			}
		}
	}
}

void Init(){
	cin>>N;
	Phi_table(N);
	rep(i,2,N){
		sum[i]=sum[i-1]+(ll)(phi[i]<<1);
	}	
}

void Solve(){
	ll ans=0;
	int n=sz(pr);
	rep(i,0,n-1){
		int k=N/pr[i];
		if(!k) break;
		ans+=1+sum[k];
	}
	cout<<ans<<endl;
}

int main(){
	Init();
	Solve();
	return 0;
}