hdu 1695(欧拉公式，容斥原理||莫比乌斯反演)

古城微笑少年丶 2022-06-10 07:42 168阅读 0赞

#                                                        GCD #

**                                                             Time Limit: 6000/3000 MS (Java/Others)    Memory Limit: 32768/32768 K (Java/Others)  
                                                                                  Total Submission(s): 12186    Accepted Submission(s): 4615**

Problem Description

Given 5 integers: a, b, c, d, k, you're to find x in a...b, y in c...d that GCD(x, y) = k. GCD(x, y) means the greatest common divisor of x and y. Since the number of choices may be very large, you're only required to output the total number of different number pairs.  
Please notice that, (x=5, y=7) and (x=7, y=5) are considered to be the same.  
  
**Yoiu can assume that a = c = 1 in all test cases.**

Input

The input consists of several test cases. The first line of the input is the number of the cases. There are no more than 3,000 cases.  
Each case contains five integers: a, b, c, d, k, 0 < a <= b <= 100,000, 0 < c <= d <= 100,000, 0 <= k <= 100,000, as described above.

Output

For each test case, print the number of choices. Use the format in the example.

Sample Input

2 1 3 1 5 1 1 11014 1 14409 9

Sample Output

Case 1: 9 Case 2: 736427

方法一：

题目分析:由于gcd(x,y)=k,我们知道gcd(x/k,y/k)=1,所以我们可以将该题目化为在1到x/k之内与1到y/k之内的互质数的对数。

我们知道欧拉公式可以解决小于等于n之内的互质数的对数，那么大于n的了?对于大于n的数值我们可以从反面思考，假设求\[1,m\]内与n互质的对数，则我们可以求\[1,m\]与n不互质的对数x,则答案是m-x.对于不互质的对数，我们可以利用容斥原理来求。

不互质，这必定含有公共质因子，假设m的质因子是p1,p2,p3,...pn.则1到m之间含有p1质因子的数目为m/p1，同理可得m/pn.但是这里面有重复的，比如6被2除了，也被3除了，所以我们就要利用容斥原理了。

![Center][]

奇数减，偶数加

#include<iostream>
    #include<cmath>
    #include<cstdio>
    #include<algorithm>
    #include<functional>
    #include<cstring>
    #define ll long long
    using namespace std;
    
    const int maxn=1e5+10;
    bool vis[maxn];
    int prime[maxn],cnt;
    void get_prime()//质数打表
    {
        memset (vis,true,sizeof (vis));
        vis[1]=false;
        for (int i=2;i<maxn;i++) {
            if (vis[i]) prime[cnt++]=i;
            for (int j=0;j<cnt&&i*prime[j]<maxn;j++) {
                vis[i*prime[j]]=false;
                if (i%prime[j]==0) break;
            }
        }
    }
    
    int euler[maxn];
    void get_euler()//利用欧拉公式打表，求1到n与n的互质对数
    {
        euler[1]=1;
        for (int i=2;i<maxn;i++)
            euler[i]=i;
        for (int i=2;i<maxn;i++) {
            if (euler[i]==i)
            for (int j=i;j<maxn;j+=i) {
                euler[j]=euler[j]/i*(i-1);
            }
        }
    }
    
    int factor[1000][2];
    int facnt;
    void get_factor(int n)//分解质因子
    {
        facnt=0;
        memset (factor,0,sizeof (factor));
        for (int i=0;i<cnt&&prime[i]*prime[i]<=n;i++) {
            if (n%prime[i]==0) {
                factor[facnt][0]=prime[i];
                while (n%prime[i]==0) {
                    factor[facnt][1]++;
                    n=n/prime[i];
                }
                facnt++;
            }
        }
        if (n!=1) {
            factor[facnt][0]=n;
            factor[facnt++][1]=1;
        }
    }
    ll get_ans(int n,int m)//容斥原理，这段代码也可以用dfs来写，这里用了状态压缩
    {
        get_factor(m);
        ll ans=0;
        for (int i=1;i<(1<<facnt);i++) {//i遍历了所有的质因子的可能，i看成二进制，每一位就表示一个质因子
            int Cnt=0,temp=1;
            for (int j=0;j<facnt;j++) {//去检查i表示那几个质因子
                if (i&(1<<j)) {
                    Cnt++;//表示有几个质因子
                    temp*=factor[j][0];//将他们乘起来
                }
            }
            if (1&(Cnt))//与1与运算，为1就是奇数
                ans+=n/temp;
            else
                ans-=n/temp;
        }
        return n-ans;
    }
    int main()
    {
        get_prime();
        get_euler();
        int t;
        scanf("%d",&t);
        for (int iCase=1;iCase<=t;iCase++) {
            int a,b,c,d,k;
            scanf("%d%d%d%d%d",&a,&b,&c,&d,&k);
            if (k==0) {
                printf("Case %d: %d\n",iCase,0);
                continue;
            }
            if (b>d) swap (b,d);//为了计算方便，我们把小的数当做范围
            b=b/k;
            d=d/k;
            ll ans=0;
            if (b==0||d==0) {
                printf("Case %d: %d\n",iCase,0);
                continue;
            }
            for (int i=1;i<=b;i++)//求出小于等于n的互质对数
                ans+=euler[i];
            for (int i=b+1;i<=d;i++)//求出大于n的数与（1到n)的互质对数
                ans+=get_ans(b,i);
            printf("Case %d: %lld\n",iCase,ans);
        }
        return 0;
    }

