「Luogu2257」YY的GCD

「Luogu2257」YY的GCD

蒟蒻的第一道莫反

跟着题解推的式子,但还是记录一下过程吧

本文可能在一定程度上存在谬误,请谨慎分析

若发现文中有错误,如您愿意,恳请您向我指出,不胜感激


problem

Solution

题目要求:

\[ans=\sum_{i=1}^N\sum_{j=1}^M[gcd(i,j)\in prime]\]

\(f(p)=\sum_{i=1}^N\sum_{j=1}^M[gcd(i,j)=p](p\in prime)\)

再令\(g(p)=\sum_{i=1}^N\sum_{j=1}^M[p|gcd(i,j)](p\in prime)\)

于是有

\[g(n)=\sum_{n|d}f(d)\]

反演后可得

\[f(n)=\sum_{n|d}\mu(\frac{d}{n})g(d)\]

又知\(g(d)=\lfloor\frac{N}{d}\rfloor\lfloor\frac{M}{d}\rfloor\)

于是有

\[ans=\sum_{n\in prime}f(n)=\sum_{n\in prime}\sum_{n|d}\mu (\frac{d}{n})\lfloor\frac{N}{d}\rfloor\lfloor\frac{M}{d}\rfloor\\=\sum_{n|d}\lfloor\frac{N}{d}\rfloor\lfloor\frac{M}{d}\rfloor\sum_{n\in prime}\mu(\frac{d}{n})\\=\sum_{d=1}^{min(N,M)}\lfloor\frac{N}{d}\rfloor\lfloor\frac{M}{d}\rfloor\sum_{n|d,n\in prime}\mu(\frac{d}{n})\]

\(sum(d)=\sum_{n|d,n\in prime}\mu(\frac{d}{n})\),预处理\(sum(d)\)

那么答案即

\[ans=\sum_{d=1}^{min(M,N)}\lfloor\frac{N}{d}\rfloor\lfloor\frac{M}{d}\rfloor sum(d)\]

\(\sum sum(d)\)仍可以利用前缀和优化,\(\sum\lfloor\frac{N}{d}\rfloor\lfloor\frac{M}{d}\rfloor\)利用整除分块优化,最终时间复杂度为\(O(T\sqrt{min(N,M)}+k)\)\(k\)为预处理复杂度

Code

#include <cstdio>
#include <iostream>
#include <cstdlib>
#include <cstring>
#include <cmath>
#include <algorithm>
#define maxn 10000005
#define N 10000000
using namespace std;
typedef long long ll;

template <typename T> void read(T &t)
{
    t=0;int f=0;char c=getchar();
    while(!isdigit(c)){f|=c=='-';c=getchar();}
    while(isdigit(c)){t=t*10+c-'0';c=getchar();}
    if(f)t=-t;
}

int T;
int n,m;
int pri[maxn],pcnt,nop[maxn];
int mu[maxn];
ll sum[maxn],up;

void GetPrime()
{
    nop[1]=1,mu[1]=1;
    for(register int i=2;i<=N;++i)
    {
        if(!nop[i])pri[++pcnt]=i,mu[i]=-1; 
        for(register int j=1;j<=pcnt && i*pri[j]<=N;++j)
        {
            nop[i*pri[j]]=1;
            if(i%pri[j]==0)break;
            else mu[i*pri[j]]=-mu[i];
        }
    }
    for(register int i=1;i<=pcnt;++i)
        for(register int j=1;pri[i]*j<=N;++j)
            sum[pri[i]*j]+=mu[j];
    for(register int i=1;i<=N;++i)
        sum[i]+=sum[i-1];
}

ll Calc()
{
    ll re=0;
    for(register int l=1,r;l<=up;l=r+1)
    {
        r=min(n/(n/l),m/(m/l));
        re+=1ll*(n/l)*(m/l)*(sum[r]-sum[l-1]);
    }
    return re;
}

int main()
{
    read(T);
    GetPrime();
    while(T--)
    {
        read(n),read(m);
        up=min(n,m);
        printf("%lld\n",Calc());
    }
    return 0;
}

转载于:https://www.cnblogs.com/lizbaka/p/10508947.html

  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值