YY的gcd

YY的GCD

原题链接

题目

神犇YY虐完数论后给傻×kAc出了一题

给定N, M,求1<=x<=N, 1<=y<=M且gcd(x, y)为质数的(x, y)有多少对

kAc这种傻×必然不会了,于是向你来请教……

多组输入

输入输出格式

输入格式:
  • 第一行一个整数T 表述数据组数
  • 接下来T行,每行两个正整数,表示N, M
输出格式:
  • T行,每行一个整数表示第i组数据的结果
说明
  • T=10000,N, M <= 10000000

解题思路

仔细读完题之后可以发现本题要求的答案ANS=\(\sum_{i=1}^{n}\sum_{j=1}^{m}[gcd(x,y)=prim]\)

那么我们开始设出两个非常套路的函数

\(f(d)=\sum_{i=1}^n \sum_{j=1}^m[gcd(i,j)=d]\),\(F(n)=\sum_{n|d}f(d)=\lfloor\frac{N}{n}\rfloor\lfloor\frac{M}{n}\rfloor\)

\(f(d)\)表示\(gcd(i,j)=d\)的方案数,\(F(n)\)表示\(gcd(i,j)=n\)\(k*n,k∈N^*\)的方案数

然后我们就可以开始用莫比乌斯函数愉快的推式子了

\(f(n)=\sum_{n|d}\mu(\lfloor\frac{d}{n}\rfloor)F(d)\)

之后开始化简ANS

\(Ans=\sum_{p\in prime}\sum_{i=1}^{n}\sum_{j=1}^{m}[gcd(i,j)=p]=\sum_{p\in prime}f(p)\)

将刚刚反演得到的式子带入一下有\(Ans=\sum_{p\in prim}\sum_{p|d}\mu(\lfloor\frac{d}{p}\rfloor)F(d)\)

由于此时存在一个\(\mu(\lfloor \frac dp \rfloor)\) 难以化简所以我们转而枚举 \(\lfloor \frac d p \rfloor\)得到

\(Ans=\sum_{p\in prim}\sum_{d=1}^{min(\lfloor\frac{n}{p}\rfloor,\lfloor\frac{m}{p}\rfloor)}\mu(d)F(dp)=\sum_{p\in prim}\sum_{d=1}^{min(\lfloor\frac{n}{p}\rfloor,\lfloor\frac{m}{p}\rfloor)}\mu(d)\lfloor\frac{n}{dp}\rfloor\lfloor\frac{m}{dp}\rfloor\)

由于现在需要枚举d和p两个数我们观察最后两项\(\lfloor\frac{n}{dp}\rfloor\lfloor\frac{m}{dp}\rfloor\)

可以直接枚举\(T=d*p\)从而得到一个更简单的式子

\(Ans=\sum_{T=1}^{min(n,m)}\sum_{t|T,t\in prime}\mu(\lfloor\frac{T}{t}\rfloor)\lfloor\frac{n}{T}\rfloor\lfloor\frac{m}{T}\rfloor\)

再次化简将最后两项\(\lfloor\frac{n}{dp}\rfloor\lfloor\frac{m}{dp}\rfloor\)提到前边即

\(Ans=\sum_{T=1}^{min(n,m)}\lfloor\frac{n}{T}\rfloor\lfloor\frac{m}{T}\rfloor(\sum_{t|T,t\in prime}\mu(\lfloor\frac{T}{t}\rfloor))\)

然后这样我们就已经化简成功了

此时显然我们可以直接O(n)求解但是对于多组数据我们必须找到一种更优秀的写法来求解这个式子。

首先我们可以对于最后面的\(\sum_{t|T,t\in prime}\mu(\lfloor\frac{T}{t}\rfloor)\)求一个前缀和设为SUM

然后我们来看这个式子\(Ans=\sum_{T=1}^{min(n,m)}\lfloor\frac{n}{T}\rfloor\lfloor\frac{m}{T}\rfloor SUM\)

那么此时我们观察前两项这不就是裸的就是整除分块吗?

然后我们预处理出SUM在对前面的两项做一个整除分块时间复杂度\(O(T\sqrt n)\)

代码

#include<bits/stdc++.h>
#define N 10000100
#define ll long long
using namespace std;
inline int read() {
    int n=1,num=0;char ch=getchar();
    while(!isdigit(ch)) {n=(ch=='-')?-1:1;ch=getchar();}
    while(isdigit(ch))  {num=(num<<1)+(num<<3)+(ch^48);ch=getchar();}
    return n*num;
}
int vis[N],prime[N],mu[N],g[N],tot;//g表示推出式子最后的Σ(t/d) ,sum为其前缀和 
ll sum[N];
void get_mu(int n)
{
    mu[1]=1;
    for(int i=2;i<=n;++i){
        if(!vis[i]) {mu[i]=-1;prime[++tot]=i;}
        for(int j=1;j<=tot&&prime[j]*i<=n;++j) {
            vis[i*prime[j]]=1;
            if(!(i%prime[j]))   break;
            else                mu[i*prime[j]]=-mu[i];
        }
    }
    for(int j=1;j<=tot;++j)
        for(int i=1;i*prime[j]<=n;++i)  g[i*prime[j]]+=mu[i];
    for(int i=1;i<=n;++i)   sum[i]=sum[i-1]+(ll)g[i];
}
int n,m,t;
int main()
{
    t=read();
    get_mu(10000000);
    while(t--) {
        n=read();m=read();
        if(n>m) swap(n,m);
        ll ans=0;
        for(int l=1,r;l<=n;l=r+1) {
            r=min(n/(n/l),m/(m/l));
            ans+=1ll*(n/l)*(m/l)*(sum[r]-sum[l-1]);
        }
        printf("%lld\n",ans);
    }
    return 0;
}

转载于:https://www.cnblogs.com/My-snowing/p/10393905.html

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

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值