bzoj3561 DZY Loves Math VI 莫比乌斯函数

       对原式运用莫比乌斯函数进行转化得到原式=

      然后直接暴力就好了。时间复杂度O(NlogN)。

AC代码如下：

#include<iostream>
#include<cstdio>
#include<cstring>
#define ll long long
#define p 1000000007
#define N 500005
using namespace std;

int m,n,cnt,c[N],mu[N],a[N],sum[N]; bool vis[N];
int ksm(int x,int y){
	int t=1; for (; y; y>>=1,x=(ll)x*x%p) if (y&1) t=(ll)t*x%p; return t;
}
void pfs(){
	int i,j; mu[1]=1;
	for (i=2; i<=n; i++){
		if (!vis[i]){ mu[i]=-1; c[++cnt]=i; }
		for (j=1; j<=cnt && i*c[j]<=n; j++){
			vis[i*c[j]]=1;
			if (i%c[j]) mu[i*c[j]]=-mu[i]; else{
				mu[i*c[j]]=0; break;
			}
		}
	}
}
int main(){
	scanf("%d%d",&m,&n); int i,j,ans=0;
	if (m>n) swap(m,n); pfs();
	for (i=1; i<=n; i++) a[i]=1;
	for (i=1; i<=m; i++){
		int x=ksm(i,i),y=0;
		for (j=1; j*i<=n; j++){
			a[j]=(ll)a[j]*j%p; sum[j]=(sum[j-1]+a[j])%p;
		}
		for (j=1; j*i<=m; j++) if (mu[j])
			y=((ll)a[j]*a[j]%p*sum[m/i/j]%p*sum[n/i/j]%p*mu[j]+y+p)%p;
		ans=(ans+(ll)x*y%p)%p;
	}
	printf("%d\n",ans);
	return 0;
}

by lych

2016.2.23