51NOD1847：奇怪的数学题

传送门

Sol

设 $f(d)$ 表示 $d$ 所有约数中第二大的，$lo{w}_d$ 表示 $d$ 的最小质因子
$$f(d)=\frac{d}{lo{w}_d}$$
那么
$$\sum _{i=1}^n\sum _{j=1}^{n}sgc{d}^k(i,j)$$
$$=\sum _{i=1}^n\sum _{j=1}^n{f}^k(gcd(i,j))$$
$$=\sum _{d=1}^n{f}^k(d)\sum _{i=1}^n\sum _{j=1}^n[gcd(i,j)=d]$$
$$=\sum _{d=1}^n{f}^k(d)\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{d}\rfloor }[gcd(i,j)=1]$$
$$=\sum _{d=1}^n{f}^k(d)(2\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\phi (i)-1)$$
可以数论分块，后面的 $\phi$ 直接杜教筛
考虑计算
$$\sum _{d=2}^n{f}^k(d)=\sum _{d=2}^n{\frac{d}{lo{w}_d}}^k$$
设 $p_j$ 表示第 $j$ 个质数
$g(x,i)$ 表示 $2$ 到 $x$ 之间最小质因子大于等于 $p_i$ 的或者质数的 $f$ 的 $k$ 次方和
$g^{\prime }(x,i)$ 表示 $2$ 到 $x$ 之间最小质因子大于等于 $p_i$ 的或者质数的 $k$ 次方和
设 $s(i)$ 表示小于等于 $p_i$ 的质数的 $k$ 次方和
那么就是要求 $g(n,1)$
$s$ 直接 $min25$ 筛
$$g(x,i)=g(x,i+1)+\sum _{e=1}^{p_i^{e+1}\le x}{p}_i^{k(e-1)}(g^{\prime }(\lfloor \frac{x}{p_i^e}\rfloor ,i+1)-s(i)+p_i^k)$$
$$g^{\prime }(x,i)=g^{\prime }(x,i+1)+\sum _{e=1}^{p_i^{e+1}\le x}{p}_i^{ke}(g^{\prime }(\lfloor \frac{x}{p_i^e}\rfloor ,i+1)-s(i)+p_i^k)$$
这里要用到自然幂数和求和，用第二类斯特林数就好了

# include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef unsigned int uint;

inline uint Pow(uint x, int y) {
	register uint ret = 1;
	for (; y; y >>= 1, x = x * x) if (y & 1) ret = ret * x;
	return ret;
}

const int maxn(1e6 + 5);

int pr[maxn], tot, k, d, id1[maxn], id2[maxn], cnt, phi[maxn];
uint n, s[55][55], f[maxn], g[maxn], val[maxn], sk[maxn];
uint s1[maxn], s2[maxn], sphi[maxn], sp[maxn];
bitset <maxn> ispr;

inline void Sieve(int mx) {
	register int i, j;
	for (phi[1] = 1, ispr[1] = 1, i = 2; i <= mx; ++i) {
		if (!ispr[i]) pr[++tot] = i, sk[tot] = sk[tot - 1] + Pow(i, k), phi[i] = i - 1;
		for (j = 1; j <= tot && pr[j] * i <= mx; ++j) {
			ispr[i * pr[j]] = 1;
			if (i % pr[j]) phi[i * pr[j]] = phi[i] * (pr[j] - 1);
			else {
				phi[i * pr[j]] = phi[i] * pr[j];
				break;
			}
		}
	}
	for (i = 1; i <= mx; ++i) sphi[i] = sphi[i - 1] + phi[i];
}

# define ID(x) ((x) <= d ? id1[x] : id2[n / (x)])

inline uint Sum(uint x) {
	register uint i, j, v = 0, t, r, tmp = k <= x ? k : x;
	for (i = 1; i <= tmp; ++i) {
		t = i + 1, r = s[k][i];
		for (j = x - i + 1; j <= x + 1; ++j)
			if (t > 1 && j % t == 0) r *= j / t, t = 1;
			else r = r * j;
		v += r;
	}
	return v;
}

uint Sumphi(uint x) {
    if (x <= d) return sphi[x];
    if (sp[ID(x)]) return sp[ID(x)];
    register uint ans = (x & 1) ? ((x + 1) >> 1) * x : (x >> 1) * (x + 1), i, j;
    for (i = 2; i <= x; i = j + 1) j = x / (x / i), ans -= Sumphi(x / i) * (j - i + 1);
    return sp[ID(x)] = ans;
}

int main() {
	register uint i, j, e, ans = 0, lst = 0, cur, r, v, tmp, now;
	scanf("%u%d", &n, &k), Sieve(d = sqrt(n));
	for (i = 1; i <= 50; ++i)
		for (s[i][1] = 1, j = 2; j <= i; ++j)
			s[i][j] = s[i - 1][j] * j + s[i - 1][j - 1];
	for (i = 1; i <= n; i = j + 1) {
		val[++cnt] = n / i, j = n / (n / i);
		val[cnt] <= d ? id1[val[cnt]] = cnt : id2[n / val[cnt]] = cnt;
		f[cnt] = val[cnt] - 1, g[cnt] = Sum(val[cnt]) - 1;
	}
	for (i = 1; i <= tot && pr[i] * pr[i] <= n; ++i)
		for (j = 1; j <= cnt && pr[i] * pr[i] <= val[j]; ++j) {
			f[j] -= f[ID(val[j] / pr[i])] - i + 1;
			g[j] -= (sk[i] - sk[i - 1]) * (g[ID(val[j] / pr[i])] - sk[i - 1]);
		}
	for (i = 1; i <= cnt; ++i) s1[i] = f[i], s2[i] = g[i];
	for (r = 1; r <= tot && pr[r] * pr[r] <= n; ++r);
	for (i = r - 1; i; --i)
		for (tmp = sk[i] - sk[i - 1], j = 1; j <= cnt && pr[i] * pr[i] <= val[j]; ++j)
			for (cur = e = 1, v = pr[i]; v <= val[j] / pr[i]; ++e, v *= pr[i], cur *= tmp)
				now = (s2[ID(val[j] / v)] - sk[i] + tmp) * cur, s1[j] += now, s2[j] += tmp * now;
	for (i = 1; i <= n; i = j + 1) {
		j = n / (n / i), cur = s1[ID(j)];
		ans += (cur - lst) * (Sumphi(n / i) * 2 - 1), lst = cur;
	}
	printf("%u\n", ans);
	return 0;
}