hdu1695：数论+容斥

题目大意：

求x属于[1,b]和 y属于[1,d]的 gcd(x,y)=k 的方案数

题解：

观察发现 gcd()=k 不好处理，想到将x=x/k,y=y/k 后 gcd(x,y)=1。。

即问题转化为求区间 [1,b/k]和 [1,d/k]的互质数对个数

由于题目规定 (x,y)和(y,x)是同一种，所以我们可以规定 x<y，，然后只需对每一个y求出比他小的即可

公共部分可以通过欧拉函数快速求出。。

非公共部分就不行了。。

所以就分解质因数，用容斥的方法求了

#include <iostream>
#include <stdio.h>
#include<string.h>
#include<algorithm>
#include<string>
#include<ctype.h>
using namespace std;
#define maxn 100000
int isnotprime[maxn+1];
int num[maxn+1];
int fac[maxn+1][30];
int prime[maxn];
int euler[maxn+1];
int np;
int a,b,c,d,k,n,m;
long long ans;
void setprime()
{
    np=0;
    memset(isnotprime,0,sizeof(isnotprime));
    for(int i=2;i<=1000;i++)
    {
        if(!isnotprime[i])
           prime[np++]=i;
        for(int j=0;j<np&&prime[j]*i<=1000;j++)
        {
            isnotprime[i*prime[j]]=1;
            if(i%prime[j]==0)
                break;
        }
    }
}
void findfac(int x)
{
    num[x]=0;
    int p=x;
    for(int i=0;i<np;i++)
    {
        if(p%prime[i]==0)
        {
            fac[x][num[x]++]=prime[i];
        }
        while(p%prime[i]==0)
        {
            p/=prime[i];
        }
    }
    if(p>1)
    {
        fac[x][num[x]++]=p;
    }
}
void setfac()
{
    for(int i=2;i<=maxn;i++)
    {
        findfac(i);
    }
}
void phi()
{
    for(int i=1;i<=maxn;i++)
        euler[i]=i;
    for(int i=2;i<=maxn;i+=2)
        euler[i]/=2;
    for(int i=3;i<=maxn;i++)
    {
        if(euler[i]==i) //未被筛到。是素数，则用此素数来筛
        {
            for(int j=i;j<=maxn;j+=i)
            {
                euler[j]=euler[j]/i*(i-1);
            }
        }
    }
    return ;
}
void ini()
{
    scanf("%d%d%d%d%d",&a,&b,&c,&d,&k);
}
int iae(int x)
{
    int res=0;
    int cur,tmp;
    for(int i=1;i<(1<<num[x]);i++)
    {
        cur=1;
        tmp=0;
        for(int j=0;j<num[x];j++)
        {
            if(i&(1<<j))
            {
                cur*=fac[x][j];
                tmp++;
            }
        }
        if(tmp&1)
        {
            res+=m/cur;
        }
        else
        {
            res-=m/cur;
        }
    }
    return m-res;
}
int cas=1;
void solve()
{
    printf("Case %d: ",cas++);
    if(k==0)
    {
        puts("0");
        return;
    }
    m=min(b/k,d/k);
    n=max(b/k,d/k);
    ans=0;
    for(int i=1;i<=m;i++)
    {
        ans+=euler[i];
    }
    for(int i=m+1;i<=n;i++)
    {
        ans+=iae(i);
    }
    printf("%I64d
",ans);
}
int main()
{
    setprime();
    setfac();
    phi();
    int t;
    scanf("%d",&t);
    while(t--)
    {
        ini();
        solve();
    }
    return 0;
}