51nod 2583 数论只会Gcd

Link
对于一组询问\((x,y)\)，若\(x>m\vee y>m\)，那么答案一定为\(0\)。
注意到\(\gcd(n,m)\)操作进行一次之后一定会有\(m<n\)，因此若\(x\le y\)，那么答案一定为\(1\)。
因此我们就只需要考虑\(x>y\)的情况了。
考虑建一棵树，每个节点代表一个二元组\((a,b)(b<a\le m)\)，根节点为\((x,y)\)，每个节点的父亲是它进行一次\(\gcd\)操作之后变成的二元组。那么当前询问的答案就是树的大小乘以二。
先将深度为偶数（根节点深度为\(1\)）的点的二元组交换两维。
然后用左儿子-右兄弟的原则将其转为二叉树，那么我们可以得到一棵这样的二叉树：

考虑用二元组\((a,b)\)来表示\(ax+b(x+y)\)。
忽略掉当前的树根\((x,y)\)，那么我们现在有\(((1,0),(0,1))\)一个四元组。
对于四元组\(((a,b),(c,d))\)，它有两个儿子：\(((a+c,b+d),(c,d)),((a,b),(a+c,b+d))\)。
不难发现这样的四元组去掉\(((1,0),(0,1))\)之后与Stern Brocot树上的相邻兄弟节点对一一对应。
对于一个四元组\(((a,b),(c,d))\)，考虑将其在\((a,b),(c,d)\)中较大的一方处统计，即对应的Stern Brocot树上的相邻兄弟节点对中的实点。
也就是说，对于Stern Brocot树中的每个实点，它都会对应两个四元组。
一个实点出现在Stern Brocot树中的充要条件是其对应二元组\((a,b)\)满足\(a\perp b\wedge a,b>0\)。
现在我们还有的限制条件就是\(ax+b(x+y)\le m\)。

\[\begin{aligned} ans&=\sum\limits_{a\ge1}\sum\limits_{b\ge1}[a\perp b][ax+b(x+y)\le m]\\ &=\sum\limits_{d\ge1}^m\mu(d)\sum\limits_{a\ge1}\sum\limits_{b\ge1}[ax+b(x+y)\le\lfloor\frac md\rfloor] \end{aligned} \]

数论分块+杜教筛即可。

#include<cstdio>
using i64=long long;
const int N=1001007;
int n,lim=1000000,cnt,mu[N],sum_mu[N],pos[N];
int read(){int x;scanf("%d",&x);return x;}
int getid(int x){return x<=lim? x:cnt+1-n/x;}
void Euler_Sieve()
{
    static int np[N],pr[N],tot;mu[1]=1;
    for(int i=2;i<=lim;++i)
    {
	if(!np[i]) pr[++tot]=i,mu[i]=-1;
	for(int j=1;j<=tot&&i*pr[j]<=lim;++j) if(np[i*pr[j]]=1,mu[i*pr[j]]=-mu[i],!(i%pr[j])) {mu[i*pr[j]]=0;break;}
    }
    for(int i=1;i<=lim;++i) sum_mu[i]=sum_mu[i-1]+mu[i];
}
void Xudyh_Sieve()
{
    for(int i=1;i<=lim;++i) pos[++cnt]=i;
    if(n>=lim) for(int l=lim+1,r;l<=n;l=r+1) r=n/(n/l),pos[++cnt]=r,sum_mu[cnt]=1;
    for(int i=lim+1;i<=cnt;++i) for(int l=2,r;l<=pos[i];l=r+1) r=pos[i]/(pos[i]/l),sum_mu[i]-=(r-l+1)*sum_mu[getid(pos[i]/l)];
}
i64 calc(i64 n,i64 a,i64 b,i64 c)
{
    i64 m=(a*n+b)/c;
    if(n<0) return 0;
    if(!a) return (b/c)*(n+1);
    if(a>=c||b>=c) return n*(n+1)/2*(a/c)+(n+1)*(b/c)+calc(n,a%c,b%c,c);
    return n*m-calc(m-1,c,c-b-1,a);
}
i64 calc(int x,int y,int n){return calc(n/x-1,x,n%x,y);}
i64 solve()
{
    int x=read(),y=read();i64 ans=1;
    if(x>n||y>n) return 0;
    if(x<=y) return 1;
    if(y+=x,y>n) return 2;
    for(int l=1,r;l<=n;l=r+1) r=n/(n/l),ans+=(sum_mu[getid(r)]-sum_mu[getid(l-1)])*calc(x,y,n/l);
    return 4*ans;
}
int main()
{
    int t=read();
    n=read(),Euler_Sieve(),Xudyh_Sieve();
    while(t--) printf("%lld\n",solve());
}