HDU 1695 GCD

题目链接：http://acm.hdu.edu.cn/showproblem.php?pid=1695

题目大意：给出区间（1,b）,(1,d)，正整数k,在(1,b)中取出一个数x，在(1,d)中取出一个数y,求满足gcd(x,y)==k的所有（x,y）的对数，（x,y）和(y,x)被当作相同的一对。

由欧几里得定理可以知道 gcd(x,y)==k可以写成gcd(x/k,y/k)==1，那么问题转化为在(1,b/k)和（1,d/k）中分别找出一个数x,y使得它们互质，

假设b<d，当x<y时就可以保证对数是不会重复的，在区间（1，b/k）中的对数可以通过求每个数的欧拉函数值加和得到，因此可以先预处理出1e5之内的欧拉函数值

而区间(b/k+1,d/k)中的数n与(1,b/k)中互质的数的个数求解不能用欧拉函数，那么可以用容斥定理求得，这里给出一个公式：

所有不与n互质的数的个数=1个因子倍数的个数-2个因子乘积的倍数的个数+3个因子乘积的倍数的个数...

最后所有与n互质的个数只需要用b/k减去上面式子的和。

而如果w是n的素因子，那么在区间(1,b/k)中w的倍数共有b/k/w个，对(b/k+1,d/k)中的每一个数枚举计算相加即可得到结果

 

代码：

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
int euler[100005];
int prime[100005];
bool isprime[100005];
int factor[1000][2];
void getEuler(){
    memset(euler,0,sizeof(euler));
    euler[1]=1;
    for(int i=2;i<=100000;i++){
        euler[i]=i;
    }
    for(int i=2;i<=100000;i++){
        if(euler[i]==i){
            for(int j=i;j<=100000;j+=i){
                euler[j]=euler[j]/i*(i-1);
            }
        }
    }
}
int getFat(int x){//求得数x的素因子以及它们的幂
    int k=0;
    memset(factor,0,sizeof(factor));
    int tmp=x;
    for(int i=0;prime[i]*prime[i]<=tmp;i++){
        if(tmp%prime[i]==0){
            factor[k][0]=prime[i];
            while(tmp%prime[i]==0){
                tmp/=prime[i];
                factor[k][1]++;
            }
            k++;
        }
    }
    if(tmp!=1){
        factor[k][0]=tmp;
        factor[k][1]=1;
        k++;
    }
    return k;
}
void getPrime(){
    int k=0;
    memset(isprime,true,sizeof(isprime));
    for(int i=2;i<=100000;i++){
        if(isprime[i]){
            prime[k++]=i;
        }
        for(int j=0;j<k&&prime[j]*i<=100000;j++){
            isprime[i*prime[j]]=false;
            if(i%prime[j]==0){
                break;
            }
        }
    }
}
int solve(int a,int b){
    int ans=0;
    int k=getFat(b);
    for(int i=1;i<(1<<k);i++){
        int cnt=0;
        int tmp=1;
        for(int j=0;j<k;j++){
            if(i&(1<<j)){
                cnt++;
                tmp*=factor[j][0];
            }
        }
        if(cnt&1)
            ans+=a/tmp;
        else
            ans-=a/tmp;
    }
    return a-ans;
}

int main()
{
    int n,p=1;
    scanf("%d",&n);
    getEuler();
    getPrime();
    int a,b,c,d,k;
    while(n--){
        scanf("%d%d%d%d%d",&a,&b,&c,&d,&k);
        if(k==0||k>b||k>d){
            printf("Case %d: 0\n",p++);
            continue;
        }
        b/=k;
        d/=k;
        if(b>d){
            swap(b,d);
        }
        ll ans=0;
        for(int i=1;i<=b;i++){
            ans+=euler[i];
        }
        for(int i=b+1;i<=d;i++){
            ans+=solve(b,i);
        }

        printf("Case %d: %lld\n",p++,ans);
    }
    return 0;
}