1. 程式人生 > >容斥+莫比烏斯反演+分塊優化-BZOJ2301

容斥+莫比烏斯反演+分塊優化-BZOJ2301

對於給出的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

題目分析:首先n就是五萬,因此每次即使是O(n)的計算總的下來也是O(n^2)了,因此對於每次操作我們要讓時間複雜度小於O(n),題意很清楚
先將原式變形:
a / k <= x / k <= b / k,c / k <= y / k <= d / k,gcd(x / k, y / k) = 1,令cal(b / k,d / k)為1 <= x / k <= b / k,1 <= y / k <= d / k時取出的滿足條件的個數,則根據容斥原理有
ans = cal(b / k,d / k) - cal((a - 1) / k,d / k) - cal((c - 1) / k,b / k) + cal((a - 1) / k,(c - 1) / k),因為1~a-1和1~c-1都不在我們所求的範圍內,又減的時候這段區間減了兩次,因此要再加上一個,接下來看cal函式,這裡要用到分塊求和,如果不做優化,就是直接列舉公約數ans += mob[i] * (l / i) * (r / i)但這樣會超時,考慮到不能整除的特性,在很大一段區間內(l / i)和(r / i)的值是相同的,舉個簡單的例子,l = 10,r = 11那麼可以看出 i從6到10,l / i和r / i的值都是1,因此考慮分塊求和,從i開始最長的相等區間長度為min(l / (l / i),r / (r / i)),這裡如果不能理解的話,比如l / (l / i),設l / i = p,p表示整除時的值,那麼l / p就是從i開始整除值為p的個數了

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
#define ll long long
using namespace std;
const int maxn=50000;
int prime[maxn+10],mu[maxn+10],isprime[maxn+10],sum[maxn];
int  cnt=0;
void init()//莫比烏斯
{
    memset(isprime,0,sizeof(isprime));
    mu[1
]=1; sum[1]=1; for(int i=2; i<maxn; i++) { if(!isprime[i]) { prime[++cnt]=i,mu[i]=-1; } for(int j=1; j<=cnt&&prime[j]*i<maxn; j++) { isprime[prime[j]*i]=1; if(i%prime[j]==0) { mu[prime[j]*i
]=0; break; } else { mu[i*prime[j]]=-mu[i]; } } sum[i]=sum[i-1]+mu[i]; } } int solve(int l,int r) { if(l>r)swap(l,r); int ans=0; for(int i=1,last=0;i<=l;i=last+1)//分塊優化 { last=min(l/(l/i),r/(r/i)); ans+=(l/i)*(r/i)*(sum[last]-sum[i-1]); } return ans; } int main() { int t; int ca=1; init(); int a,b,c,d,k; scanf("%d",&t); while(t--) { scanf("%d %d %d %d %d",&a,&b,&c,&d,&k); int ans=0; ans+=solve(b/k,d/k); ans-=solve((a-1)/k,d/k); ans-=solve((c-1)/k,b/k); ans+=solve((a-1)/k,(c-1)/k); printf("%d\n",ans); } return 0; }