# bzoj4407 于神之怒加强版

T行，每行一个结果

1 2
3 3

20

### Hint

T<=2000,1<=N,M,K<=5000000

$Ans=\sum_{i=1}^{n}\sum_{j=1}^{m}\gcd(i,j)^{k}$

$Ans=\sum_{d=1}^{min(n,m)}d^{k}\sum_{i=1}^{n}\sum_{j=1}^{m}[\gcd(i,j)==d]$

 1 //It is made by wfj_2048~
2 #include <algorithm>
3 #include <iostream>
4 #include <complex>
5 #include <cstring>
6 #include <cstdlib>
7 #include <cstdio>
8 #include <vector>
9 #include <cmath>
10 #include <queue>
11 #include <stack>
12 #include <map>
13 #include <set>
14 #define inf (1<<30)
15 #define N (5000010)
16 #define il inline
17 #define RG register
18 #define ll long long
19 #define rhl (1000000007)
20 #define File(s) freopen(s".in","r",stdin),freopen(s".out","w",stdout)
21
22 using namespace std;
23
24 int vis[N],prime[N],T,k,n,m,cnt,pos;
25 ll bin[N],mu[N],f[N],ans;
26
27 il int gi(){
28     RG int x=0,q=1; RG char ch=getchar(); while ((ch<'0' || ch>'9') && ch!='-') ch=getchar();
29     if (ch=='-') q=-1,ch=getchar(); while (ch>='0' && ch<='9') x=x*10+ch-48,ch=getchar(); return q*x;
30 }
31
32 il ll qpow(RG ll a,RG ll b){
33     RG ll ans=1;
34     while (b){
35     if (b&1) ans=ans*a; a=a*a,b>>=1;
36     if (ans>=rhl) ans%=rhl;
37     if (a>=rhl) a%=rhl;
38     }
39     return ans;
40 }
41
42 il void sieve(){
43     vis[1]=mu[1]=f[1]=bin[1]=1;
44     for (RG int i=2;i<N;++i){
45     if (!vis[i]) vis[i]=1,mu[i]=rhl-1,bin[i]=qpow(i,k),f[i]=(bin[i]+rhl-1)%rhl,prime[++cnt]=i;
46     for (RG int j=1,k=i*prime[j];j<=cnt && k<N;++j,k=i*prime[j]){
47         vis[k]=1;
48         if (i%prime[j]) mu[k]=rhl-mu[i],f[k]=f[i]*f[prime[j]]%rhl;
49         else{ f[k]=f[i]*bin[prime[j]]%rhl; break; }
50     }
51     }
52     for (RG int i=2;i<N;++i){
53     mu[i]+=mu[i-1],f[i]+=f[i-1];
54     if (mu[i]>=rhl) mu[i]-=rhl;
55     if (f[i]>=rhl) f[i]-=rhl;
56     }
57     return;
58 }
59
60 il void work(){
61     n=gi(),m=gi(),ans=0; if (n>m) swap(n,m);
62     for (RG int i=1;i<=n;i=pos+1){
63     pos=min(n/(n/i),m/(m/i));
64     ans+=(ll)(n/i)*(ll)(m/i)%rhl*(f[pos]-f[i-1]+rhl)%rhl;
65     if (ans>=rhl) ans-=rhl;
66     }
67     printf("%lld\n",ans); return;
68 }
69
70 int main(){
71     File("4407");
72     T=gi(),k=gi(),sieve();
73     while (T--) work();
74     return 0;
75 }

posted @ 2017-03-16 11:36  wfj_2048  阅读(146)  评论(0编辑  收藏  举报