1. 程式人生 > >最容易理解的莫比烏斯反演

最容易理解的莫比烏斯反演

對於給出的n個詢問,每次求有多少個數對(x,y),滿足axbcyd,且gcd(x,y) = kgcd(x,y)函式為xy的最大公約數

1≤n≤50000,1≤a≤b≤50000,1≤c≤d≤50000,1≤k≤50000

Sample Input
2
2 5 1 5 1
1 5 1 5 2
Sample Output
14
3




上敘對於區間最遠的端點的理解:

若a為10,x為 4 ,因為 10除以4有餘數,所以會存在x在擴大的過程中儘量減少這個餘數,很容易知道當x為5的時候剛好沒有餘數,同時他們的商都是2

最後發現,其實上面敘述的P 其實就是莫比烏斯函式

而莫比烏斯函式反演有兩種形式,一般是第二種:


從這兩個形式中可以發現莫比烏斯函式的兩個應用:

對於一些函式,我們很難直接求他的值,但是很容易求這個函式的倍數或者約數的值

此時我們可以用莫比烏斯反演得到此函式,本題就是這麼個道理

#include<stdio.h>
#include<string.h>
#include<algorithm>
using namespace std;

#define MAXN 100000
#define LL long long
bool check[MAXN+10];
int primer[MAXN+10];
int mu[MAXN+10];
int sum[MAXN+10];

void Moblus()
{
    memset(check,0,sizeof(check));
    mu[1]=1;
    int tot=0;
    for(int i=2;i<=MAXN;i++){
        if(!check[i]){
            primer[tot++]=i;
            mu[i]=-1;
        }
        for(int j=0;j<tot&&i*primer[j]<=MAXN;j++){
            check[i*primer[j]]=true;
            if(i%primer[j]==0){
                mu[i*primer[j]]=0;
                break;
            }
            mu[i*primer[j]]=-mu[i];
        }
    }
    sum[0]=0;
    for(int i=1;i<=MAXN;i++)
        sum[i]=sum[i-1]+mu[i];
}

LL solve(int n,int m)
{
    LL ans=0;
    if(n>m) swap(n,m);
    for(int i=1,last;i<=n;i=last+1){
        last=min(n/(n/i),m/(m/i));
        ans+=(LL)(sum[last]-sum[i-1])*(n/i)*(m/i);
    }
    return ans;
}

int main()
{
    int T;
    Moblus();
    int a,b,c,d,k;
    //freopen("in.txt","r",stdin);
    scanf("%d",&T);
    while(T--)
    {
        scanf("%d%d%d%d%d",&a,&b,&c,&d,&k);

        LL ans=solve(b/k,d/k)-solve((a-1)/k,d/k)
        -solve(b/k,(c-1)/k)+solve((a-1)/k,(c-1)/k);
        printf("%lld\n",ans);
    }
    return 0;
}

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).
        
 

本題和上題不同的是,需要去重

(1,2)和(2,1)是一樣的

#include<stdio.h>
#include<string.h>
#include<algorithm>
using namespace std;

#define MAXN 1000000
#define LL long long
bool check[MAXN+10];
LL primer[MAXN+10];
LL mu[MAXN+10];
LL sum[MAXN+10];

void Moblus()
{
    memset(check,0,sizeof(check));
    mu[1]=1;
    int tot=0;
    for(int i=2;i<=MAXN;i++){
        if(!check[i]){
            primer[tot++]=i;
            mu[i]=-1;
        }
        for(int j=0;j<tot&&i*primer[j]<=MAXN;j++){
            check[i*primer[j]]=true;
            if(i%primer[j]==0){
                mu[i*primer[j]]=0;
                break;
            }
            mu[i*primer[j]]=-mu[i];
        }
    }
    sum[0]=0;
    for(int i=1;i<=MAXN;i++)
        sum[i]=sum[i-1]+mu[i];
}

LL solve(int n,int m)
{
    LL ans=0;
    if(n>m) swap(n,m);
    for(int i=1,last;i<=n;i=last+1){
        last=min(n/(n/i),m/(m/i));
        ans+=(LL)(sum[last]-sum[i-1])*((LL)(n/i)*(m/i));
    }
    return ans;
}

int main()
{
    int T;
    Moblus();
    int a,b,c,d,k,cases=1;
   // freopen("in.txt","r",stdin);
    scanf("%d",&T);
    while(T--)
    {
        printf("Case %d: ",cases++);
        scanf("%d%d%d%d%d",&a,&b,&c,&d,&k);
        if(k==0){
            puts("0");
            continue;
        }
        LL ans=solve(b/k,d/k)-solve((a-1)/k,d/k)
        -solve(b/k,(c-1)/k)+solve((a-1)/k,(c-1)/k);
        int x=max(a,c);
        int y=min(b,d);
        if(y>=x){
            ans-=solve(y/k,y/k)/2;
            ans+=solve(x/k,x/k)/2;
        }
        printf("%I64d\n",ans);
    }
    return 0;
}