题解 P2257 【YY的GCD】

· · 题解

在这里提供一个新思路

想必大家都做过这个题吧

\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=1]

有一种做法是利用\mu函数的一个性质

\sum_{d\mid n}\mu(d)=[n=1]

从而可以O(\sqrt n+\sqrt m)时间内求解答案

那么这个题让我们求

\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)\ is\ prime]

怎么办？？

如果有个函数(设为f(n))它如果能有以下性质就好了

\sum_{d\mid n}f(d)=[n\ is \ prime]

没有还不能自己造嘛？

应该可以通过莫比乌斯反演找到另外一个数论函数求得我们的这个f(n),但找不到怎么办？？(其实是我不会)

有一种方法可以通过一些加加减减来获得这个函数，就是我先让f(n)等于我们想要的那个值，比如我们要构造的这个函数就这样附上初值

f(n)=\begin{cases}1&n\ is\ prime\\0& n\ is\ composite \end{cases}

其核心就是再让

f(n)=f(n)-\sum_{d\mid n,d\neq n}f(d)

这样的话不就有

\sum_{d\mid n}f(d)=f(n)+\sum_{d\mid n,d\neq n}f(d)=\begin{cases}1& n\ is\ prime\\0&n\ is\ composite\end{cases}

（除了f(n)的初值以外其他前面的函数值的和都被抵消掉了）

然后附上实现第二步的代码(赋初值的代码就省略掉了)

for (int i=1;i<=MAX;++i)
    for (int j=i<<1;j<=MAX;j+=i)
        f[j]-=f[i];

这两层for循环第二维是一个调和级数，总复杂度约为nln(n)

这种方法也是比较常用的吧，有些题目里面也需要自己用这个方法构造函数，而且它的常数也比直接用莫比乌斯反演还小(我们学校的同学竟然出了卡莫比乌斯反演的数据，让我们只能用这种方法求)

那么这个题就愉快地结束了，除了复杂度是O(nln(n)+T(\sqrt N+\sqrt M))在4秒时限下跑了3.6秒勉强卡过这道题外，其他都还是很好的

最后感谢DT_Kang神犇的指导，以及附上AC代码

#include<cstdio>
#define N 10230040
#define U 10000000
#define ll long long
int min(int a,int b) {return a<b?a:b;}
long long f[N];int prime[N/10],cnt;bool ip[N];//其实10000000内只有660000个质数
int read(){//丑陋的快读
    int n=0;char a;bool z=false;
    while(a=getchar())
    {
        if (a>'9'||a<'0')
            if (z) break;
            else continue;
        if (!z) z=true;
        n=(n<<1)+(n<<3)+(a^48);
    }
    return n;
}
int main()
{
    int n,T,m,up,R;ll ans=0;
    for (int i=2;i<=U;++i)
    {
        if (!ip[i])
        {
            prime[++cnt]=i;
            f[i]=1;
        }
        for (int j=1;j<=cnt&&i*prime[j]<=U;++j)
        {
            ip[i*prime[j]]=true;
            if (i%prime[j]==0) break;
        }
    }
    for (int i=2;i<=U;++i)
        for (int j=i<<1;j<=U;j+=i)
            f[j]-=f[i];
    for (int i=2;i<=U;++i) f[i]+=f[i-1];
    T=read();
    while(T--)
    {
        ans=0;
        n=read();m=read();up=min(n,m);
        for (int i=1;i<=up;i=R+1)
        {
            R=min(n/(n/i),m/(m/i));
            ans+=1ll*(f[R]-f[i-1])*(n/i)*(m/i);
        }
        printf("%lld\n",ans);
    }
    return 0;
}