[bzoj 3994] sdoi2015 约数个数和

需要用到一个神奇的式子：

d(nm)=sigma [ gcd(i,j)==1 ] 其中i|n&&j|m

把这个做前缀和应用到最终答案中去可得

ans = sigma（ {n/i} * {m/j}  * [ gcd(i,j)==1 ] ） 其中{x}表示x向下取整

应用莫比乌斯反演可得

ans = sigma（ {n/i} * {m/j}  *  sigma mu(d) ）其中d|i&&d|j

ans =  sigma （ mu(d) * d[n/i] * d[m/i] ) 其中d[i]为i的约束个数

后面就不用我说了吧。。。。

#include <cstdio>
#include <cstring>
#include <algorithm>
 
using namespace std;
const int Maxn=500005;
typedef long long LL;
LL mu[Maxn],d[Maxn],ans;
int prime[Maxn],cnt[Maxn];
int n,m,i,j,tot,Case;
bool check[Maxn];
 
void init(){
  mu[1] = 1; d[1] = 1;
  for (i=2;i<Maxn;i++){
    if (!check[i]){
      prime[++tot]=i; mu[i]=-1;
      d[i]=2; cnt[i]=1;
    }
    for (j=1;j<=tot;j++){
      if (i*prime[j]>=Maxn) break;
      check[i*prime[j]]=1;
      if (i%prime[j]){
        mu[i*prime[j]] = -mu[i];
        d[i*prime[j]] = (d[i]<<1);
        cnt[i*prime[j]] = 1;
      }else
      {
        mu[i*prime[j]] = 0;
        d[i*prime[j]] = d[i]/(cnt[i]+1)*(cnt[i]+2);
        cnt[i*prime[j]] = cnt[i]+1;
        break;
      }
    }
  }
  for (i=1;i<Maxn;i++){
    mu[i] += mu[i-1];
    d[i] += d[i-1];
  }
}
 
int main(){
  //freopen("3994.in","r",stdin);
  //freopen("3994.out","w",stdout);
  init();
  scanf("%d",&Case);
  while (Case--){
    scanf("%d%d",&n,&m);
    ans = 0;
    for (i=1;i<=n&&i<=m;i=j+1){
      j = min( n/(n/i), m/(m/i) );
      ans += (mu[j]-mu[i-1])*d[n/i]*d[m/i];
    }
    //for (i=1;i<=n;i++)
    //  for (j=1;j<=m;j++)
    //  if (gcd(i,j)==1) ans += (n/i)*(m/j);
      //  for (int i1=1;i1<=i;i1++) if (i%i1==0)
      //    for (int j1=1;j1<=j;j1++) if (j%j1==0)
      //    if (gcd(i1,j1)==1) ans++;
    printf("%lld\n",ans);
  }
  return 0;
}