bzoj 3518 Dirichlet卷积
详情见代码，回头再填坑。。。
1 #include<iostream> 2 #include<cstdio> 3 #include<algorithm> 4 #include<cstring> 5 #define int long long 6 #define p 1000000007 7 using namespace std; 8 int n,m; 9 int phi[50005],su[50005],pr[50005],cnt; 10 void shai() 11 { 12 phi[1]=1; 13 for(int j=2;j<=50000;j++) 14 { 15 if(!phi[j])phi[j]=j-1,su[++cnt]=j,pr[j]=j; 16 for(int i=1;su[i]<=pr[j]&&i<=cnt&&su[i]*j<=50000;i++) 17 { 18 pr[su[i]*j]=su[i]; 19 if(!phi[su[i]*j])phi[su[i]*j]=su[i]*j; 20 if(su[i]==pr[j])phi[su[i]*j]=phi[j]*su[i]; 21 else phi[su[i]*j]=phi[j]*(su[i]-1); 22 } 23 } 24 } 25 signed main() 26 { 27 shai(); 28 scanf("%lld%lld",&n,&m); 29 int ans=0; 30 ans+=n*(n-1)*(n-2)*m/6;ans%=p; 31 ans+=m*(m-1)*(m-2)*n/6;ans%=p; 32 int tmp=0; 33 for(int i=1;i<=min(n-1,m-1);i++)tmp=(tmp+m*n%p*phi[i]*((n-1)/i)%p*((m-1)/i))%p; 34 for(int i=1;i<=min(n-1,m-1);i++)tmp=(tmp+phi[i]*i%p*i%p*(((n-1)/i)*(1+(n-1)/i)/2)%p*(((m-1)/i)*(1+(m-1)/i)/2)%p)%p; 35 for(int i=1;i<=min(n-1,m-1);i++)tmp=(tmp-n*phi[i]%p*i%p*((n-1)/i)%p*((((m-1)/i)*(1+(m-1)/i)/2)%p)%p+p)%p; 36 for(int i=1;i<=min(n-1,m-1);i++)tmp=(tmp-m*phi[i]%p*i%p*((m-1)/i)%p*((((n-1)/i)*(1+(n-1)/i)/2)%p)%p+p)%p; 37 tmp-=((1+n-1)*(n-1)/2)%p*((1+m-1)*(m-1)/2)%p;tmp=(tmp+p)%p; 38 ans=(ans+tmp*2)%p; 39 printf("%lld\n",ans); 40 return 0; 41 }
