BZOJ_P3561 DZY Loves Math VI(数论+莫比乌斯反演)

来源:互联网 发布:java图形界面时间控件 编辑:程序博客网 时间:2024/04/30 14:42

BZOJ传送门
Time Limit: 10 Sec Memory Limit: 256 MB
Submit: 90 Solved: 66
[Submit][Status][Discuss]
Description
给定正整数n,m。求
这里写图片描述
Input
一行两个整数n,m。

Output
一个整数,为答案模1000000007后的值。

Sample Input
5 4

Sample Output
424

HINT
数据规模:
1<=n,m<=500000,共有3组数据。

Source
By Jcvb

#include<cstdio>#include<iostream>using namespace std;#define N 500005#define P 1000000007typedef long long ll;bool b[N];int mu[N],p[N];int n,m,cnt,ans;int sum[N],a[N];inline int Pow(int a,int b,int res=1){while(b){if(b&1) res=(ll)res*a%P;a=(ll)a*a%P,b>>=1;}return res;}void GetMu(){    mu[1]=1;    for(int i=2;i<=n;i++){        if(!b[i]) mu[i]=-1,p[cnt++]=i;        for(int j=0;j<cnt&&i*p[j]<=n;j++){            b[i*p[j]]=1;            if(i%p[j]) mu[i*p[j]]=-mu[i];            else break;        }    }}int main(){    scanf("%d%d",&n,&m);if(n<m) swap(n,m);GetMu();    for(int i=1;i<=n;i++) a[i]=1;    for(int d=1;d<=m;d++){        int x=Pow(d,d),y=0,limm=m/d,limn=n/d;        for(int j=1;j<=limn;j++) a[j]=(ll)a[j]*j%P,sum[j]=(sum[j-1]+a[j])%P;        for(int k=1;k<=limm;k++)if(mu[k])            y=((ll)a[k]*a[k]%P*sum[limm/k]%P*sum[limn/k]%P*mu[k]+y+P)%P;        ans=(ans+(ll)x*y%P)%P;    }    printf("%d",ans);    return 0;}
1 0
原创粉丝点击