莫比乌斯进阶：bzoj 4407 于神之怒加强版（Mobius）

传送门

Tip：那两对很高的括号也是向下取整
P.S.现在才真正理解到原创的辛苦orz

#include<bits/stdc++.h>using namespace std;typedef long long ll;const int MAXN=5e6+2;const ll MOD=1e9+7;int n,m,kase,prime[MAXN],tot=0;ll f[MAXN],g[MAXN],k;bool vis[MAXN];ll fpow(ll a,ll b,ll p) {    ll ret=1;    while (b) {        if (b&1) ret=ret*a%p;        b>>=1,a=a*a%p;    }    return ret;}inline void linear_shaker() {    memset(vis,false,sizeof(vis));    f[1]=1;    for (register int i=2;i<MAXN;++i) {        if (!vis[i]) prime[++tot]=i,f[i]=(g[i]=fpow((ll)i,k,MOD))-1;//如果每次不记录fpow的值运行时间会变成1.5倍，能过        for (int j=1;j<=tot&&i*prime[j]<MAXN;++j) {            vis[i*prime[j]]=true;            if (i%prime[j]==0) {f[i*prime[j]]=f[i]*g[prime[j]]%MOD;break;}//g[x]=x^k%MOD            f[i*prime[j]]=f[i]*f[prime[j]]%MOD;        }    }    for (register int i=2;i<MAXN;++i) f[i]=(f[i]+f[i-1])%MOD;}int main() {//  freopen("bzoj 4407.in","r",stdin);    scanf("%d%d",&kase,&k);    linear_shaker();    while (kase--) {        scanf("%d%d",&n,&m);        if (n>m) n^=m^=n^=m;        int last=0;        ll ans=0;        for (int i=1;i<=n;i=last+1) {            last=min(n/(n/i),m/(m/i));            ans=(ans+(f[last]-f[i-1]+MOD)*(n/i)%MOD*(m/i)%MOD)%MOD;//小心负数        }        printf("%lld\n",ans);    }    return 0;}