简要题意:
给定一个固定的 k k k, T T T 组询问求:
∑ i = 1 n ∑ j = 1 m gcd ( i , j ) k \sum_{i=1}^n \sum_{j=1}^m \gcd(i,j)^k i=1∑nj=1∑mgcd(i,j)k
其中 n , m , k ≤ 5 × 1 0 6 n,m,k \leq 5 \times 10^6 n,m,k≤5×106, T ≤ 2 × 1 0 3 T \leq 2 \times 10^3 T≤2×103.
没什么好说的,没有部分分。根据 莫比乌斯反演性质 开始推式子!
∑ i = 1 n ∑ j = 1 m gcd ( i , j ) k \sum_{i=1}^n \sum_{j=1}^m \gcd(i,j)^k i=1∑nj=1∑mgcd(i,j)k
= ∑ i = 1 min ( n , m ) d k ∑ i = 1 ⌊ n d ⌋ ∑ j = 1 ⌊ n d ⌋ [ gcd ( i , j ) = = 1 ] = \sum_{i=1}^{\min(n,m)} d^k \sum_{i=1}^{\lfloor \frac{n}{d} \rfloor} \sum_{j=1}^{\lfloor \frac{n}{d} \rfloor} [\gcd(i,j)==1] =i=1∑min(n,m)dki=1∑⌊dn⌋j=1∑⌊dn⌋[gcd(i,j)==1]
= ∑ i = 1 min ( n , m ) d k ∑ i = 1 ⌊ n d ⌋ ∑ j = 1 ⌊ n d ⌋ ∑ k ∣ gcd ( i , j ) μ k = \sum_{i=1}^{\min(n,m)} d^k \sum_{i=1}^{\lfloor \frac{n}{d} \rfloor} \sum_{j=1}^{\lfloor \frac{n}{d} \rfloor} \sum_{k|\gcd(i,j)} \mu_k =i=1∑min(n,m)dki=1∑⌊dn⌋j=1∑⌊dn⌋k∣gcd(i,j)∑μk
= ∑ i = 1 min ( n , m ) d k ∑ k = 1 ⌊ min ( n , m ) d ⌋ ⌊ n k d ⌋ ⌊ m k d ⌋ = \sum_{i=1}^{\min(n,m)} d^k \sum_{k=1}^{\lfloor \frac{\min(n,m)}{d} \rfloor} \lfloor \frac{n}{kd} \rfloor \lfloor \frac{m}{kd} \rfloor =i=1∑min(n,m)dkk=1∑⌊dmin(n,m)⌋⌊kdn⌋⌊kdm⌋
= ∑ i = 1 T ⌊ n T ⌋ ⌊ m T ⌋ ∑ d ∣ T d k × μ T d = \sum_{i=1}^T \lfloor \frac{n}{T} \rfloor \lfloor \frac{m}{T} \rfloor \sum_{d|T} d^k \times \mu_{\frac{T}{d}} =i=1∑T⌊Tn⌋⌊Tm⌋d∣T∑dk×μdT
显然 O ( T n ) O(T \sqrt{n}) O(Tn) 的时间足以通过,前面的部分整除分块足够。但是对于后面这一块较难解决,则设:
g x = ∑ d ∣ x d k × μ x d g_x= \sum_{d|x} d^k \times \mu_{\frac{x}{d}} gx=d∣x∑dk×μdx
通过观察可知这是个 积性函数。那么显然,对于:
x = ∏ i = 1 q p i a i x = \prod_{i=1}^q p_i^{a_i} x=i=1∏qpiai
存在:
g x = ∏ i = 1 q g ( p i a i ) g_x = \prod_{i=1}^q g(p_i^{a_i}) gx=i=1∏qg(piai)
= ∏ i = 1 q ( p i k × ( a i − 1 ) × μ p i + p i k × a i × μ 1 ) = \prod_{i=1}^q (p_i^{k \times (a_i-1)} \times \mu_{p_i} + p_i^{k \times a_i} \times \mu_1) =i=1∏q(pik×(ai−1)×μpi+pik×ai×μ1)
= ∏ i = 1 q p i k × ( a i − 1 ) × ( p i k − 1 ) = \prod_{i=1}^q p_i^{k \times (a_i-1)} \times (p_i^k-1) =i=1∏qpik×(ai−1)×(pik−1)
首先,上面对于 g g g 的几步需要解释。那个关于 x x x 的式子实际上是 x x x 质因数分解以后的结果, p i ∈ prime p_i \in \text{prime} pi∈prime.
然后第一步是线性筛积性函数的套路。
第二步是一个变形,考虑展开 g ( p i a i ) g(p_i^{a_i}) g(piai) 的过程.
第三步就是个乘法分配律,算出 μ \mu μ 而已。
好了,看看这个式子,每个部分都可以筛。结束。
时间复杂度: O ( n + T n + p log k ) O(n + T \sqrt{n} + p \log k) O(n+Tn+plogk). p p p 的值之后解释。
实际得分: 100 p t s 100pts 100pts.
#pragma GCC optimize(2)
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=5e6+1;
const ll MOD=1e9+7;
inline int read(){char ch=getchar(); int f=1; while(ch<'0' || ch>'9') {if(ch=='-') f=-f; ch=getchar();}
int x=0; while(ch>='0' && ch<='9') x=(x<<3)+(x<<1)+ch-'0',ch=getchar(); return x*f;}
inline void write(ll x) {
if(x<0) {putchar('-');write(-x);return;}
if(x<10) {putchar(char(x%10+'0'));return;}
write(x/10);putchar(char(x%10+'0'));
}
int prime[N],mu[N];
ll g[N],f[N]; bool h[N];
int T,k,cnt=0; ll ans=0;
inline ll pw(ll x,ll y) {
ll ans=1; while(y) {
if(y&1) ans=ans*x%MOD;
x=x*x%MOD; y>>=1;
} return ans;
}
inline void Euler(int n) {
f[1]=1; for(int i=2;i<=n;i++) {
if(!h[i]) prime[++cnt]=i,g[cnt]=pw(i,k),f[i]=(g[cnt]-1+MOD)%MOD;
//f 表示 分解质因数后后面的一堆东西
for(int j=1;j<=cnt && i*prime[j]<=n;j++) {
h[i*prime[j]]=1;
if(i%prime[j]==0) {f[i*prime[j]]=1ll*f[i]*g[j]%MOD;break;}
f[i*prime[j]]=1ll*f[i]*f[prime[j]]%MOD; //线性筛模板
}
}
for(int i=1;i<=n;i++) f[i]=(f[i]+f[i-1])%MOD; //前缀和
}
int main() {
T=read(); k=read(); Euler(N-1);
while(T--) {
int n=read(),m=read();
ll ans=0; int r;
for(int i=1;i<=min(n,m);i=r+1) { //整除分块
r=min(n/(n/i),m/(m/i));
ans=(ans+1ll*(f[r]-f[i-1]+MOD)*(n/i)%MOD*(m/i)%MOD)%MOD;
} write(ans); putchar('\n');
}
return 0;
}
显然,你发现对于素数我们需要暴力计算 p k p^k pk. 那么,首先 5 × 1 0 6 5 \times 10^6 5×106 以内的素数有多少个呢?代码测试:
得到结果约 3.4 × 1 0 5 3.4 \times 10^5 3.4×105,也即 p ≤ 3.4 × 1 0 5 p \leq 3.4 \times 10^5 p≤3.4×105. 而 log k = log ( 5 × 1 0 6 ) \log k = \log (5 \times 10^6) logk=log(5×106),大概是 22 22 22,不会出问题,显然 3.4 × 1 0 5 × 22 3.4 \times 10^5 \times 22 3.4×105×22 是 1 0 7 10^7 107 级的,这也是程序较慢的原因(不过不影响通过,照样 2.55 s 2.55s 2.55s)