【SPOJ7001】Visible Lattice Points-莫比乌斯反演+分块

测试地址：Visible Lattice Points
题目大意：在三维空间中，我们说一个点是可见的当且仅当它与点$(0,0,0)$连成的线段不经过任何其他坐标为整数的点。有$T(T\le50)$组询问，每组询问给出一个参数$N$，意为询问在所有点$(x,y,z)(0\le x,y,z\le N,1\le N\le 1*10^6)$中，有多少个可见的点，其中$x,y,z$为整数。特别地，$(0,0,0)$不算作可见的点。
做法：这一道题是POJ3090的加强版，从二维扩展到了三维，我写的POJ3090题解在此。而这题就不能简单地使用欧拉函数的性质来解决问题了，而需要用莫比乌斯反演来解决。
将点分为三种情况：1.$x,y,z$均不为0。2.$x,y,z$中只有一个是0。3.$x,y,z$中只有一个不是0。对于第一种情况，很容易看出当$\gcd(x,y,z)=1$时，点$(x,y,z)$是可见的，那么我们就是要求下面这个式子的值（注意，下文中方括号$[]$表示如果括号内的式子为真，值为1，否则值为0）：

$$\sum_{x=1}^{N}\sum_{y=1}^{N}\sum_{z=1}^{N}[\gcd(x,y,z)=1]$$
显然暴力求这个式子会炸，这时候我们就要用到一个莫比乌斯函数的性质： $\sum_{d|n}\mu(d)=[n=1]$，这个性质使得我们可以把上式化成：
$\sum_{x=1}^{N}\sum_{y=1}^{N}\sum_{z=1}^{N}\sum_{d|\gcd(x,y,z)}\mu(d)$。
分析第四个求和号下面的条件，显然 $d|\gcd(x,y,z)$的充要条件是 $d|x且d|y且d|z$，那么上式就可以化成：
$\sum_{d=1}^{N}\mu(d)\sum_{1\le x\le N且d|x}\sum_{1\le y\le N且d|y}\sum_{1\le z\le N且d|z}1$
这个式子等价于：
$\sum_{d=1}^{N}\mu(d)(\sum_{1\le x\le N且d|x} 1)\times(\sum_{1\le y\le N且d|y}1)\times(\sum_{1\le z\le N且d|z}1)$
那么这个式子显然等于：
$\sum_{d=1}^{N}\mu(d)\lfloor\frac{N}{d}\rfloor^3$
这就是第一种情况的答案了。第二种情况就是要分别求 $gcd(x,y)=1,gcd(y,z)=1,gcd(x,z)=1$的点数，式子的推导和上面比较类似，这里就不再赘述了。第三种情况显然只有 $(0,0,1),(0,1,0),(1,0,0)$三个点。这样我们就找到了一个可以 $O(N)$计算单个询问答案的方法。但是 $O(TN)$的总复杂度对于严苛的时间限制好像有点拙计，需要想办法优化。注意到 $\lfloor\frac{N}{d}\rfloor$的值不会超过 $2\sqrt N$种，而且对于同一个 $\lfloor\frac{N}{d}\rfloor$的值，满足条件的 $d$是一个连续的区间，这启发我们利用分块思想优化，问题转变为如何求$\lfloor\frac{N}{d}\rfloor$值相等的区间。思考对于任意一个 $d$，求满足$\lfloor\frac{N}{d}\rfloor=\lfloor\frac{N}{p}\rfloor$的最大的 $p$。令$q=\lfloor\frac{N}{d}\rfloor$，那么 $N=dq+r_0=pq+r_1$，我们知道在 $q$相等的情况下，$r_1$越逼近 $0$，$p$就越大，所以最大的 $p=\lfloor\frac{N}{q}\rfloor=\lfloor\frac{N}{\lfloor\frac{N}{d}\rfloor}\rfloor$。那么我们再预处理出 $\mu(i)$的前缀和就可以把单个询问的时间复杂度加速到 $O(\sqrt N)$，处理询问的总时间复杂度为 $O(T\sqrt N)$，可以通过此题。
以下是本人代码：

#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <iostream>
#include <algorithm>
#define ll long long
using namespace std;
int T;
ll N,mu[1000010],sum[1000010];
bool prime[1000010]={0};

void calc_mu()
{
  N=1000000;
  for(int i=1;i<=N;i++) mu[i]=1;
  for(int i=2;i<=N;i++)
    if (!prime[i])
    {
      for(int j=1;i*j<=N;j++)
      {
        mu[i*j]*=-1;
        if (j>1) prime[i*j]=1;
        if (!(j%i)) {mu[i*j]=0;continue;}
      }
    }
  sum[0]=0;
  for(int i=1;i<=N;i++) sum[i]=sum[i-1]+mu[i];
}

int main()
{
  scanf("%d",&T);
  calc_mu();
  while(T--)
  {
    scanf("%lld",&N);
    ll ans=3,last;
    for(ll d=1;d<=N;d=last+1)
    {
      last=N/(N/d);
      ans+=((N/d)+3)*(N/d)*(N/d)*(sum[last]-sum[d-1]);
    }
    printf("%lld\n",ans);
  }

  return 0;
}