●BZOJ 4407 于神之怒加强版

http://www.lydsy.com/JudgeOnline/problem.php?id=4407

\begin{align*}ANS&=\sum_{i=1}^{n}\sum_{j=1}^{m}gcd(i,j)^k(接下来先枚举gcd的值g，然后求出有多少对(i,j)的gcd=g) \\&=\sum_{g=1}^{min(n,m)}g^k\sum_{d=1}{\frac{n}{g}\mu(d)\lfloor \frac{n}{gd} \rfloor \lfloor \frac{m}{gd} \rfloor}\\ &=\sum_{D=gd=1}^{min(n,m)} \lfloor \frac{n}{D} \rfloor \lfloor \frac{m}{D} \rfloor \sum_{g|D} g^k \mu(\frac{D}{g})\;(令f(D)=\sum_{g|D} g^k \mu(\frac{D}{g}))\\ &=\sum_{D=1}^{min(n,m)} \lfloor \frac{n}{D} \rfloor \lfloor \frac{m}{D} \rfloor f(D) \end{align*}

\begin{align*} f(D)&=\sum_{g|D} g^k \mu(\frac{D}{g}) \\&=\sum_{g|D} w(g) \mu(\frac{D}{g}) \end{align*}

#include<cstdio>
#include<cstring>
#include<iostream>
#define MAXN 5000050
using namespace std;
const int mod=1000000007;
int f[MAXN],k;
int Pow(int a,int b){
int ret=1; a=(a%mod+mod)%mod;
while(b){
if(b&1) ret=1ll*ret*a%mod;
b>>=1; a=1ll*a*a%mod;
}
return ret;
}
void Sieve(){
static bool np[MAXN];
static int prime[MAXN],pnt;
f[1]=1;
for(int i=2,tmp,d;i<=5000000;i++){
if(!np[i]) prime[++pnt]=i,f[i]=(Pow(i,k)-1+mod)%mod;
for(int j=1;j<=pnt&&i<=5000000/prime[j];j++){
np[i*prime[j]]=1; tmp=i; d=prime[j];
while(tmp%prime[j]==0) tmp/=prime[j],d*=prime[j];
if(tmp!=1) f[tmp*d]=1ll*f[tmp]*f[d]%mod;
else f[d]=(1ll*Pow(d,k)-Pow(d/prime[j],k)+mod)%mod;
if(i%prime[j]==0) break;
}
}
for(int i=1;i<=5000000;i++) f[i]=(1ll*f[i]+f[i-1])%mod;
}
int main(){
int Case,n,m,mini,ans;
scanf("%d%d",&Case,&k);
Sieve();
while(Case--){
scanf("%d%d",&n,&m);
mini=min(n,m); ans=0;
for(int D=1,last;D<=mini;D=last+1){
last=min(n/(n/D),m/(m/D));
ans=(1ll*ans+1ll*(n/D)*(m/D)%mod*(((1ll*f[last]-f[D-1]+mod)%mod))%mod)%mod;
}
printf("%d\n",ans);
}
return 0;
}


Do not go gentle into that good night.
Rage, rage against the dying of the light.
————Dylan Thomas
posted @ 2018-01-17 12:36  *ZJ  阅读(107)  评论(0编辑  收藏  举报