计算过程:
令(S(n)=sum^n_{a_1=1}sum^n_{a_2=1}cdotssum^n_{a_n=1}(prod^x_{j=1}a_j^k)=sum^n_{a_1=1}{a_1}^ksum^n_{a_2=1}{a_2}^kcdots sum^n_{a_n=1}{a_n}^k=({sum^n_{i=1}i^k})^x)
(sum^n_{a_1=1}sum^n_{a_2=1}cdotssum^n_{a_n=1}(prod^x_{j=1}a_j^k)f(gcd(a_1,a_2,cdots,a_n))gcd(a_1,a_2,cdots,a_n))
(=sum^n_{d=1}df(d)S(n)[gcd(a_1,a_2,cdots,a_n)=d])
(=sum^n_{d=1}df(d)S(lfloor{frac{n}{d}} floor)[gcd(a_1,a_2,cdots,a_n)=1]d^{kx})
(=sum^n_{d=1}d^{kx+1}f(d)S(lfloor{frac{n}{d}} floor)[gcd(a_1,a_2,cdots,a_n)=1])
(=sum^n_{d=1}d^{kx+1}f(d)S(lfloor{frac{n}{d}} floor)sum_{t|gcd(a_1,a_2,cdots,a_n)}mu(t))
枚举(t)
(=sum^n_{d=1}d^{kx+1}f(d)sum^{lfloor{frac{n}{d}} floor}_{t=1}mu(t)t^{kx}(sum^{lfloor{frac{n}{dt}} floor}_{i=1}i^k)^x)
(=sum^n_{d=1}d^{kx+1}f(d)sum^{lfloor{frac{n}{d}} floor}_{t=1}mu(t)t^{kx}S(lfloor{frac{n}{dt}} floor))
令(T=dt)
(sum^n_{T=1}S(lfloor{frac{n}{T}} floor)sum_{d|T}d^{kx+1}(frac{T}{d})^{kx}mu(frac{T}{d})f(d))
(sum^n_{T=1}S(lfloor{frac{n}{T}} floor)sum_{d|T}d^{kx+1}(frac{T}{d})^{kx}mu(frac{T}{d})f(d))
(sum^n_{T=1}S(lfloor{frac{n}{T}} floor)sum_{d|T}dT^{kx}mu(frac{T}{d})f(d))
具体实现
#include <bits/stdc++.h>
#include <ext/pb_ds/assoc_container.hpp>
#include <ext/pb_ds/hash_policy.hpp>
#include <ext/pb_ds/tree_policy.hpp>
#include <ext/pb_ds/trie_policy.hpp>
using namespace __gnu_pbds;
using namespace std;
/* freopen("k.in", "r", stdin);
freopen("k.out", "w", stdout); */
// clock_t c1 = clock();
// std::cerr << "Time:" << clock() - c1 <<"ms" << std::endl;
//#pragma comment(linker, "/STACK:1024000000,1024000000")
#define de(a) cout << #a << " = " << a << endl
#define rep(i, a, n) for (int i = a; i <= n; i++)
#define per(i, a, n) for (int i = n; i >= a; i--)
#define ls ((x) << 1)
#define rs ((x) << 1 + 1)
typedef long long ll;
typedef unsigned long long ull;
typedef pair<int, int> PII;
typedef pair<double, double> PDD;
typedef pair<ll, ll> PLL;
typedef vector<int> VI;
#define inf 0x3f3f3f3f
const ll INF = 0x3f3f3f3f3f3f3f3f;
const ll MAXN = 2e5 + 7;
const ll MAXM = 4e5 + 7;
const ll MOD = 1e9 + 7;
const double eps = 1e-7;
const double pi = acos(-1.0);
ll k, x;
ll quick_pow(ll a, ll b)
{
ll ans = 1;
while (b)
{
if (b & 1)
ans = (1LL * ans * a) % MOD;
a = (1LL * a * a) % MOD;
b >>= 1;
}
return ans;
}
ll mu[MAXN], pri[MAXN], vis[MAXN], tot = 0;
ll sum[MAXN];
void init()
{
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 && pri[j] * i < MAXN; j++)
{
vis[i * pri[j]] = 1;
if (i % pri[j] == 0)
mu[i * pri[j]] = 0;
else
mu[i * pri[j]] = -mu[i];
}
}
for (int i = 1; i < MAXN; i++)
sum[i] = sum[i - 1] + mu[i];
}
ll f[MAXN];
void calf()
{
for (int i = 0; i < MAXN; i++)
f[i] = 1;
for (ll i = 2; i * i < MAXN; i++)
for (ll j = 1; j * i * i < MAXN; j++)
f[j * i * i] = 0;
}
ll s[MAXN] = {0};
void cals()
{
ll now = 0;
for (ll i = 1; i < MAXN; i++)
{
(now += quick_pow(i, k)) %= MOD;
s[i] = quick_pow(now, x);
}
}
ll g[MAXN] = {0};
ll sumg[MAXN] = {0};
void calg()
{
for (ll d = 1; d < MAXN; d++)
for (ll t = 1; t * d < MAXN; t++)
{
ll temp = (d * mu[t]) % MOD * f[d] % MOD;
(g[t * d] += temp % MOD) %= MOD;
}
for (int i = 1; i < MAXN; i++)
(sumg[i] = sumg[i - 1] + g[i] * quick_pow(i, k * x)) %= MOD;
}
int main()
{
int t;
scanf("%d%lld%lld", &t, &k, &x);
init(), calf(), cals(), calg();
while (t--)
{
ll n;
scanf("%lld", &n);
ll ans = 0;
for (ll l = 1, r; l <= n; l = r + 1)
{
r = n / (n / l);
ll temp = s[n / l];
(temp *= (sumg[r] - sumg[l - 1])) %= MOD;
(temp += MOD) %= MOD;
(ans += temp) %= MOD;
}
printf("%lld
", (ans % MOD + MOD) % MOD);
}
return 0;
}