题目大意
给定(k,x)。(t)次询问，每次给出一个(n)。请你求：
(gcd(a_1,a_2,ldots ,a_n))是(a_1,a_2,...,a_n)的最大公约数。
函数(f(x))定义为：如果存在一个整数(k) ((k>1))使得(k^2)是(x)的约数，则(f(x)=0)，否则(f(x)=1)。
答案对(10^9+7)取模。
数据范围：(1leq k,xleq 10^9)，(1leq tleq 10^4)，(1leq nleq 2 imes 10^5)。
本题题解
先考虑单组数据的情况。
可以计算每个(gcd)值的贡献。不过精确考虑“(gcd)等于”几有些困难，我们可以先对每个(d)求出“(gcd)是(d)的倍数”的贡献，再容斥回去：
// F[d] 表示gcd是d的倍数,对答案的贡献
// 容斥成 G[d] 表示gcd就是d,对答案的贡献
int ans=0;
for(int d=n;d>=1;--d){
G[d]=F[d];
for(int j=d+d;j<=n;j+=d){
G[d]-=G[j];
}
ans+=G[d]*f[d]*d;
}
考虑钦定了(gcd)是(d)的倍数，如何求出(sum_{a_1=1}^{n}sum_{a_2=1}^{n}ldots sum_{a_x=1}^{n}left (prod_{j=1}^{x}a_j^k
ight ))，也就是代码片段里的F[d]。此时的限制，就是所有(a_i)都要是(d)的倍数。可以先令所有(a_i)都除以(d)，那剩下部分的取值，就只有([1,lfloorfrac{n}{d}
floor])，每个(a_i)都能取这些值，那它们的乘积之和，就等于：((1^k+2^k+dots+lfloorfrac{n}{d}
floor^k)^x)。再把所有的(d)乘回去，则(F[d]=d^{kx}(1^k+2^k+dots+lfloorfrac{n}{d}
floor^k)^x)。求出F[]后，再用上面的代码容斥求出答案。
但这只是单组数据时的做法。
对多组数据，每次独立求答案是不现实的。我们考虑从(n-1)向(n)递推答案，每次只需要计算这两者之间微小的变化。
具体来说，要使(n)对新答案有影响，则(d)一定要是(n)的约数（否则(lfloorfrac{n}{d} floor=lfloorfrac{n-1}{d} floor)，求出来的(F)值显然还是一样的）。那么我们每次枚举(n)的约数，时间复杂度相当于(1dots n)的约数数量之和，约等于(O(nlog n))（可以考虑每个数能作为哪些数的约数，显然是(i,2i,3idots ,lfloorfrac{n}{i} floorcdot i)，所以数量是(frac{n}{i})。(1dots n)的总数就是(sum_{i=1}^{n}frac{n}{i})，也就是调和级数）。
现在，我们虽然能快速求出哪些(F)有变化，也能得到新的(F)值，但如果每次重新做一遍容斥，复杂度还是无法接受的。我们考虑把这个容斥转化一下，改成求出每个(F[d])在容斥里的贡献系数(c[d])。这样不管(F[d])怎么变化，我们只需要往答案里加上变化值( imes c[d])即可。
(c[d])怎么求呢？可以先初始化，令(c[d]=f(d)cdot d)（这里的(f)是题面中定义的函数，千万不要与我们定义的大(F)搞混）。然后：
for(int i=1;i<=n;++i)
c[i]=f[i]*i;
for(int i=1;i<=n;++i){
for(int j=i+i;j<=n;j+=i){
c[j]-=c[i];
}
}
这是因为每个(F)值会在它的约数处被减掉。可以把这个代码对照上面的容斥理解，上面的容斥（求出的ans）就等价于(sum_{i=1}^{n}F[i]cdot c[i])。这就是我们求出这个系数(c)数组的妙处。
于是我们就能(O(1))修改并更新答案了。
总时间复杂度(O(nlog n+t))。
参考代码：
//problem:
#include <bits/stdc++.h>
using namespace std;
#define pb push_back
#define mk make_pair
#define lob lower_bound
#define upb upper_bound
#define fi first
#define se second
#define SZ(x) ((int)(x).size())
typedef unsigned int uint;
typedef long long ll;
typedef unsigned long long ull;
typedef pair<int,int> pii;
template<typename T>inline void ckmax(T& x,T y){x=(y>x?y:x);}
template<typename T>inline void ckmin(T& x,T y){x=(y<x?y:x);}
const int MOD=1e9+7;
inline int mod1(int x){return x<MOD?x:x-MOD;}
inline int mod2(int x){return x<0?x+MOD:x;}
inline void add(int& x,int y){x=mod1(x+y);}
inline void sub(int& x,int y){x=mod2(x-y);}
inline int pow_mod(int x,int i){int y=1;while(i){if(i&1)y=(ll)y*x%MOD;x=(ll)x*x%MOD;i>>=1;}return y;}
const int MAXN=2e5;
int K,X,n;
int pw[MAXN+5],sum[MAXN+5],pwX_of_pwK[MAXN+5],pwX_of_sum[MAXN+5];
int F[MAXN+5],coef[MAXN+5],ans[MAXN+5];
vector<int>factors[MAXN+5];
int p[MAXN+5],cnt_p,f[MAXN+5];
bool v[MAXN+5];
void sieve(int lim){
f[1]=1;
for(int i=2;i<=lim;++i){
if(!v[i]){
p[++cnt_p]=i;
f[i]=1;
}
for(int j=1;j<=cnt_p && i*p[j]<=lim;++j){
v[i*p[j]]=1;
if(i%p[j]==0){
break;
}
if(f[i])
f[i*p[j]]=1;
}
}
}
void init(){
sieve(MAXN);
for(int i=1;i<=MAXN;++i){
pw[i]=pow_mod(i,K);
sum[i]=mod1(pw[i]+sum[i-1]);
pwX_of_pwK[i]=pow_mod(pw[i],X);
pwX_of_sum[i]=pow_mod(sum[i],X);
}
}
int main() {
int T;cin>>T>>K>>X;
init();
/*
//单组数据的简单容斥:
cin>>n;
for(int g=1;g<=n;++g){
F[g]=(ll)pow_mod(sum[n/g],X)*pow_mod(pw[g],X)%MOD;
}
int ans=0;
for(int i=n;i>=1;--i){
G[i]=F[i];
for(int j=i+i;j<=n;j+=i){
sub(G[i],G[j]);
}
ans=((ll)ans + (ll)G[i]*f[i]%MOD*i)%MOD;
}
cout<<ans<<endl;
*/
for(int i=1;i<=MAXN;++i)
coef[i]=f[i]*i;
for(int i=1;i<=MAXN;++i){
for(int j=i+i;j<=MAXN;j+=i){
sub(coef[j],coef[i]);
}
}
ans[1]=F[1]=1;
for(int i=1;i<=MAXN;++i){
for(int j=i;j<=MAXN;j+=i){
factors[j].pb(i);
}
}
for(int i=2;i<=MAXN;++i){
ans[i]=ans[i-1];
for(int j=0;j<SZ(factors[i]);++j){
int g=factors[i][j];
sub(ans[i],(ll)F[g]*coef[g]%MOD);
F[g]=(ll)pwX_of_sum[i/g]*pwX_of_pwK[g]%MOD;
add(ans[i],(ll)F[g]*coef[g]%MOD);
}
}
while(T--){
cin>>n;
cout<<ans[n]<<endl;
}
return 0;
}