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

测试地址：Visible Lattice Points
题目大意：在三维空间中，我们说一个点是可见的当且仅当它与点(0,0,0)连成的线段不经过任何其他坐标为整数的点。有T(T≤50)组询问，每组询问给出一个参数N，意为询问在所有点(x,y,z)(0≤x,y,z≤N,1≤N≤1∗106)中，有多少个可见的点，其中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）：

∑x=1N∑y=1N∑z=1N[gcd(x,y,z)=1]

显然暴力求这个式子会炸，这时候我们就要用到一个莫比乌斯函数的性质：∑d|nμ(d)=[n=1]，这个性质使得我们可以把上式化成：
∑Nx=1∑Ny=1∑Nz=1∑d|gcd(x,y,z)μ(d)。
分析第四个求和号下面的条件，显然d|gcd(x,y,z)的充要条件是d|x且d|y且d|z，那么上式就可以化成：
∑Nd=1μ(d)∑1≤x≤N且d|x∑1≤y≤N且d|y∑1≤z≤N且d|z1
这个式子等价于：
∑Nd=1μ(d)(∑1≤x≤N且d|x1)×(∑1≤y≤N且d|y1)×(∑1≤z≤N且d|z1)
那么这个式子显然等于：
∑Nd=1μ(d)⌊Nd⌋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)的总复杂度对于严苛的时间限制好像有点拙计，需要想办法优化。注意到⌊Nd⌋的值不会超过2N−−√种，而且对于同一个⌊Nd⌋的值，满足条件的d是一个连续的区间，这启发我们利用分块思想优化，问题转变为如何求⌊Nd⌋值相等的区间。思考对于任意一个d，求满足⌊Nd⌋=⌊Np⌋的最大的p。令q=⌊Nd⌋，那么N=dq+r0=pq+r1，我们知道在q相等的情况下，r1越逼近0，p就越大，所以最大的p=⌊Nq⌋=⌊N⌊Nd⌋⌋。那么我们再预处理出μ(i)的前缀和就可以把单个询问的时间复杂度加速到O(N−−√)，处理询问的总时间复杂度为O(TN−−√)，可以通过此题。
以下是本人代码：

#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;
}