BZOJ2301 [HAOI2011]Problem b

本文版權歸ljh2000和博客園共有，歡迎轉載，但須保留此聲明，並給出原文鏈接，謝謝合作。

  

本文作者：ljh2000
作者博客：http://www.cnblogs.com/ljh2000-jump/
轉載請注明出處，侵權必究，保留最終解釋權！

 

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

 

 

正解：線性篩+莫比烏斯函數+分塊

解題報告：

　　 這道題是一道莫比烏斯函數的入門題，我從昨晚一直看題解和分析，直到今天才完全弄懂。系統性地寫一份解題報告吧。

　　 ps：我在學習的時候參考了這些論文與博客，貼出來以供參考吧：

　　 http://wenku.baidu.com/view/fbe263d384254b35eefd34eb.html

　　 http://wenku.baidu.com/view/fbec9c63ba1aa8114431d9ac.html

　　 http://www.cnblogs.com/zhsl/p/3269288.html

 

　　 可以發現題目要求的是∑ ∑ [gcd(i,j)==k] (1<=i<=n && 1<=j<=m)，這樣一個東西。（至於區間的話我們很快可以發現根據容斥原理[a,b]和[c,d]組合起來的答案就等於[1,b]組合[1,d] - [1,a-1]組合[1][d] - [1][b]組合[1][c-1] + [1][a-1] 組合 [1][c-1]，問題從而轉化。）但是直接求的復雜度很高，還是多組詢問，肯定會TLE。

　　題目所求的可以變成

　　∑ ∑ [gcd(i,j)==1] (1<=i<=n/k && 1<=j<=m/k) ，相當於是同時除以k；

　　考慮再把式子變形：

　　∑ ∑ ∑ μ(d) (1<=i<=n && 1<=j<=m && d|gcd(i,j))；

　　

      根據莫比烏斯函數的性質，當n=1時，∑μ(i)=1(1<=i<=n)；否則∑μ(i)=0(1<=i<=n)，所以實際上只有在gcd(i,j)==1時上式才有貢獻。可是復雜度還是沒有改變啊。我看了一篇解題報告，里面對於這個問題講的很詳細，大概是說根據容斥原理，可以把1<=i<=n && 1<=j<=m且滿足gcd(i,j)==k的方案數轉化成1<=i<=n && 1<=j<=m且滿足gcd(i,j)>=k的方案數。二者的聯系就在於莫比烏斯函數。論文里面畫了一個圖，一目了然地展示了二者的聯系。

 

　　可以根據莫比烏斯函數的性質很快反應過來二者的聯系。

 　  所以得到新式：

    　　 ans=∑μ(i)*[n/i]*[m/i] (1<=i<=min(n,m),[n/i]表示n/i下取整)

　　容易發現這個式子可以使我的復雜度降低到O(n)，但是詢問次數也是5w級別，還是會TLE。考慮繼續優化上式，考慮到[n/i]、[m/i]都會有大量的完全相等的部分，我們可以把[n/i]、[m/i]都相等的部分放在一起算，也就是一個分塊的思想。容易發值為[n/i]的起點是i，終點是n/([n/i])，而取值不會超過根號n種，也就是說我只要算根號次。但是因為是兩個都要相等，在數軸上畫線段可以發現，這樣的話總取值不會超過2×根號n種（論文里面也有圖說明）。所以復雜度降低到可以接受的范圍內。

　　值得一提的是，由於我們是分塊計算的，所以我們必須要預處理出莫比烏斯函數和他的前綴和。預處理的方式是根據莫比烏斯函數的性質，進行一次線性篩。在篩出質數的同時隨便求一下莫比烏斯函數，比較巧妙，然后直接構前綴和就可以了。

 

//It is made by ljh2000
#include <iostream>
#include <cstdlib>
#include <cstring>
#include <cstdio>
#include <cmath>
#include <algorithm>
#include <ctime>
#include <vector>
#include <queue>
#include <map>
#include <set>
using namespace std;
typedef long long LL;
const int inf = (1<<30);
const int MAXN = 50011;
const int MAXL = 50000;
LL a,b,c,d,k,ans;
int prime[MAXN],cnt;
bool ok[MAXN];
int mobius[MAXN],sum[MAXN];//莫比烏斯函數及其前綴和

inline int getint()
{
    int w=0,q=0; char c=getchar();
    while((c<'0' || c>'9') && c!='-') c=getchar(); if(c=='-') q=1,c=getchar(); 
    while (c>='0' && c<='9') w=w*10+c-'0', c=getchar(); return q ? -w : w;
}

inline void init(){//遞推得到莫比烏斯函數
    //1的莫比烏斯函數是1；質數的莫比烏斯函數為1；含有相同質因子的數莫比烏斯函數為0；
    //不含有相同質因子的數的莫比烏斯函數函數為（-1）^k，k為質因子個數
    mobius[1]=1;
    for(int i=2;i<=MAXL;i++) {
    if(!ok[i]) {
        mobius[i]=-1;
        prime[++cnt]=i;
    }
    for(int j=1;j<=cnt && i*prime[j]<=MAXL;j++) {
        ok[i*prime[j]]=1;//標記合數
        if(i%prime[j]) mobius[i*prime[j]]=-mobius[i];//互質的兩個數乘起來得到一個不含有相同質因子的數，質因子個數奇偶性改變，莫比烏斯函數變號
        else { mobius[i*prime[j]]=0; break; }//留到后面再篩，此處已經可以break
    }
    }
    for(int i=1;i<=MAXL;i++) sum[i]=sum[i-1]+mobius[i];//求莫比烏斯函數的前綴和
}

inline LL solve(LL n,LL m){//計算a在[1,n]且b在[1,m]中的gcd(a,b)==1的數目
    n/=k; m/=k;
    if(n>m) swap(n,m); if(n==0) return 0;
    LL i,next1,next2,next;//把相等的部分直接分塊一起計算
    LL tot=0;
    for(i=1;i<=n;i=next) {
    next1=n/(n/i); next2=m/(m/i);
    next=min(next1,next2);
    tot+=(n/i)*(m/i)*(sum[next]-sum[i-1]);
    next++;
    }
    return tot;
}

inline void work(){
    int T=getint(); init();
    while(T--) {
    a=getint(); b=getint(); c=getint(); d=getint(); k=getint();
    ans=solve(b,d)-solve(a-1,d)-solve(b,c-1)+solve(a-1,c-1);//容斥原理
    printf("%lld\n",ans);
    }
}

int main()
{
    work();
    return 0;
}