BZOJ 2301 Problem b (莫比乌斯反演+容斥)

这道题和 HDU-1695不同的是,a,c不一定是1了。还是莫比乌斯的套路,加上容斥求结果。
\(F(n,m,k)\)为满足\(gcd(i,j)=k(1\leq i\leq n,1\leq j\leq m)\)的对数。则\(ans = F(b,d,k)-F(a-1,d,k)-F(c-1,b,k)+F(a-1,c-1,k)\)
预处理莫比乌斯函数的前缀和,分块加速求和即可

#include<bits/stdc++.h>
using namespace std;
typedef long long LL;
const int maxn=1e5+5;
bool vis[maxn];
int prime[maxn],mu[maxn];
int sum[maxn];
void init(){
    memset(vis,false,sizeof(vis));
    mu[1] = 1;
    prime[0] = 0;
    int cnt=0;
    for(int i=2;i<maxn;++i){
        if(!vis[i]){
            mu[i] = -1;
            sum[i] = 1;
            prime[++cnt] = i;
        }
        for(int j=1;j<=cnt;++j){
            if(i*prime[j] >= maxn)  break;
            vis[i*prime[j]] = true;
            if(i % prime[j]){
                mu[i*prime[j]] = -mu[i];
                sum[i*prime[j]] = mu[i] - sum[i];
            }
            else{
                mu[i*prime[j]] = 0;
                sum[i*prime[j]] = mu[i];
                break;
            }
        }
    }
    for(int i =2;i<maxn;++i) sum[i]+=sum[i-1];
}
void prepare(){
    int i,j,cnt=0;
    mu[1]=sum[1]=1;
    for(i=2;i<maxn;i++){
        if(!vis[i])
            prime[++cnt]=i,mu[i]=-1;
        for(j=1;prime[j]*i<maxn;j++){
            vis[prime[j]*i]=1;
            if(i%prime[j]==0){
                mu[prime[j]*i]=0;
                break;
            }
            mu[prime[j]*i]=-mu[i];
        }
        sum[i]=sum[i-1]+mu[i];
    }
}

LL gao(LL n,LL m,LL k)
{
    if(n>m) swap(n,m);
    n/=k,m/=k;
    LL ans = 0;
    for(LL i = 1,j;i<=n;i=j+1){
        j = min(n/(n/i),m/(m/i));
        ans += (sum[j]-sum[i-1]) *(n/i) *(m/i);
    }   
    return ans; 
}

int main()
{
    #ifndef ONLINE_JUDGE
        freopen("in.txt","r",stdin);
        freopen("out.txt","w",stdout);
    #endif
    prepare();
    LL a,b,c,d,k;
    int T; scanf("%d",&T);
    while(T--){
        scanf("%lld %lld %lld %lld %lld",&a,&b,&c,&d,&k);
        LL res=0;
        res += gao(b,d,k);
        res -= gao(a-1,d,k);
        res -= gao(c-1,b,k);
        res += gao(a-1,c-1,k);
        printf("%lld\n",res);
    } 
    return 0;
}
posted @ 2018-08-26 22:58  xiuwenL  阅读(127)  评论(0编辑  收藏  举报