【UR #5】怎样跑得更快

\[b_i\equiv \sum_{j=1}^n \gcd(i,j)^c\cdot \operatorname {lcm(i,j)^d}\cdot x_j\pmod p \]

\[b_i\cdot i^{-d}\equiv \sum_{j=1}^n \gcd(i,j)^{c-d}\cdot j^d\cdot x_j\pmod p \]

换元，设 \(B_i=b_i\cdot i^{-d},X_i=x_i\cdot i^d\)

\[B_i\equiv \sum_{j=1}^n \gcd(i,j)^{c-d}\cdot X_j\pmod p \]

换一下元。

\[f_i=\sum_{i|j}X_j\\ X_{i}=\sum_{T=1}^{\frac{n}{i}}X_{iT}\cdot [T==1] \\ =\sum_{j=1}^{\frac{n}{i}}\mu(j)\sum_{k=1}^{\frac{n}{ij}}X_{ijk} \\ =\sum_{j=1}^{\frac{n}{i}}\mu(j)f_{i\cdot j}\\ \]

对于质数 \(x\)，\(f_x\) 很好求，考虑非质数

\[f_i=\sum_{j=1}^{\frac{n}{i}}X_{ij} \]

\[B_i\equiv \sum_{k|i} k^{c-d}\sum_{j=1}^{\frac{n}{k}}\cdot X_{jk}\pmod p \]

注意到这里需要容斥。

\[k_i\cdot f_{i}=B_i-\sum_{j|i\&j\ne i}k_j\cdot f_j \\ k_i=i^{c-d}-\sum_{j|i\And j\ne i}k_j,k_1=1 \]

#include <vector>
#include <cstdio>
#define LL long long
using namespace std;
const int Rea=1e5+3;
struct Rin
{
    char c;
    inline char gc()
    {
        static char rea[Rea];
        static char *head,*tail;
        return head==tail&&(tail=(head=rea)+fread(rea,1,Rea,stdin),head==tail)?EOF:*head++;
    }
    Rin&operator >>(int &x)
    {
        bool tag=false;x=0;
        for(c=gc();c>'9'||c<'0';c=gc())if(c=='-'){c=gc();tag=true;break;}
        for(;c>='0'&&c<='9';c=gc())x=(x<<1)+(x<<3)+(c^'0');if(tag)x=-x;return *this;
    }
    Rin&operator >>(LL &x)
    {
        bool tag=false;x=0;
        for(c=gc();c>'9'||c<'0';c=gc())if(c=='-'){c=gc();tag=true;break;}
        for(;c>='0'&&c<='9';c=gc())x=(x<<1)+(x<<3)+(c^'0');if(tag)x=-x;return *this;
    }
}rin;
inline void jh(int &x,int &y){if(x^y)x^=y^=x^=y;return;}
inline void jh(LL &x,LL &y){if(x^y)x^=y^=x^=y;return;}
inline int min(int x,int y){return x<y?x:y;}
inline int max(int x,int y){return x>y?x:y;}
inline LL min(LL x,LL y){return x<y?x:y;}
inline LL max(LL x,LL y){return x>y?x:y;}

const int N=3e5+3;
const int M=998244353;
#define lmod(x) if(x<0)x+=M
#define rmod(x) if(x>=M)x-=M
inline int prpr(int x,int y){return 1LL*x*y%M;}
inline int ksm(int x,int y){int ans=1;for(;y;y>>=1){if(y&1)ans=prpr(ans,x);x=prpr(x,x);}return ans;}
inline int inside(int x,int M){return (x%M+M)%M;}

int mu[N];
bool pri[N];
vector<int>prime;

int n,c,d,q;
int B[N];
int X[N];
int f[N];
int k[N];
void work()
{
    for(int i=1;i<=n;i++)B[i]=X[i]=f[i]=0;
    for(int i=1;i<=n;i++)rin>>B[i],B[i]=prpr(B[i],ksm(i,inside(-d,M-1)));
    for(int i=1;i<=n;i++)f[i]=B[i];
    for(int i=1;i<=n;i++)for(int j=2;i*j<=n;j++){f[i*j]-=f[i];lmod(f[i*j]);}
    for(int i=1;i<=n;i++){if(f[i]>0&&k[i]==0){puts("-1");return;}f[i]=prpr(f[i],ksm(k[i],M-2));}
    for(int i=1;i<=n;i++)for(int j=1;i*j<=n;j++){X[i]+=mu[j]*f[i*j];lmod(X[i]);rmod(X[i]);}
    for(int i=1;i<=n;i++)printf("%d ",prpr(X[i],ksm(ksm(i,d),M-2)));puts("");
    return;
}

void Read(){rin>>n>>c>>d>>q;return;}
void Init()
{
    int s=inside(c-d,M-1);mu[1]=1;
    for(int i=1;i<=n;i++)k[i]=ksm(i,s);
    for(int i=1;i<=n;i++)for(int j=2;i*j<=n;j++){k[i*j]-=k[i];lmod(k[i*j]);}
    for(int i=2;i<=n;i++)
    {
        if(!pri[i])prime.push_back(i),mu[i]=-1;
        for(auto j:prime){int now=i*j;if(now>n)break;pri[now]=true;if(i%j==0)break;mu[i*j]=-mu[i];}
    }
    return;
}
void Work(){for(;q;q--)work();return;}
int main()
{
    Read();
    Init();
    Work();
    return 0;
}