BZOJ 2301-Problem b（莫比乌斯反演）

2301: [HAOI2011]Problem b

Time Limit: 50 Sec  Memory Limit: 256 MB
Submit: 5747  Solved: 2621
[ Submit][ Status][ Discuss]

Description

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

Input

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

 

Output

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

 

Sample Input

2

2 5 1 5 1

1 5 1 5 2

Sample Output

14

3

HINT

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

Source

        学了将近一天的莫比乌斯，感觉还是没吃透，好难呀，主要还是不明白原理和证明，做这道题也是超级牵强，还是有很多不懂的地方，翻了官方的解题报告也是一脸懵逼，然后就是翻来覆去的看，思路倒是理解的差不多了，但是证明还是难懂，觉得还是多写几道较好。。。本题的思路好像和一道GCD很像，然而直接杠了这道题，遇到了很多瓶颈，首先便是问题的转化，听大佬们说，对于有下界的区域来说，可以考虑用容斥原理来做。比如本题的x,y的取值范围来说，我们可以将其转化为1<=x<=n和1<=y<=m;减去重复的话就是容斥了。。

       我们先定义一个函数g(a,b,k)表示1<=x<=a && 1<=y<=b的所有满足gcd(x,y)=k的数的对数，呢么由容斥原理得：

                                            ans=g(c,d,k)−g(a−1,d,k)−g(b,c−1,k)+g(a−1,c−1,k)

1<=x<=n,1<=y<=m，满足gcd(x,y)是k的(x,y)的对数也等价于1<=x<=n/k,1<=y<=m/k,(x,y)互质的对数，即 

g(n,m,k)=g(n/k,m/k,1)

                                                      g ( n , m , k ) = g ( n / k , m / k , 1 )

设f(i)表示满足gcd(x,y)=i的（x,y）de 对数，F(i)表示满足i|gcd(x,y)的(x,y)的对数，则F(i)=[n/i][m/i]:[]表示向下取整

则根据莫比乌斯反演定理我们可以得出

                

然后将i=1带入即可。。。最后的难点在于此时实际上复杂度是O(n)的，然而我们还有100000次查询，故还是会超时，因此我们要找一个小于O（n）的解决方法，回到式子中，我们发现[n/d]这部分，你会发现其实他只有不超过2sqrt（n）个不同的取值，呢就分块咯，剩下的就是对μ函数取个前缀和，毕竟你分块了。。。。。。

还是推荐这个PPT：https://wenku.baidu.com/view/fbec9c63ba1aa8114431d9ac.html

#include<map>   
#include<stack>          
#include<queue>          
#include<vector>  
#include<string>
#include<math.h>          
#include<stdio.h>          
#include<iostream>          
#include<string.h>          
#include<stdlib.h>  
#include<algorithm> 
#include<functional>  
using namespace std;          
typedef long long  ll;         
#define inf  1000000000     
#define mod 1000000007           
#define maxn  160050
#define lowbit(x) (x&-x)          
#define eps 1e-9
ll a[maxn]={1,1},b[maxn],mu[maxn],sum[maxn],cnt;
void init()
{
	ll i,j;mu[1]=1;
	for(i=2;i<maxn;i++)
	{
		if(a[i]==0)
			b[++cnt]=i,mu[i]=-1;
		for(j=1;j<=cnt && i*b[j]<=maxn;j++)
		{
			a[b[j]*i]=1;
			if(i%b[j]==0)
			{
				mu[b[j]*i]=0;
				break;
			}
			else
				mu[b[j]*i]=-mu[i];
		}
	}
}
ll solve(ll n,ll m)
{
	ll i,j,last,res=0;
	if(n>m) swap(n,m);
	for(i=1,last=0;i<=n;i=last+1)
	{
		last=min(n/(n/i),m/(m/i));
		res+=(sum[last]-sum[i-1])*(n/i)*(m/i);
	}
	return res;
}
int  main(void)
{
	init();
	ll T,i,j,x,y,xx,yy,k;
	for(i=1;i<=50005;i++)
		sum[i]=sum[i-1]+mu[i];
	scanf("%lld",&T);
	while(T--)
	{
		scanf("%lld%lld%lld%lld%lld",&x,&y,&xx,&yy,&k);
		ll ans=solve(y/k,yy/k)-solve((x-1)/k,yy/k)-solve((xx-1)/k,y/k)+solve((x-1)/k,(xx-1)/k);
		printf("%lld\n",ans);
	}
	return 0;
}