前往小程序,Get更优阅读体验!
立即前往
首页
学习
活动
专区
工具
TVP
发布
社区首页 >专栏 >专题练习---(数论)莫比乌斯反演

专题练习---(数论)莫比乌斯反演

作者头像
Gxjun
发布2018-03-26 16:37:31
8190
发布2018-03-26 16:37:31
举报
文章被收录于专栏:mlml

            GCD

Time Limit: 6000/3000 MS (Java/Others)    Memory Limit: 32768/32768 K (Java/Others) Total Submission(s): 7026    Accepted Submission(s): 2584

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

Hint

For the first sample input, all the 9 pairs of numbers are (1, 1), (1, 2), (1, 3), (1, 4), (1, 5), (2, 3), (2, 5), (3, 4), (3, 5).

Source

2008 “Sunline Cup” National Invitational Contest

题意:  给出a,b,c,d,k, 使得  a<=x<=b , c<=y<=d  &&gcd(x,y)=k的数对(x,y)的对数。

  限制:  a=c=1 ,  0<b,c<=100000 ;(n1,n2)和(n2 ,n1)算为同一种情况、

---------------------------------分割线-----------------------------------------------

个人理解:  莫比乌斯反演就是一种类似的容斥原理,首先来看看莫比乌斯函数。

   莫比乌斯函数:

                     mu(d) =  1      d==1

        mu(d) =  (-1)^k         d=p1p2p3....pk (pk为素数)

        mu(d) =  0                其他        

 莫比乌斯反演:

         求和函数:    

   其中  d|n 表示  n%d==0   1<=d<=n

         对于上面的公式,举例:

         F(1) = f(1)                         F(6) = f(1) +f(2) +f(3) + f(6)

         F(2) = f(1) +f(2)                 F(7) = f(1) + f(7)

    F(3) = f(1) + f(3)               F(8) = f(1) +f(2)+f(4)+f(8)

    F(4) = f(1) +f(2) +f(4)        F(9) = f(1)+f(3)+f(9)

    F(5) = f(1) + f(5)      F(10) = f(1)+f(2)+f(5)+f(10)

       实际中,往往是,我们求解F(N)比较容易,而求解f(x) 比较困难。但是我们亦反推求解f(x)

        f(1) = F(1)         f(6) = F(6) -F(3) -F(2) +F(1)

   f(2) = F(2) -F(1)                 f(7) = F(7) -F(1)

        f(3) = F(3) - F(1)      f(8) = F(8) - F(4)

        f(4) = F(4) -F(2)       f(9) = F(9) - F(3)

        f(5) = F(5) -F(1)       f(10) = F(10) - F(5) -F(2) +F(1)

  对于上述,我们不难总结出,已知F(x), 求解 f(x)的公式:

       d|n 意思是   n%d==0

  而且我们可以很欣喜的发现,U(d)其实就是我们上上面求解的mu(d)(即莫比乌斯函数)

  其实对于上面的公式,我们如果仔细器观察,也不难发现,这个和容斥原理还是有那么点类似的。

   对于容斥原理

         P(AUBUC) =P(A)+P(B)+P(C)-P(AUB)-P(BUC)-P(AUC)+P(ANBNC) ;

  而如  f(6) = F(6) -F(3) -F(2) +F(1)  ,这样来看,还是是挺相似的。

----------------------------------------分割线-----------------------------------------------

  那么对于这道题,我们可以这样分析:

代码语言:javascript
复制
1  设f(k)为gcd(x,y)=k的数对(x,y)的对数,我们要求的是f(1)
2 设F(k)为gcd(x,y)为k的倍数的数对(x,y)的对数,可以想到F(k)=floor(b/k)*floor(d/k),
   F(k)=E[b/k]*[d/k]*mu[k];

3 由莫比乌斯反演得:
4 令lim=min(b/k,d/k)
5 f(1)=mu[1]*F(1) + mu[2]*F[2] + ... + mu[lim]*F(lim)
6 因为(n1,n2)和(n2,n1)算为同一种情况,所以最后结果还要减掉重复的情况。

 代码:

代码语言:javascript
复制
 1 #include<cstdio>
 2 #include<cstdlib>
 3 #include<cstring>
 4 #include<iostream>
 5 #include<algorithm>
 6 #define _llint long long
 7 using namespace std;
 8 const int maxn=100010;
 9 bool jud[maxn];
10 int mu[maxn] , prime[maxn];
11 
12 void Mobius(){
13   int cnt=0;
14   memset(jud , 0 , sizeof jud);
15   mu[1]=1;   //d=1
16   for(int i=2 ; i<maxn ;i++){
17 
18       if(!jud[i]){
19           mu[i]=-mu[1];  //1*p
20           prime[cnt++]=i;
21       }
22     for(int j=0 ; j<cnt&&i*prime[j]<maxn ; j++ ){
23 
24         jud[i*prime[j]]=1;
25         if(i%prime[j]==0){
26             mu[i*prime[j]]=0 ; //²»ÊÇ»¥ÒìµÄËØÊý
27             break;
28         }else{
29              mu[i*prime[j]]=  -1*mu[i];
30         }
31     }
32   }
33   return ;
34 }
35 
36 _llint solve(int n, int m){
37 
38    _llint res=0;
39    for(int i=1; i<=n ;i++){
40      res+=(_llint)(m/i)*(n/i)*mu[i];
41    }
42    return res;
43 }
44 
45 int main()
46 {
47   int T,a,b,c,d,k;
48   Mobius();
49   scanf("%d",&T);
50   for(int i=1 ;i<=T ;i++){
51     scanf("%d%d%d%d%d",&a,&b,&c,&d,&k);
52     if(0==k)
53     {
54         printf("Case %d: 0\n",i);
55         continue;
56     }
57      b/=k  , d/=k;
58      if(b>d)    swap( b , d );
59     _llint ans1=solve(b,d);
60     _llint ans2=solve(b,b);
61      printf("Case %d: %lld\n",i,ans1-ans2/2);
62   }
63  return 0;
64 }
本文参与 腾讯云自媒体分享计划,分享自作者个人站点/博客。
原始发表:2015-06-16 ,如有侵权请联系 cloudcommunity@tencent.com 删除

本文分享自 作者个人站点/博客 前往查看

如有侵权,请联系 cloudcommunity@tencent.com 删除。

本文参与 腾讯云自媒体分享计划  ,欢迎热爱写作的你一起参与!

评论
登录后参与评论
0 条评论
热度
最新
推荐阅读
目录
  •             GCD
领券
问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档