方法二: 莫比乌斯反演

由于gcd(x,y)=k,我们知道gcd(x/k,y/k)=1,所以我们可以将该题目化为在1到x/k之内与1到y/k之内的互质数的对数。

![SouthEast][]

其中根据了莫比乌斯函数和的性质

![SouthEast 1][]

记得去重

代码如下:

#include<iostream>
    #include<cmath>
    #include<cstring>
    #include<cstdio>
    #define ll long long
    using namespace std;
    
    const int MAXN=1e6+100;
    bool check[MAXN+10];
    int prime[MAXN+10];
    int mu[MAXN+10];
    void Moblus()
    {
        memset(check,false,sizeof(check));
        mu[1] = 1;
        int tot = 0;
        for(int i = 2; i <= MAXN; i++)
        {
            if( !check[i] ){
                prime[tot++] = i;
                mu[i] = -1;
            }
            for(int j = 0; j < tot; j++)
            {
                if(i * prime[j] > MAXN) break;
                check[i * prime[j]] = true;
                if( i % prime[j] == 0){
                    mu[i * prime[j]] = 0;
                    break;
                }else{
                    mu[i * prime[j]] = -mu[i];
                }
            }
        }
    }
    
    int main()
    {
        Moblus();
        int t;
        scanf("%d",&t);
        for (int icase=1;icase<=t;icase++) {
            int a,b,c,d,k;
            scanf("%d%d%d%d%d",&a,&b,&c,&d,&k);
            if(k == 0){
                printf("Case %d: 0\n",icase);
                continue;
            }
            b /= k;
            d /= k;
            if(b > d)swap(b,d);
            ll ans1 = 0;
            for(int i = 1; i <= b;i++)
                ans1 += (ll)mu[i]*(b/i)*(d/i);
            ll ans2 = 0;
            for(int i = 1;i <= b;i++)
                ans2 += (ll)mu[i]*(b/i)*(b/i);
            ans1 -= ans2/2;
            printf("Case %d: %I64d\n",icase,ans1);
        }
        return 0;
    }

[Center]: /images/20220610/bddf252f56de448f84cf827067ddc27e.png
[SouthEast]: /images/20220610/b6bb36c6d4c645f4899f789a24f99095.png
[SouthEast 1]: /images/20220610/25c25a6057764287928a695c34c8fec9.png