HDU1695 GCD 数论之 莫比乌斯反演

做了一段时间的DP，继续回头啃数论了，这是一道莫比乌斯反演的题目，比较繁琐

给你a,b,c,d,k五个数，其中a=c=1固定的，让你从[a,b]中找出x，[c,d]中找出y,是的gcd(x,y) == k，注意gcd(x,y) 与gcd(y,x)归为同一种，问一共能找到多少组x,y;

分析：

因为gcd(x,y) = k，充要条件gcd(x/k,y/k) = 1,所以区间可以缩小为[1/k,b/k]，[1/k,d/k]，注意此处的d是题目中的d，不是数论中默认的d最大为公约数含义,这样酒吧问题转换为 去两个区间中的元素x,y使得gcd(x,y) = 1，因为gcd(y,x)和gcd(x,y)是相同的。假设x<y，可以先取小区间进行操作，由于gcd(x,y) = 1，那么可以利用欧拉函数
区间计算为 ans+=phi[1]+phi[2]+……+phi[b]

对于[b/k+1,d/k]，设y为区间的一个元素，对y进行素因子分解，得到集合｛p1,p2,p3……｝pi为素数，求的是gcd(x,y)=1的组合，但是反过来可以求cgd(x,y)!=1的组合，这样就是利用了容斥原理进行统计能被这些素数整出的数的个数，最后相减，求补数即可，然后再加上欧拉函数得到的ans值 就是最后的答案了，

#include<iostream>
#include<cstdio>
#include<list>
#include<algorithm>
#include<cstring>
#include<string>
#include<queue>
#include<stack>
#include<map>
#include<vector>
#include<cmath>
#include<memory.h>
#include<set>

#define ll long long

#define eps 1e-8

#define inf 0xfffffff
//const ll INF = 1ll<<61;

using namespace std;

//vector<pair<int,int> > G;
//typedef pair<int,int > P;
//vector<pair<int,int> > ::iterator iter;
//
//map<ll,int >mp;
//map<ll,int >::iterator p;

//#define IN freopen("c:\\Users\\linzuojun\\desktop\\input.txt", "r", stdin)
//#define OUT freopen("c:\\Users\\linzuojun\\desktop\\output.txt", "w", stdout)

//筛选法 同时求出了n以内的素数
const int N = 100000 + 5;
int num
,cnt;
int a,b,c,d,k,sum;

int phi
,prime
;
bool vis
;

int len;

void init()
{
ll i,j;
memset(vis,false,sizeof(vis));
vis[0]=vis[1]=true;
len = 0;
for(i=2;i<=N;i++){
if(!vis[i])
prime[len++] = i;
for(int j=0;j<len && i*prime[j] < N;j++) {
vis[i*prime[j]] = true;
if(!(i%prime[j]))break;
}
} //这段求出了N内的所有素数
for(i=1;i<=N;i++) {
phi[i]=i;
}
//sum[0]=0;
//sum[1]=phi[1];
memset(phi, 0, sizeof(phi));
phi[1] = 1;
for(int i = 2; i < N; i++) if(!phi[i]){
for(int j = i; j < N; j+=i){
if(!phi[j]) phi[j] = j;
phi[j] = phi[j] / i * (i-1);
}
}
}//这段求出了N内所有数的欧拉函数值

void dfs(int i,int u,int x,int v) {
if(u == x) {
sum += b/v;
return;
}
if(i == cnt)return;
dfs(i+1,u+1,x,v*num[i]);
dfs(i+1,u,x,v);
return;
}

int cal() {//容斥原理
int s = 0;
for(int i=1;i<=cnt;i++) {
sum = 0;
dfs(0,0,i,1);
if(i&1)
s += sum;//奇数加
else
s -= sum;//偶数减
}
return b - s;
}

int main() {
init();
int t;
int Case = 0;
scanf("%d",&t);
while(t--) {
scanf("%d %d %d %d %d",&a,&b,&c,&d,&k);
if(k == 0) {
printf("Case %d: %d\n",++Case,0);
continue;
}
b /= k;
d /= k;
if(b > d) {
a = b;
b = d;
d = a;
}
ll ans = 0;
for(int i=1;i<=b;i++)
ans += phi[i];
for(int i=b+1;i<=d;i++) {
cnt = 0;
a = i;
for(int j=0;j<len && prime[j]*prime[j] <= a;j++)
if(!(a%prime[j])) {
num[cnt++] = prime[j];
do
a /= prime[j];
while(a%prime[j] == 0);
}
if(a > 1)
num[cnt++] = a;
ans += cal();
}
printf("Case %d: %I64d\n",++Case,ans);
}
return EXIT_SUCCESS;
}