简单题（确信）
不难发现：
[f(n)=mu^2(n)
]
把这个代进原式，然后开始推式子：
[egin{aligned}
ext{Ans}&=sum_{i=1}^nsum_{j=1}^n(i+j)^kmu^2({
m{gcd}}(i,j)){
m{gcd}}(i,j)\
&=sum_{d=1}^nsum_{i=1}^nsum_{j=1}^n(i+j)^kmu^2(d)d[{
m{gcd}}(i,j)=d]\
&=sum_{d=1}^nd^{k+1}mu^2(d)sum_{i=1}^{lfloorfrac{n}{d}
floor}sum_{j=1}^{lfloorfrac{n}{d}
floor}(i+j)^ksum_{x|i,x|j}mu(x)\
&=sum_{d=1}^nd^{k+1}mu^2(d)sum_{x=1}^{lfloorfrac{n}{d}
floor}mu(x)x^ksum_{i=1}^{lfloorfrac{n}{dx}
floor}sum_{j=1}^{lfloorfrac{n}{dx}
floor}(i+j)^k
end{aligned}
]
令 (S(n)=sum_{i=1}^nsum_{j=1}^n(i+j)^k) ，代入：
[egin{aligned}
ext{Ans}&=sum_{d=1}^nd^{k+1}mu^2(d)sum_{x=1}^{lfloorfrac{n}{d}
floor}mu(x)x^kS(lfloorfrac{n}{dx}
floor)\
&=sum_{T=1}^nS(lfloorfrac{n}{T}
floor)sum_{d|T}d^{k+1}mu^2(d)mu(frac{T}{d})(frac{T}{d})^k\
&=sum_{T=1}^nS(lfloorfrac{n}{T}
floor)sum_{d|T}dmu^2(d)mu(frac{T}{d})T^k
end{aligned}
]
都是套路。
令 (f(n)=sum_{d|n}dmu^2(d)mu(frac{n}{d})) 。
如果我们能求得 (S(n)) 和 (sum_{i=1}^nf(i)) ，再用一下数论分块就做完了。
先考虑如何求 (S(n)) 。
令 (F(n)=sum_{i=1}^ni^k,G(n)=sum_{i=1}^nF(i)) ，则有：
[S(n)=sum_{i=n+1}^{2n}F(i)-sum_{i=1}^nF(i)=G(2n)-2 imes G(n)
]
证明：
考虑枚举 (x=i+j) 。将 (S(n)=sum_{i=1}^nsum_{j=1}^n(i+j)^k) 和式展开，找找规律，发现：
- 若 (x leq n+1) ，则有 (x-1) 个数对 ((i,j)(1 leq i,j leq n)) 能表示出 (x) ，则此时 (x) 对 (S(n)) 的贡献是 (x-1) 。
- 若 (n+1 <x leq 2n) ，设 (x=n+i (2 leq i leq n)) 则有 (n-i+1) 个数对 ((i,j)(1 leq i,j leq n)) 能表示出 (x) ，则此时 (x) 对 (S(n)) 的贡献是 (n-i+1) 。
则有：
[egin{aligned}
S(n)&=sum_{i=2}^{n+1}(i-1)*i^k+sum_{i=2}^{n}(n-i+1)(n+i)^k\
&=2^k+2 imes 3^k+cdots+n imes (n+1)^k+(n-1) imes(n+2)^k+(n-2) imes(n+3)^k+cdots+(2n)^k
end{aligned}
]
我们再把要证明的式子展开：
[egin{aligned}
S(n)&=sum_{i=n+1}^{2n}sum_{j=1}^ij^k-sum_{i=1}^nsum_{j=1}^ij^k\
&=n imes1^k+n imes 2^k+n imes3^k+cdots+n imes(n+1)^k+(n-1) imes(n+2)^k+cdots+(2n)^k-\
&(n imes 1^k+(n-1) imes 2^k+(n-2) imes 3^k+cdots+n^k)\
&=2^k+2 imes 3^k+cdots+n imes (n+1)^k+(n-1) imes(n+2)^k+(n-2) imes(n+3)^k+cdots+(2n)^k
end{aligned}
]
所以 (S(n)=sum_{i=n+1}^{2n}F(i)-sum_{i=1}^nF(i)) 。
然后发现这东西能用数学归纳法证明，上面这东西太暴力了。
当 (n=1) 时， (S(1)=G(2)-2G(1)) 显然成立。
假设 (S(n)=G(2n)-2 imes G(n)) 成立，证明 (S(n+1)=G(2n+2)-2 imes G(n+1)) 成立：
[egin{aligned}
S(n+1)&=S(n)+2 imessum_{i=1}^n(i+n+1)^k+(2n+2)^k\
&=S(n)+2 imes F(2n+1)-2 imes F(n+1)+F(2n+2)-F(2n+1)\
&=S(n)+F(2n+2)+F(2n+1)-2 imes F(n+1)\
&=G(2n)-2 imes G(n)+F(2n+2)+F(2n+1)-2 imes F(n+1)\
&=G(2n)+F(2n+1)+F(2n+2)-2 imes G(n)-2 imes F(n+1)\
&=G(2n+2)-2 imes G(n+1)
end{aligned}
]
证毕。
线性筛出 ({ m{id}}_k) 再做两次前缀和，就能做到 (O(1)) 询问 (S(n)) 。
[f(n)=sum_{d|n}dmu^2(d)mu(frac{n}{d})
]
这东西就是几个积性函数卷起来，所以 (f) 也是积性函数，可以线性筛。
对于质数 (p) ，有 (f(p)=mu(p)+pmu^2(p)=p-1) 。
对于 (p^k) ，对 (k) 进行分类讨论：
- 若 (k=2) ，有 (f(p^2)=mu(p^2)+pmu^2(p)mu(p)+p^2mu^2(p^2)=-p) 。
- 若 (k>2) ，对于任意 (d) ，都有 (d) 或 (frac{n}{d}) 有平方因子，所以 (f(p^k)=0) 。
然后就可以线性筛筛 (f) 了：
inline void sieve()
{
f[1]=1;
for(int i=2;i<maxn;++i)
{
if(!flag[i]) prime[++cnt]=i,f[i]=i-1;
for(int j=1;j<=cnt&&i*prime[j]<maxn;++j)
{
flag[i*prime[j]]=true;
if(i%prime[j]) f[i*prime[j]]=f[i]*f[prime[j]]%mod;// 积性函数
else
{
if((i/prime[j])%prime[j])// 互质则用积性函数的性质
f[i*prime[j]]=(mod-prime[j])*f[i/prime[j]]%mod;
// 不互质则说明prime[j]的次数大于2
break;
}
}
}
}
( ext{Code}:)
#include <iostream>
#include <cstring>
#include <cstdio>
#include <algorithm>
#include <cmath>
#define maxn 10000005
#define Rint register int
#define INF 0x3f3f3f3f
using namespace std;
typedef long long lxl;
const lxl mod=998244353;
template <typename T>
inline T read()
{
T x=0,f=1;char ch=getchar();
while(ch<'0'||ch>'9') {if(ch=='-') f=-1;ch=getchar();}
while(ch>='0'&&ch<='9') {x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
return x*f;
}
inline lxl fmi(lxl a,lxl b)
{
lxl ans=1;
a%=mod;
while(b>0)
{
if(b&1) ans=(ans*a)%mod;
a=(a*a)%mod;
b>>=1;
}
return ans;
}
int n;
lxl k;
int prime[maxn],cnt;
bool flag[maxn];
lxl f[maxn],F[maxn];
inline void sieve()
{
f[1]=F[1]=1;
for(int i=2;i<maxn;++i)
{
if(!flag[i]) prime[++cnt]=i,f[i]=i-1,F[i]=fmi(i,k);
for(int j=1;j<=cnt&&i*prime[j]<maxn;++j)
{
flag[i*prime[j]]=true;
F[i*prime[j]]=F[i]*F[prime[j]]%mod;
if(i%prime[j]) f[i*prime[j]]=f[i]*f[prime[j]]%mod;
else
{
if((i/prime[j])%prime[j])
f[i*prime[j]]=(mod-prime[j])*f[i/prime[j]]%mod;
break;
}
}
}
for(int i=1;i<maxn;++i) f[i]=(f[i-1]+f[i]*F[i]%mod)%mod,F[i]=(F[i]+F[i-1])%mod;
for(int i=1;i<maxn;++i) F[i]=(F[i]+F[i-1])%mod;
}
inline lxl S(int n)
{
return (F[n<<1]-2*F[n]%mod+mod)%mod;
}
inline lxl calcu(lxl n)
{
lxl res=0;
for(lxl l=1,r=0;l<=n;l=r+1)
{
r=n/(n/l);
res=(res+S(n/l)*(f[r]-f[l-1]+mod)%mod)%mod;
}
return res;
}
int main()
{
// freopen("P6156.in","r",stdin);
n=read<int >(),k=read<lxl >();
sieve();
printf("%lld
",calcu(n));
return 0;
}