1. 程式人生 > >UOJ #62. 【UR #5】怎樣跑得更快

UOJ #62. 【UR #5】怎樣跑得更快

tint ret d+ void getchar() return isdigit clas 套路

題目分析

顯然不可能高斯消元。

考慮反演。

\(b_i=\sum\limits_{j=1}^n\gcd(i,j)^C\cdot \text{lcm}(i,j)^D\cdot x_j\)

\(b_i=\sum\limits_{j=1}^n\gcd(i,j)^C\cdot \frac{i^D\cdot j^D}{\gcd(i,j)^D}\cdot x_j\)

\(b_i=\sum\limits_{j=1}^n\gcd(i,j)^{C-D}\cdot i^D\cdot j^D\cdot x_j\)

實際上形如\(b_i=\sum\limits_{j=1}^nf(\gcd(i,j))\cdot g(i)\cdot h(j)\cdot x_j\)

都可以做。

我們按照套路化一下式子。

\(b_i=\sum\limits_{d|i}\sum\limits_{d|j}[\gcd(i,j)=d]\cdot f(d)\cdot g(i)\cdot h(j)\cdot x_j\)

\([\gcd(i,j)=d]\)換成\(\sum\limits_{k|\frac{\gcd(i,j)}{d}}\mu(k)\)

\(b_i=\sum\limits_{d|i}\sum\limits_{d|j}\sum\limits_{k|\frac{\gcd(i,j)}{d}}\mu(k)\cdot f(d)\cdot g(i)\cdot h(j)\cdot x_j\)

\(b_i=\sum\limits_{d|i}\sum\limits_{d|j}\sum\limits_{k \cdot d|\gcd(i,j)}\mu(k)\cdot f(d)\cdot g(i)\cdot h(j)\cdot x_j\)

\(b_i=\sum\limits_{T|i}\sum\limits_{T|j}\sum\limits_{d|T}\mu(\frac{T}{d})\cdot f(d)\cdot g(i)\cdot h(j)\cdot x_j\)

\(\frac{b_i}{g(i)}=\sum\limits_{T|i}\sum\limits_{T|j}\sum\limits_{d|T}\mu(\frac{T}{d})\cdot f(d)\cdot h(j)\cdot x_j\)

\(fr(T)=\sum\limits_{d|T}\mu(\frac{T}{d})\cdot f(d)\)

\(\frac{b_i}{g(i)}=\sum\limits_{T|i}\sum\limits_{T|j}fr(T)\cdot h(j)\cdot x_j\)

\(\frac{b_i}{g(i)}=\sum\limits_{T|i}fr(T)\sum\limits_{T|j}h(j)\cdot x_j\)

\(q(T)=\sum\limits_{T|j}h(j)\cdot x_j\)

\(\frac{b_i}{g(i)}=\sum\limits_{T|i}fr(T)\cdot q(T)\)

\(fr(i)\cdot q(i)=\sum\limits_{T|i}\mu(\frac{i}{T})\cdot \frac{b_T}{g(T)}\)

所以就可以求出\(q(i)\)了。

求出\(q(i)\)後,再次反演,

\(h(i)\cdot x_i=\sum\limits_{i|j}\mu(\frac{j}{i})\cdot q(j)\)

那麽就很容易求出\(x_i\)了。

註意一下無解的情況即可。

#include <bits/stdc++.h>
using namespace std;
inline int Getint(){
    int x=0,f=1;char ch=getchar();
    while(!isdigit(ch))ch!='-'?:f=-1,ch=getchar();
    while(isdigit(ch))x=x*10+ch-'0',ch=getchar();
    return x*f;
}
typedef long long ll;
const int Maxn=100005,mod=998244353,pmod=mod-1;
int n,q,c,d,b[Maxn],ans[Maxn],h[Maxn],fr[Maxn],g[Maxn];
int mu[Maxn],Prime[Maxn];
bool vis[Maxn];
ll Pow(ll x,ll k){
    ll ret=1;
    while(k){
        if(k&1)ret=ret*x%mod;
        k>>=1;x=x*x%mod;
    }
    return ret;
}
void init(){
    mu[1]=1;
    for(int i=2;i<=100000;i++){
        if(!vis[i]){Prime[++Prime[0]]=i;mu[i]=-1;}
        for(int j=1;j<=Prime[0]&&i*Prime[j]<=100000;j++){
            vis[i*Prime[j]]=1;
            if(i%Prime[j]==0){mu[i*Prime[j]]=0;break;}
            mu[i*Prime[j]]=-mu[i];
        }
    }
    int mi=((c-d)%pmod+pmod)%pmod;
    for(int i=1;i<=100000;i++){
        int tmp=Pow(i,mi);
        for(int j=1;i*j<=100000;j++)
            fr[i*j]=(fr[i*j]+(ll)mu[j]*tmp)%mod;
    }
    for(int i=1;i<=100000;i++)g[i]=Pow(i,d);
}
void solve(){
    memset(h,0,sizeof(h));
    memset(ans,0,sizeof(ans));
    int mi=((-d)%pmod+pmod)%pmod;
    for(int i=1;i<=n;i++)b[i]=(ll)b[i]*Pow(i,mi)%mod;
    for(int i=1;i<=n;i++)
        for(int j=1;i*j<=n;j++)
            h[i*j]=(h[i*j]+(ll)mu[j]*b[i])%mod;
    for(int i=1;i<=n;i++){
        if(fr[i]==0&&h[i]!=0){puts("-1");return;}
        h[i]=(ll)h[i]*Pow(fr[i],mod-2)%mod;
    }
    for(int i=1;i<=n;i++)
        for(int j=1;i*j<=n;j++)
            ans[i]=(ans[i]+(ll)mu[j]*h[i*j])%mod;
    for(int i=1;i<=n;i++){
        if(g[i]==0&&ans[i]!=0){puts("-1");return;}
        if(g[i])ans[i]=(ll)ans[i]*Pow(g[i],mod-2)%mod;
        else ans[i]=0;
    }
    for(int i=1;i<=n;i++)cout<<(ans[i]+mod)%mod<<" \n"[i==n];
}
int main(){
    n=Getint();c=Getint();d=Getint();q=Getint();
    init();
    while(q--){
        for(int i=1;i<=n;i++)b[i]=Getint();
        solve();
    }
}

UOJ #62. 【UR #5】怎樣跑得更快