bzoj4407 于神之怒加强版

题目链接：bzoj4409
题目大意：
给N,M,K。求

$$\sum_{i=1}^{n}\sum_{j=1}^{m}gcd(i,j)^k \mod (10^9+7)$$
多组数据
$1<=N,M,K<=5000000,1<=T<=2000$
时限80s，空间512M

题解：
莫比乌斯反演
设$gcd(i,j)=d$
则原式看成$\sum_{i=1}^{n}\sum_{j=1}^{m}d^k$
交换枚举约数和倍数

$$\sum_{d=1}^{min(n,m)}\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}d^k[gcd(i,j)=1]$$
莫比乌斯反演一下
$$\sum_{d=1}^{min(n,m)}d^k\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}\sum_{t|(i,j)}\mu(t)$$
继续交换枚举约数和倍数 $$\sum_{d=1}^{min(n,m)}d^k\sum_{t=1}^{min(\lfloor\frac{n}{d}\rfloor,\lfloor\frac{m}{d}\rfloor)}\sum_{i=1}^{\lfloor\frac{n}{dt}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{dt}\rfloor}\mu(t)$$
即 $$\sum_{d=1}^{min(n,m)}d^k\sum_{t=1}^{min(\lfloor\frac{n}{d}\rfloor,\lfloor\frac{m}{d}\rfloor)}\lfloor\frac{n}{dt}\rfloor\times \lfloor\frac{m}{dt}\rfloor\times \mu(t)$$
too naive的我看到80s时限化到这里就直接敲代码了
预处理 $d^k$的前缀和、 $\mu$的前缀和然后分两次块什么的
于是完美TLE。估计排队的人骂死我
[后知后觉：。。。原来我当成T≤20来算时间了

瞄了一眼别人的题解
还要继续套路下去啊。
换元，设$T=dt$
则原式化为

$$\sum_{T=1}^{min(n,m)}\lfloor\frac{n}{T}\rfloor\times \lfloor\frac{m}{T}\rfloor\sum_{d|T}d^k\mu(\frac{T}{d})$$
设函数 $f(T)=\sum_{d|T}d^k\mu(\frac{T}{d})$,即 $f=Id_k*\mu$
单位幂函数 $Id_k$与莫比乌斯函数 $\mu$都是积性函数
根据积性函数的性质， $f$也是积性函数
考虑线性筛把$f$搞出来

• 若$i$为质数，则$f[i]=\mu(1)*i^k+\mu(i)*1=i^k-1$
• 若$i$为合数，
  
  • 若$i与prime互质$，则$f[i*prime]=f[i]*f[prime]$
  • 若$i与prime不互质$，据$\mu$的定义，含平方因子的$\mu$值为0可得，$f[i*prime]=f[i]*prime^k$…当我在口胡= =看dalao的吧戳我

  然后预处理埋它的前缀和
  那么

  $$Ans=\sum_{T=1}^{min(n,m)}\lfloor\frac{n}{T}\rfloor\times \lfloor\frac{m}{T}\rfloor\times f(T)$$
  此时再分块就好了求就好啦

  #include<cstdio>
  #include<cstdlib>
  #include<cstring>
  #include<iostream>
  #include<algorithm>
  using namespace std;
  typedef long long LL;
  #define N 5000100
  
  const LL mod=1000000007;
  LL f[N];bool ispri[N];int cnt,K,pri[N/10];
  int mymin(int x,int y){return (x<y)?x:y;}
  LL qpow(LL x,int t)
  {
      LL ret=1;
      while (t)
      {
          if (t&1) ret=(ret*x)%mod;
          x=(x*x)%mod;t>>=1;
      }return ret;
  }
  void pre(int lim)
  {
      cnt=0;f[1]=1;
      for (int i=2;i<=lim;i++)
      {
          if (!ispri[i]) {pri[++cnt]=i;f[i]=(qpow((LL)i,K)-1+mod)%mod;}
          for (int j=1;j<=cnt && i*pri[j]<=lim;j++)
          {
              ispri[i*pri[j]]=true;
              if (i%pri[j]==0)
              {
                  f[i*pri[j]]=f[i]*qpow((LL)pri[j],K)%mod;
                  break;
              }f[i*pri[j]]=(f[i]*f[pri[j]])%mod;
          }
      }
      for (int i=2;i<=lim;i++) f[i]=(f[i-1]+f[i])%mod;
  }
  int main()
  {
      int T,lim,n,m,i,r;LL ans;
      scanf("%d%d",&T,&K);
      lim=5000000;pre(lim);
      while (T--)
      {
          scanf("%d%d",&n,&m);ans=0;
          if (n>m){int tt=n;n=m;m=tt;}
          for (i=1;i<=n;i=r+1)
          {
              r=mymin(n/(n/i),m/(m/i));
              ans+=(LL)(n/i)*(LL)(m/i)%mod*(f[r]-f[i-1])%mod;
              ans%=mod;
          }
          if (ans<0) ans+=mod;
          printf("%lld\n",ans);
      }
      return 0;
  }