UOJ#62【UR #5】怎样跑得更快(反演)

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

题目描述

大力水手问禅师：“大师，我觉得我光有力气是不够的。比如我吃菠菜可以让力气更大，但是却没有提升跑步的速度。请问怎样才能跑得更快？我试过吃白菜，没有效果。”
禅师浅笑，答：“方法很简单，不过若想我教你，你先看看这道UOJ Round的C题。”
令 \(p = 998244353\)（\(7 \times 17 \times 2^{23} + 1\)，一个质数）。
给你整数 \(n, c, d\)。现在有整数 \(x_1, \dots, x_n\) 和 \(b_1, \dots, b_n\) 满足 \(0 \leq x_1, \dots, x_n, b_1, \dots, b_n \lt p\)，且对于 \(1 \leq i \leq n\) 满足：

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

其中 \(v \equiv u \pmod{p}\) 表示 \(v\) 和 \(u\) 除以 \(p\) 的余数相等。\(\gcd(i, j)\) 表示 \(i\) 和 \(j\) 的最大公约数，\(lcm(i, j)\) 表示 \(i\) 和 \(j\) 的最小公倍数。
有 \(q\) 个询问，每次给出 \(b_1, \dots, b_n\)，请你解出 \(x_1, \dots, x_n\) 的值。

数据范围

\[1 \le n \le 10^5, ~~~~q \le 3 \]

解题思路

好题 + 1

先说一个有启发的 \(\Theta(n^3+qn^2)\) 做法

直接将系数暴力算出直接跑高斯消元是 \(\Theta(n^3q)\) 的

有一个技巧，因为系数矩阵相同，所以对系数矩阵求逆即可

\[A * V = Ans\\ A^T*A*V = A^T*Ans\\ V =A^T*Ans \]

向量乘矩阵是 \(\Theta(n^2)\) 的

正解（参考vfk解法）

转化一波式子，变为

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

类似于

\[\sum_{j = 1}^{n} F(\gcd(i, j))\cdot h(i)*g(j) \cdot x_j \equiv b_i \pmod{p}\\ \]

的式子都可以做

同时设 \(F(x) = \sum_{i|d}f(i)\)，有

\[f(x) = \sum_{d|x}F(d)*\mu(\frac xd)\\ f(x) = F(x) - \sum_{d|x且d != x} f(d)\\ \]

两种做法均可在 \(\Theta(NlogN)\) 的时间复杂度内求出

将 F 换为 f 得到

\[\sum_{j =1}^n \sum_{d}[d|i][d|j]f(d)\cdot h(i)*g(j) \cdot x_j \equiv b_i \pmod{p}\\ \sum_{d|i}f(d)\sum_{j =1}^n [d|j]\cdot g(j) \cdot x_j \equiv b_i*h(i)^{-1} \pmod{p}\\ let ~~H_i = b_i\cdot h(i)^{-1},~~G_i= \sum_{i|j}g(j) \cdot x_j\\ \sum_{d|i}f(d)G_d=H_i\\ \]

这个东西，很熟悉，和刚才设 \(F(x) = \sum_{i|d}f(i)\) 一样可以递推过去，然后再对 \(G_i\) 进行一次递推

总复杂度 \(\Theta(qNlogN)\)

代码

#pragma GCC optimize(3, "inline")
#include <iostream>
#include <cstdio>
#include <cstring>
#define ll long long
using namespace std;

template <typename T>
void read(T &x) {
    x = 0; bool f = 0;
    char c = getchar();
    for (;!isdigit(c);c=getchar()) if (c=='-') f=1;
    for (;isdigit(c);c=getchar()) x=x*10+(c^48);
    if (f) x=-x;
}

const int P = 998244353;
int n, c, d, q;

ll fpw(ll x, ll mi) {
	ll res = 1;
	while (mi) {
		if (mi & 1) res = res * x % P;
		x = x * x % P, mi >>= 1;
	}
	return res;
}

const int N = 200500;
ll f[N], h[N], G[N], b[N], X[N], B[N];
void getf(void) {
	for (int i = 1;i <= n; i++) {
		f[i] += fpw(i, (P - 1 + c - d) % (P - 1));
		for (int j = i + i;j <= n; j += i) f[j] -= f[i];
		f[i] = (f[i] % P + P) % P, f[i] = fpw(f[i], P - 2);
	}
}

int main() {
	read(n), read(c), read(d), read(q);
	c %= (P - 1), d %= (P - 1), getf();
	for (int i = 1;i <= n; i++) 
		h[i] = fpw(i, P - 1 - d);
	while (q--) {
		for (int i = 1;i <= n; i++) 
			read(G[i]), G[i] = G[i] * h[i] % P;
		bool fl = 0;
		for (int i = 1;i <= n; i++) {
			for (int j = i + i;j <= n; j += i)
				G[j] -= G[i];
			G[i] = (G[i] % P + P) % P;
			if (G[i] && !f[i]) { fl = 1; break; }
			G[i] = G[i] * f[i] % P;
		}
		if (fl) { puts("-1"); continue; }
		for (int i = n;i >= 1; i--) {
			for (int j = i + i;j <= n; j += i)
				G[i] -= G[j];
			G[i] = (G[i] % P + P) % P;
		}
		for (int i = 1;i <= n; i++) printf ("%lld ", G[i] * h[i] % P);
		puts("");
	}
	return 0;
}

给个最短代码?

#include <bits/stdc++.h>
#define ll long long
#define rep() for(int i=1;i<=n;i++)
#define Rep() for(int j=i+i;j<=n;j+=i)
const int P=998244353;
const int N=200500;
ll f[N],h[N],G[N];
int n,c,d,q;
ll w(ll x,ll m) {
	ll r=1;
	for (;m;m>>=1,x=x*x%P)(m&1)&&(r=r*x%P);
	return r;
}
int main() {
	scanf ("%d%d%d%d",&n,&c,&d,&q);c%=(P-1),d%=(P-1);
	rep(){
		h[i]=w(i,P-1-d),f[i]+=w(i,(P-1+c-d)%(P-1));
		Rep()f[j]-=f[i];f[i]=(f[i]%P+P)%P,f[i]=w(f[i],P-2);
	}
	while (q--) {
		rep() scanf("%lld",G+i),G[i]=G[i]*h[i]%P;ll fl=0;
		rep(){
			Rep() G[j]-=G[i];G[i]=(G[i]%P+P)%P;
			if(G[i]&&!f[i]){fl=1;break;}G[i]=G[i]*f[i]%P;}
		if(fl){puts("-1");continue;}
		for (int i=n;i>=1;G[i]=(G[i]%P+P)%P,i--) Rep()G[i]-=G[j]; 
		rep()printf("%lld ",G[i]*h[i] % P);puts("");
	}
	return 0;
}