洛谷 P2522 [HAOI2011]Problem b 莫比乌斯反演

题目描述

对于给出的n个询问，每次求有多少个数对(x,y)，满足a≤x≤b，c≤y≤d，且gcd(x,y) = k，gcd(x,y)函数为x和y的最大公约数。

输入输出格式

输入格式：
第一行一个整数n，接下来n行每行五个整数，分别表示a、b、c、d、k

输出格式：
共n行，每行一个整数表示满足要求的数对(x,y)的个数

输入输出样例

输入样例#1：
2
2 5 1 5 1
1 5 1 5 2
输出样例#1：
14
3
说明

100%的数据满足：1≤n≤50000，1≤a≤b≤50000，1≤c≤d≤50000，1≤k≤50000

分析：刚刚学了莫比乌斯反演，打道版题。
我们设ans(n,m)表示满足x<=n,y<=m且gcd(x,y)=k的数对(x,y)的数量，那么答案就是ans(b,d)-ans(a-1,d)-ans(b,c-1)+ans(a-1,c-1)
问题就转换成了如何求ans(n,m)
我们设F(i)表示x<=n,y<=m且i|gcd(x,y)的数对(x,y)的数量，f(i)表示x<=n,y<=m且gcd(x,y)=i的数对数量。然后就比较显然了，其实纠结最久的是最后用sqrt(n)求解过程。
先跑一个前缀和，其实想快速统计miu(i)×(n/i/k)×(m/i/k)，就是要后面的(n/i/k)*(m/i/k)相等，就可以分配律搞一搞，显然要解一个方程n/k/i=n/k/j，其中j最大，i最小，最小可以用上一个j+1得到，那么j怎么求呢？
设p=n/k,所以有c=p/i。这个好像是一个定理，不管它。有p/j=c,则p>=j*c,所以j<=p/c，所以j的最大值就是p/c,把p=n/k带入，即n/k/(n/k/i),m同理，取最小值，就可以了。

代码：

#include <iostream>
#include <cstdio>
#include <cmath>

const int maxn=50007;

using namespace std;

int mu[maxn],prime[maxn],not_prime[maxn],sum[maxn];
int n,i,A,B,C,D,k,cnt;

void get_mu(int n)//筛莫比乌斯函数
{
    int i,j;
    mu[1]=1;
    for (i=2;i<=n;i++)
    {
        if (!not_prime[i])
        {
            prime[++cnt]=i;
            mu[i]=-1;
        }
        for (j=1;j<=cnt;j++)
        {
            if (prime[j]*i>n) break;
            not_prime[i*prime[j]]=1;
            if (i%prime[j]==0)
            {
                mu[i*prime[j]]=0;
                break;
            }
            mu[i*prime[j]]=-mu[i];
        }
    }
    sum[1]=mu[1];
    for (i=2;i<=n;i++) sum[i]=sum[i-1]+mu[i];//前缀和
}

int get(int n,int m,int k)
{
    if (n>m) swap(n,m);//保证n<=m
    int i=1,j=0,ans=0;//i为左边界，j为右边界
    for (i=1;i<=n/k;i=j+1)//显然i最大为n/k,因为n/k/i>=1
    {
        j=min(n/k/(n/k/i),m/k/(m/k/i));//这个很明显啦，解释过了
        ans+=(n/i/k)*(m/i/k)*(sum[j]-sum[i-1]);//求和
    }
    return ans;
}

int main()
{
    scanf("%d",&n);
    get_mu(50000);
    for (i=1;i<=n;i++)
    {
        scanf("%d%d%d%d%d",&A,&B,&C,&D,&k);     
        int ans=get(B,D,k)-get(A-1,D,k)-get(B,C-1,k)+get(A-1,C-1,k);//容斥一下
        printf("%d\n",ans);
    }
}