#include"bits/stdc++.h"
#define last last_
using namespace std;
typedef long long LL;
const int P=1000000007;
struct node{
int n,m;
inline LL hash(){return (LL)n<<30|m;}
node(int _n=0,int _m=0){n=_n,m=_m;}
};
const int N=5000005,M=500005,hmod=1000003,L=5000000,X=100005;
struct HashMap{
node w[M];
int val[M],flagnull,h[M],nxt[M],ndc;
HashMap(){flagnull=-1,ndc=0;}
int&find(node x){
for(int i=h[x.hash()%hmod];i;i=nxt[i])
if(x.n==w[i].n&&x.m==w[i].m)return val[i];
return flagnull;
}
int&operator[](node x){
int&r=find(x);
if(r!=-1)return r;
w[++ndc]=x,nxt[ndc]=h[x.hash()%hmod],h[x.hash()%hmod]=ndc;
return val[ndc];
}
}q;
int n,m,mu[X],eu[N],seu[N],ed[X],pri[N>>3],pcnt,ans;
bool fl[N];
vector<int>di[X];
int G(int n){
if(n<=L)return seu[n];
node x(1,n);
int r=q.find(x);
if(r!=-1)return r;
r=((LL)n*(n+1)>>1)%P;
for(int i=1,j;i<=n;i=j+1){
j=n/(n/i);
if(n/i==1)continue;
r=(r-(LL)(n/i-1)*(G(j)-G(i-1)))%P;
}
return q[x]=(r+P)%P;
}
int F(int n,int m){
if(n==1)return G(m);
if(m==0)return 0;
if(mu[n]==0)return (LL)F(ed[n],m)*(n/ed[n])%P;
node x(n,m);
int r=q.find(x);
if(r!=-1)return r;
r=0;
for(int i=0,d;i<di[n].size();i++){
d=di[n][i];
r=(r+(LL)eu[n/d]*F(d,m/d))%P;
}
return q[x]=r;
}
int main(){
#ifndef ONLINE_JUDGE
freopen("bzoj3512.in","r",stdin);
freopen("bzoj3512.out","w",stdout);
#endif
scanf("%d%d",&n,&m);
mu[1]=ed[1]=eu[1]=seu[1]=1;
for(int i=2;i<=L;i++){
if(i<=n)di[i].push_back(1);
if(!fl[i]){
pri[++pcnt]=i;
if(i<=n){
ed[i]=i;
mu[i]=-1;
}
eu[i]=i-1;
}
for(int j=1;i*pri[j]<=L;j++){
fl[i*pri[j]]=1;
if(i%pri[j]==0){
if(i*pri[j]<=n){
mu[i*pri[j]]=0;
ed[i*pri[j]]=ed[i];
}
eu[i*pri[j]]=eu[i]*pri[j];
break;
}
if(i*pri[j]<=n){
mu[i*pri[j]]=-mu[i];
ed[i*pri[j]]=ed[i]*pri[j];
}
eu[i*pri[j]]=eu[i]*(pri[j]-1);
}
for(int j=1;i*j<=n;j++)di[i*j].push_back(i);
seu[i]=(seu[i-1]+eu[i])%P;
}
for(int i=1;i<=n;i++)ans=(ans+F(i,m))%P;
cout<<ans<<endl;
return 0;
}
[BZOJ3512] DZY Loves Math IV -- 杜教筛,搜索
最新推荐文章于 2020-01-21 11:42:19 发布