BZOJ-2301 [HAOI2011]Problem b(莫比乌斯反演+容斥)

题目描述

  计算：

\[\sum\limits_{x=a}^{b}\sum\limits_{y=c}^{d}[\gcd(x,y)=k] \]

  数据范围：\(1\leq T,k,a,b,c,d\leq 5\times 10^4\)。

分析

  设 \(f(a,b)=\displaystyle\sum_{i=1}^{a}\sum_{j=1}^{b}[\gcd(i,j)=k]\)。

  则：

\[\begin{aligned}&f(a,b)\\ =&\sum_{i=1}^{a}\sum_{j=1}^{b}[\gcd(i,j)=k]\\ =&\sum_{i=1}^{\frac{a}{k}}\sum_{j=1}^{\frac{b}{k}}[\gcd(i,j)=1]\\ =&\sum_{i=1}^{\frac{a}{k}}\sum_{j=1}^{\frac{b}{k}}\sum_{d\mid \gcd(i,j)}\mu(d)\\ =&\sum_{i=1}^{\frac{a}{k}}\sum_{j=1}^{\frac{b}{k}}\sum_{d=1}^{\frac{a}{k}}\mu(d)·[d\mid \gcd(i,j)]\\ =&\sum_{d=1}^{\frac{a}{k}}\mu(d)\sum_{i=1}^{\frac{a}{k}}\sum_{j=1}^{\frac{b}{k}}[d\mid \gcd(i,j)]\\ =&\sum_{d=1}^{\frac{a}{k}}\mu(d)\Big\lfloor\frac{a}{kd}\Big\rfloor\Big\lfloor\frac{b}{kd}\Big\rfloor \end{aligned} \]

  答案为：

\[\sum\limits_{x=a}^{b}\sum\limits_{y=c}^{d}[\gcd(x,y)=k]\\ =\sum\limits_{x=1}^{b}\sum\limits_{y=1}^{d}[\gcd(x,y)=k]-\sum\limits_{x=1}^{b}\sum\limits_{y=1}^{c-1}[\gcd(x,y)=k]-\sum\limits_{x=1}^{a-1}\sum\limits_{y=1}^{d}[\gcd(x,y)=k]+\sum\limits_{x=1}^{a-1}\sum\limits_{y=1}^{c-1}[\gcd(x,y)=k]\\ =f(b,d)-f(b,c-1)-f(a-1,d)+f(a-1,c-1) \]

代码

#include<bits/stdc++.h>
using namespace std;
const int N=1e5+10;
bool vis[N+10];
long long cnt,prime[N+10],mu[N+10],sum[N+10];
void init()
{
    mu[1]=1;
    for(int i=2;i<=N;i++)
    {
        if(!vis[i])
        {
            prime[++cnt]=i;
            mu[i]=-1;
        }
        for(int j=1;j<=cnt&&i*prime[j]<=N;j++)
        {
            vis[i*prime[j]]=1;
            if(i%prime[j]==0)
                break;
            else
                mu[i*prime[j]]=-mu[i];
        }
    }
    for(int i=1;i<=N;i++)
        sum[i]=sum[i-1]+mu[i];
}
long long solve(long long b,long long d,long long k)
{
    b=b/k,d=d/k;
    long long ans=0,n=min(b,d);
    for(int l=1,r;l<=n;l=r+1)
    {
        r=min(b/(b/l),d/(d/l));
        ans=ans+(sum[r]-sum[l-1])*(b/l)*(d/l);
    }
    return ans;
}
int main()
{
    init();
    int T;
    cin>>T;
    while(T--)
    {
        long long a,b,c,d,k;
        scanf("%lld %lld %lld %lld %lld",&a,&b,&c,&d,&k);
        printf("%lld\n",solve(b,d,k)-solve(b,c-1,k)-solve(a-1,d,k)+solve(a-1,c-1,k));
    }
    return 0;
}