(color{#0066ff}{ 题目描述 })
给定n,m,k,计算
(sum_{i=1}^n sum_{j=1}^m mathrm{gcd}(i,j)^k)
对1000000007取模的结果
(color{#0066ff}{输入格式})
多组数据。 第一行是两个数T,K; 之后的T行，每行两个整数n，m；
(color{#0066ff}{输出格式})
K行，每行一个结果
(color{#0066ff}{输入样例})
1 2
3 3
(color{#0066ff}{输出样例})
20
(color{#0066ff}{数据范围与提示})
T<=2000,1<=N,M,K<=5000000
(color{#0066ff}{ 题解 })
就是要求
[sum_{i=1}^n sum_{j=1}^m gcd(i,j)^k
]
枚举gcd
[sum_{d=1}^{min(n,m)}sum_{i=1}^n sum_{j=1}^m [gcd(i,j)==d]d^k
]
把(d^k)提出来，d再除上去，就是一个基本模型了
[sum_{d=1}^{min(n,m)}d^ksum_{i=1}^{lfloorfrac{n}{d}
floor} sum_{j=1}^{lfloorfrac{m}{d}
floor} [gcd(i,j)==1]
]
[sum_{d=1}^{min(n,m)}d^ksum_{i=1}^{lfloorfrac{n}{d}
floor} sum_{j=1}^{lfloorfrac{m}{d}
floor} sum_{k|gcd(i,j)} mu(k)
]
[sum_{d=1}^{min(n,m)}d^ksum_{k=1}^{min(lfloorfrac{n}{d}
floor,lfloorfrac{m}{d}
floor)}mu(k)sum_{i=1}^{lfloorfrac{n}{kd}
floor} sum_{j=1}^{lfloorfrac{m}{kd}
floor} 1
]
后面好像空了。。。
[sum_{d=1}^{min(n,m)}d^ksum_{k=1}^{min(lfloorfrac{n}{d}
floor,lfloorfrac{m}{d}
floor)}mu(k) * lfloorfrac{n}{kd}
floor * lfloorfrac{m}{kd}
floor
]
来一发kd换q
[sum_{q=1}^{min(n,m)} lfloorfrac{n}{q}
floor * lfloorfrac{m}{q}
floor sum_{d|q}mu(frac q d)*d^k
]
“额，这怎么处理”
“暴力了解一下”
线性筛出(mu) 然后(O(nlogn))求出(d^k)
之后枚举倍数(O(nlogn)把后面的)sum$搞出来，数列分块就行了
#include<bits/stdc++.h>
#define LL long long
LL in() {
char ch; LL x = 0, f = 1;
while(!isdigit(ch = getchar()))(ch == '-') && (f = -f);
for(x = ch ^ 48; isdigit(ch = getchar()); x = (x << 1) + (x << 3) + (ch ^ 48));
return x * f;
}
const int mod = 1e9 + 7;
const int maxn = 5e6 + 100;
int k;
int mu[maxn], pri[maxn], tot, mi[maxn];
bool vis[maxn];
LL h[maxn];
LL ksm(LL x, LL y) {
LL re = 1LL;
while(y) {
if(y & 1) re = re * x % mod;
x = x * x % mod;
y >>= 1;
}
return re;
}
void predoit() {
mu[1] = 1;
for(int i = 2; i < maxn; i++) {
if(!vis[i]) pri[++tot] = i, mu[i] = -1;
for(int j = 1; j <= tot && (LL)i * pri[j] < maxn; j++) {
vis[i * pri[j]] = true;
if(i % pri[j] == 0) break;
else mu[i * pri[j]] = -mu[i];
}
}
for(int i = 1; i < maxn; i++) mi[i] = ksm(i, k);
for(int i = 1; i < maxn; i++)
for(int j = i; j < maxn; j += i)
(h[j] += (1LL * mu[j / i] * mi[i] % mod)) %= mod;
for(int i = 2; i < maxn; i++) (h[i] += h[i - 1]) %= mod;
}
LL work(LL n, LL m) {
LL ans = 0;
for(LL l = 1, r; l <= std::min(n, m); l = r + 1) {
r = std::min(n / (n / l), m / (m / l));
LL tot1 = (n / l) * (m / l) % mod;
tot1 = (tot1 * ((h[r] - h[l - 1]) % mod + mod) % mod) % mod;
(ans += tot1) %= mod;
}
return ans;
}
int main() {
int T;
for(T = in(), k = in(), predoit(); T --> 0;)
printf("%lld
", work(in(), in()));
return 0;
}