原文链接https://www.cnblogs.com/zhouzhendong/p/UOJ62.html
题解
太久没更博客了，该拯救我的博客了。
$$
\sum_{1\leq j \leq n} \gcd(i,j) ^{c-d} i^dj^dx_j = b_i\\
A_i = i^d x_i, B_i = \frac{b_i}{i^d}, f(x) = x^{c-d}\\
f(x) = \sum_{d|x} g(d) \\
\begin{eqnarray*}
\sum_{1\leq j \leq n} f(\gcd(i,j)) A_j&=&B_i\\
\sum_{1\leq j \leq n} \sum_{1\leq d \leq n} [d|i]\cdot[d|j]\cdot g(d) A_j&=&B_i\\
\sum_{1\leq d \leq n} [d|i]g(d) \sum_{1\leq j\leq n}[d|j]A_j &=& B_i\\
\sum_{1\leq d\leq n}[d|i]h(d) &=& B_i\\
\sum_{d|i} h(d)&=&B_i\\
(h(x) &=& g(x) \sum_{d|j} A_j)
\end{eqnarray*}
$$
于是只需要 2 次因数反演，1 次倍数反演，三次莫比乌斯反演就好了。
代码
#include <bits/stdc++.h>
using namespace std;
typedef long long LL;
LL read(){
LL x=0,f=1;
char ch=getchar();
while (!isdigit(ch)&&ch!='-')
ch=getchar();
if (ch=='-')
f=0,ch=getchar();
while (isdigit(ch))
x=(x<<1)+(x<<3)+(ch^48),ch=getchar();
return f?x:-x;
}
const int N=100005,mod=998244353;
void Del(int &x,int y){if ((x-=y)<0) x+=mod; }
void Add(int &x,int y){if ((x+=y)>=mod) x-=mod; }
int n,q,c,d;
int b[N],f[N],g[N],h[N],x[N];
int Pow(int x,int y){
if (y<0)
return Pow(x,y+mod-1);
int ans=1;
for (;y;y>>=1,x=(LL)x*x%mod)
if (y&1)
ans=(LL)ans*x%mod;
return ans;
}
void Mobius(int *f,int *g,int n,int flag){
for (int i=1;i<=n;i++)
g[i]=f[i];
if (flag==0)
for (int i=1;i<=n;i++)
for (int j=i<<1;j<=n;j+=i)
Del(g[j],g[i]);
else
for (int i=n;i>=1;i--)
for (int j=i<<1;j<=n;j+=i)
Del(g[i],g[j]);
}
void solve(){
for (int i=1;i<=n;i++)
b[i]=(LL)read()*Pow(i,mod-1-d)%mod;
Mobius(b,h,n,0);
for (int i=1;i<=n;i++)
if (g[i])
h[i]=(LL)h[i]*Pow(g[i],-1)%mod;
else if (h[i])
return (void)(puts("-1"));
Mobius(h,x,n,1);
for (int i=1;i<=n;i++)
x[i]=(LL)Pow(i,mod-1-d)*x[i]%mod;
for (int i=1;i<=n;i++)
printf("%d ",x[i]);
puts("");
}
int main(){
n=read(),c=read(),d=read(),q=read();
for (int i=1;i<=n;i++)
f[i]=Pow(i,c-d);
Mobius(f,g,n,0);
while (q--)
solve();
return 0;
}