[SDOI2015] 約數個數和 (莫比烏斯反演)

[SDOI2015]約數個數和

題目描述

設d(x)爲x的約數個數,給定N、M,求 \(\sum^N_{i=1}\sum^M_{j=1}d(ij)\)html

輸入輸出格式

輸入格式:c++

輸入文件包含多組測試數據。第一行,一個整數T,表示測試數據的組數。接下來的T行,每行兩個整數N、M。函數

輸出格式:測試

T行,每行一個整數,表示你所求的答案。優化

輸入輸出樣例

輸入樣例#1:spa

2
7 4
5 6code

輸出樣例#1:htm

110
121blog

說明

\(1<=N, M<=50000\)get

\(1<=T<=50000\)

Solution

莫比烏斯反演套路題

作這道題須要一個關於約數函數的公式
\[d(ij)=\sum_{x|i}\sum_{y|j}[gcd(x,y)=1]\]
關於這個公式,先感性理解一下,能夠理解成由於不能重複,因此gcd(x,y)=1,求這種數的總個數
下面是證實

設ij的質因子p
左邊\(i=ii\times p^a,j=jj\times p^b,ij=ii\times jj\times p^{a+b}\),後面的p有0~a+b共(a+b+1)種取值,即每一個因子p對d(ij)有(a+b+1)的貢獻
右邊\((xx\times p^a,yy),(xx \times p^{a-1},yy)....(xx,yy)....(xx,yy\times p^{b-1}),(xx,yy\times p^{b})\),也是(a+b+1)個數對,故上式成立
由於\(gcd(x,y)=1\),因此兩個因子之間是互不影響的

而後給我講莫比烏斯反演的shenlao告訴我只須要記住下面這個公式,作題賊爽(好像是的)
\[[gcd(i,j)=1]=\sum_{d|gcd(i,j)}\mu(d)\]
必定要記住!!!(敲重點)

如今式子變成了
\[Ans=\sum_{i=1}^n\sum_{j=1}^m\sum_{x|i}\sum_{y|j}\sum_{d|gcd(x,y)}\mu(d)\]
更換枚舉項,由枚舉gcd(i,j)的約數改成直接枚舉d
\[Ans=\sum_{i=1}^n\sum_{j=1}^m\sum_{x|i}\sum_{y|j}\sum_{d=1}^{min(n,m)}\mu(d)\times[d|gcd(x,y)]\]
\(\mu(d)\)\(i,j\)無關,提出來
\[Ans=\sum_{d=1}^{min(n,m)}\mu(d)\sum_{i=1}^n\sum_{j=1}^m\sum_{x|i}\sum_{y|j}[d|gcd(x,y)]\]
再次更換枚舉項,由分別枚舉i,j和及i,j的因子改成枚舉i,j的因子,再直接乘上每一個因子對它們產生的貢獻
\[Ans=\sum_{d=1}^{min(n,m)}\mu(d)\sum_{x=1}^n\sum_{y=1}^m[d|gcd(x,y)]\lfloor\frac{n}{x}\rfloor\lfloor\frac{m}{y}\rfloor\]
繼續更換枚舉項,將枚舉x,y換成枚舉dx,dy,這樣\([d|gcd(x,y)]\)這個條件就能夠省去,由於保證了此時x=dx,y=dy,d必定是x,y的因子,那麼枚舉範圍也要相應縮小
\[Ans=\sum_{d=1}^{min(n,m)}\mu(d)\sum_{x=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{y=1}^{\lfloor\frac{m}{d}\rfloor}\lfloor\frac{n}{x}\rfloor\lfloor\frac{m}{y}\rfloor\]
\(\lfloor\frac{n}{x}\rfloor\)對枚舉y沒有影響,提早就行了,因此最後答案就變成了
\[Ans=\sum_{d=1}^{min(n,m)}\mu(d)\sum_{x=1}^{\lfloor\frac{n}{d}\rfloor}\lfloor\frac{n}{x}\rfloor\sum_{y=1}^{\lfloor\frac{m}{d}\rfloor}\lfloor\frac{m}{y}\rfloor\]
由於多組詢問,預處理以後\(O(\sqrt n)\)處理就可了
不會預處理的先去作一作這道題[AHOI2005]約數研究
怎麼\(O(\sqrt n)\)優化?整除分塊

Code

#include<bits/stdc++.h>
#define il inline
#define rg register
#define lol long long
#define NN 50000
#define Min(a,b) (a)<(b)?(a):(b)
#define Max(a,b) (a)?(b)?(a):(b)

using namespace std;
const int N=5e4+10;

void in(int &ans) {
    ans=0; char i=getchar();
    while(i<'0' || i>'9') i=getchar();
    while(i>='0' && i<='9') ans=(ans<<1)+(ans<<3)+i-'0',i=getchar();
}

int tot,T;
int mu[N],prime[N];
lol f[N];
bool vis[N];

il void init() {
    mu[1]=1;
    for(rg int i=2;i<=NN;i++) {
        if(!vis[i]) prime[++tot]=i,mu[i]=-1;
        for(rg int j=1;j<=tot && prime[j]*i<=NN;j++) {
            vis[i*prime[j]]=1;
            if(i%prime[j]==0) break;
            else mu[i*prime[j]]=-mu[i];
        }
    }
    for(rg int i=1;i<=NN;i++) {
        lol ans=0; mu[i]+=mu[i-1];
        for(rg int l=1,r,len;l<=i;l=r+1) {
            r=i/(i/l),len=r-l+1;
            ans+=1ll*len*(i/l);
        }f[i]=ans;
    }
}

int main()
{
    in(T); init();
    while(T--) {
        int n,m,maxn; lol ans=0;
        in(n),in(m);  maxn=Min(n,m);
        for(rg int l=1,r;l<=maxn;l=r+1) {
            r=Min(n/(n/l),m/(m/l));
            ans+=1ll*(mu[r]-mu[l-1])*f[n/l]*f[m/l];
        }
        printf("%lld\n",ans);
    }
}

博主蒟蒻,隨意轉載.但必須附上原文連接

http://www.cnblogs.com/real-l/

相關文章
相關標籤/搜索