bzoj4407 于神之怒加强版 莫比乌斯函数

首先运用莫比乌斯函数对原式进行化简，得到：

$$原式=\sum_{i=1}^{m}\sum_{j=1}^{n} gcd(i,j)^{k} \rightarrow \sum_{d} d^{k} \sum_{i=1}^{\left \lfloor m/d \right \rfloor}\sum_{j=1}^{\left \lfloor n/d \right \rfloor} [gcd(i,j)==1]$$

$$\rightarrow \sum_{d} d^{k} \sum_{i=1}^{\left \lfloor m/d \right \rfloor}\sum_{j=1}^{\left \lfloor n/d \right \rfloor}\sum_{p|i,p|j}\mu(p) \rightarrow \sum_ {k} \left \lfloor m/k \right \rfloor \left \lfloor n/k \right \rfloor \sum_ {d|k} d^{k} \mu(k/d)$$

然后，令$f(x)=\sum_ {d|k} d^{k} \mu(k/d)$，显然可以得到$f(x)$是一个积性函数，如果令辅助数组$g(x)$为x的最小质因数的最高次幂，就可以利用线性筛$O(N)$得到$f$数组。
然后就可以分段计算了。时间复杂度$O(Max+T\sqrt{N})$。
AC代码如下：

#include<iostream>
#include<cstdio>
#define ll long long
#define N 5000005
#define mod 1000000007
using namespace std;

int f[N],g[N],c[N]; bool vis[N];
int ksm(int x,int y){
    int t=1; for (; y; y>>=1,x=(ll)x*x%mod) if (y&1) t=(ll)t*x%mod;
    return t;
}
int main(){
    int cas,k,i,j; f[1]=1;
    scanf("%d%d",&cas,&k);
    for (i=2; i<=5000000; i++){
        if (!vis[i]){ c[++c[0]]=g[i]=i; f[i]=(ksm(i,k)+mod-1)%mod; }
        for (j=1; j<=c[0] && i*c[j]<=5000000; j++){
            int x=i*c[j]; vis[x]=1;
            if (i%c[j]){
                g[x]=c[j]; f[x]=(ll)f[i]*f[c[j]]%mod;
            } else{
                g[x]=g[i]*c[j];
                if (g[i]==i) f[x]=(ll)f[i]*(f[c[j]]+1)%mod; else f[x]=(ll)f[i/g[i]]*f[g[i]*c[j]]%mod;
                break;
            }
        }
    }
    for (i=2; i<=5000000; i++) f[i]=(f[i]+f[i-1])%mod;
    while (cas--){
        int m,n,ans=0; scanf("%d%d",&m,&n);
        for (i=1; i<=m && i<=n; i=j+1){
            j=min(m/(m/i),n/(n/i));
            ans=((ll)(f[j]-f[i-1]+mod)*(m/i)%mod*(n/i)%mod+ans)%mod;
        }
        printf("%d\n",ans);
    }
    return 0;
}

                                        by lych
                                       2016.2.16