以下用等号代替同余
这个式子是$\sum\limits_{j=1}^n(i,j)^{c-d}i^dj^dx_j=b_i$
令$g(n)=\sum\limits_{e|n}\mu\left(\frac ne\right)e^{c-d}$,那么原式变为$\sum\limits_{e|i}g(e)\sum\limits_{\substack{e|j\\1\leq j\leq n}}j^dx_j=\frac{b_i}{i^d}$
令$z_e=\sum\limits_{\substack{e|j\\1\leq j\leq n}}j^dx_j$,先处理原式
$\sum\limits_{e|i}g(e)z_e=\frac{b_i}{i^d}\Leftrightarrow g(i)z_i=\sum\limits_{e|i}\mu\left(\frac ie\right)\frac{b_e}{e^d}$
$g(i)$和上式右边都可以$O(n\log n)$预处理,此时如果存在$i$使$g(i)=0$且右边$\ne0$那么无解,两个都是$0$就多解,随便给$z_i$赋值即可
最后求$x_i$:$x_ii^d=\sum\limits_{\substack{i|j\\1\leq j\leq n}}\mu\left(\frac ji\right)z_j$
正如vfk的题解所说:三个莫比乌斯反演掷地有声==
写快速幂时要注意负指数的处理...
#include<stdio.h>
#include<string.h>
char s[4000000];
int ns;
#define NUM(x) ('0'<=x&&x<='9')
inline int rd(){
while(!NUM(s[ns]))ns++;
int x=0;
while(NUM(s[ns]))x=x*10+s[ns++]-'0';
return x;
}
typedef long long ll;
const int mod=998244353,T=100000;
int mul(int a,int b){return(ll)a*b%mod;}
void inc(int&a,int b){(a+=b)%=mod;}
int pow(int a,int b){
int s=1;
while(b<0)b+=mod-1;
while(b){
if(b&1)s=mul(s,a);
a=mul(a,a);
b>>=1;
}
return s;
}
int pr[T+10],mu[T+10],nd[T+10],d;
bool np[T+10];
void sieve(){
int i,j,M=0;
mu[1]=1;
nd[1]=1;
for(i=2;i<=T;i++){
if(!np[i]){
pr[++M]=i;
mu[i]=-1;
nd[i]=pow(i,-d);
}
for(j=1;j<=M&&i*pr[j]<=T;j++){
np[i*pr[j]]=1;
nd[i*pr[j]]=mul(nd[i],nd[pr[j]]);
if(i%pr[j]==0)break;
mu[i*pr[j]]=-mu[i];
}
}
}
int b[T+10],g[T+10],h[T+10],n,c;
int main(){
fread(s,1,4000000,stdin);
int q,i,j,t;
n=rd();
c=rd();
d=rd();
q=rd();
sieve();
for(i=1;i<=n;i++){
t=pow(i,c-d);
for(j=1;i*j<=n;j++)inc(g[i*j],mu[j]*t);
}
for(i=1;i<=n;i++)g[i]=pow(g[i],mod-2);
while(q--){
for(i=1;i<=n;i++)b[i]=mul(rd(),nd[i]);
memset(h,0,sizeof(h));
for(i=1;i<=n;i++){
t=b[i];
for(j=1;i*j<=n;j++){
if(mu[j])inc(h[i*j],mu[j]*t);
}
}
for(i=1;i<=n;i++){
if(g[i]==0&&h[i]!=0){
i=-1;
break;
}
h[i]=mul(h[i],g[i]);
}
if(i==-1){
puts("-1");
continue;
}
for(i=1;i<=n;i++){
t=0;
for(j=1;i*j<=n;j++)inc(t,mu[j]*h[i*j]);
t=mul(nd[i],t);
inc(t,mod);
printf("%d ",t);
}
putchar('\n');
}
}