【UOJ#62】【UR #5】怎样跑得更快（莫比乌斯反演）
【UOJ#62】【UR #5】怎样跑得更快（莫比乌斯反演）
题面
题解
众所周知，\(lcm(i,j)=\frac{ij}{gcd(i,j)}\)，于是原式就变成了：
\[\sum_{j=1}^n gcd(i,j)^{c-d}i^dj^dx_j\equiv b_i
\]
于是我们就可以写成函数的形式：
\[\sum_{j=1}^n f(gcd(i,j))h(i)h(j)x_j\equiv b_i
\]
然后就开始枚举\(gcd\)。
\[\begin{aligned}
b_i&=\sum_{d=1}^n f(d)\sum_{j=1}^n [gcd(i,j)=d]h(i)h(j)x_j\\
&=\sum_{d=1}^n f(d) \sum_{j=1}^{n}[\frac{gcd(i,j)}{d}=1]h(i)h(j)x_j\\
&=\sum_{d|i} f(d) \sum_{d|j}h(i)h(j)\sum_{k|\frac{gcd(i,j)}{d}}\mu(k)x_j
\end{aligned}\]
条件等价于\(kd|gcd(i,j)\)，令\(T=kd\)。
\[\begin{aligned}
b_i&=h(i)\sum_{d|i}f(d)\sum_{d|j}h(j)\sum_{T|gcd(i,j)}\mu(k)x_j\\
&=h(i)\sum_{T|i}\sum_{T|j}\sum_{d|T}f(d)x_j\mu(\frac{T}{d})h(j)\\
&=h(i)\sum_{T|i}\sum_{T|j}g(j)x_j\sum_{d|T}\mu(\frac{T}{d})f(d)
\end{aligned}\]
后半部分可以提前预处理出来，记做\(fr(T)\)。
继续往下就是：
\[\begin{aligned}
b_i&=h(i)\sum_{T|i}\sum_{T|j}h(j)x_j fr(T)\\
&=h(i)\sum_{T|i}fr(T)\sum_{T|j}h(j)x_j
\end{aligned}\]
考虑把后半部分的那个东西也给提前算出来，记做\(g(T)\)。
那么要求的就变成了\(\displaystyle b_i=h(i)\sum_{T|i}fr(T)g(T)\)
令\(gr(T)=fr(T)g(T)\)，于是\(\displaystyle b_i=h(i)\sum_{T|i}gr(T)\)。
根据莫比乌斯反演可以得到：
\[\frac{b_i}{h(i)}=\sum_{T|i}gr(T)
\]
\[gr(i)=\sum_{T|i}\mu(\frac{i}{T})\frac{b_T}{h(T)}
\]
而\(gr\)是可以算出来的，\(fr\)也是可以算出来的，所以\(g\)也是可以算出来的。
而\(\displaystyle g(T)=\sum_{T|j}h(j)x_j\)，那么通过莫比乌斯反演可以算出\(h(j)x_j\)，直接除掉之后就能算出\(x_j\)了。
如果无解就是存在除法的时候出现了非零数除以\(0\)。
#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<cmath>
#include<algorithm>
#include<vector>
using namespace std;
#define ll long long
#define MOD 998244353
#define MAX 100100
inline int read()
{
int x=0;bool t=false;char ch=getchar();
while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();
if(ch=='-')t=true,ch=getchar();
while(ch<='9'&&ch>='0')x=x*10+ch-48,ch=getchar();
return t?-x:x;
}
int fpow(int a,int b){int s=1;while(b){if(b&1)s=1ll*s*a%MOD;a=1ll*a*a%MOD;b>>=1;}return s;}
int n,C,D,Q;
int pri[MAX],tot,mu[MAX];
bool zs[MAX];
void Sieve()
{
mu[1]=1;
for(int i=2;i<=n;++i)
{
if(!zs[i])pri[++tot]=i,mu[i]=-1;
for(int j=1;i*pri[j]<=n;++j)
{
zs[i*pri[j]]=true;
if(i%pri[j]==0)break;
mu[i*pri[j]]=-mu[i];
}
}
}
int b[MAX],gr[MAX],h[MAX],invh[MAX],f[MAX],g[MAX],fr[MAX],w[MAX],invfr[MAX];
int main()
{
n=read();C=read();D=read();Q=read();Sieve();
for(int i=1;i<=n;++i)h[i]=fpow(i,D),invh[i]=fpow(h[i],MOD-2);;
for(int i=1,p=(C-D+MOD-1)%(MOD-1);i<=n;++i)f[i]=fpow(i,p);
for(int i=1;i<=n;++i)
for(int j=i;j<=n;j+=i)
fr[j]=(0ll+fr[j]+mu[j/i]*f[i]+MOD)%MOD;
for(int i=1;i<=n;++i)invfr[i]=fpow(fr[i],MOD-2);
while(Q--)
{
for(int i=1;i<=n;++i)b[i]=read();
bool fl=true;
for(int i=1;i<=n;++i)w[i]=1ll*b[i]*invh[i]%MOD;
for(int i=1;i<=n;++i)if(b[i]&&!h[i])fl=false;
for(int i=1;i<=n;++i)gr[i]=0;
for(int i=1;i<=n;++i)
for(int j=i;j<=n;j+=i)
gr[j]=(0ll+gr[j]+mu[j/i]*w[i]+MOD)%MOD;
for(int i=1;i<=n;++i)g[i]=1ll*gr[i]*invfr[i]%MOD;
for(int i=1;i<=n;++i)if(gr[i]&&!invfr[i])fl=false;
for(int i=1;i<=n;++i)w[i]=0;
for(int i=1;i<=n;++i)
for(int j=i;j<=n;j+=i)
w[i]=(0ll+w[i]+mu[j/i]*g[j]+MOD)%MOD;
for(int i=1;i<=n;++i)if(w[i]&&!h[i])fl=false;
for(int i=1;i<=n;++i)w[i]=1ll*w[i]*invh[i]%MOD;
if(!fl){puts("-1");continue;}
for(int i=1;i<=n;++i)printf("%d ",w[i]);
puts("");
}
return 0;
}
