[CF947E]Perpetual Subtraction

题目

题目背景
数学：一门严谨地说明 $\text{crashed}$ 是神的理论。

题目描述
传送门 to CF

思路

首先变换是一个线性变换。问题相当于是求矩阵的幂，乘一个向量。

有两种思考方向：一是解构这个矩阵，常用相似变换，比如变成置换；二是说明这个矩阵代表了特殊线性变换，比如卷积矩阵。后者一般可以抛开矩阵的背景，直接从别的角度刻画问题。

因此我试了试从 $\mathrm{GF}$ 角度推导。直接从过程出发，得到的结果涉及 $[x^m]$ 远处值，完蛋。当然后来 $\text{Tiw-Air-OAO}$ 告诉我可以用 $\mathrm{GF}$ 求解。话说这真的不算 $\mathrm{constructive}$ 吗 😱

考虑解构这个矩阵。该矩阵是三角阵，特征值一眼即知；手算几个特征向量，规律就很容易发现了。

若 $f_x=(-1{)}^{n-x}(\genfrac{}{}{0ex}{}{n}{x})(x⩽n)$ 则
$$\begin{array}{rl}{g}_x=\sum _{y⩾x}\frac{f_y}{y+1}& =\sum _{y⩾x}\frac{(-1{)}^{n-y}(\genfrac{}{}{0ex}{}{n}{y})}{y+1}\\ & =\sum _{y⩾x}\left(\genfrac{}{}{0ex}{}{n+1}{y+1}\right)(n+1{)}^{-1}(-1{)}^{n-y}\\ & =\frac{1}{n+1}\sum _{y⩾x}(-1{)}^{n-y}(\genfrac{}{}{0ex}{}{n}{y})+\sum _{y⩾x}(-1{)}^{n-y}\left(\genfrac{}{}{0ex}{}{n}{y+1}\right)\\ & =\frac{1}{n+1}\left(\genfrac{}{}{0ex}{}{n}{x}\right)(-1{)}^{n-x}=\frac{f_x}{n+1}\end{array}$$

所以这就是特征向量。对角化，等价于解 { η i } \{\eta_i\} {ηi​} 使得
$$f_i=\sum _{j⩾i}(-1{)}^{j-i}(\genfrac{}{}{0ex}{}{j}{i}){\eta }_j$$

我今天才知道这是 二项式反演 的标准形式 😓

或者你数学直觉足够好，可以指出：记 $F(x)$ 为 $f_i$ 的 $\mathrm{OGF}$，事实上 ${\eta }_i$ 的生成函数是 $F(x-1)$ 。因此可以很容易得到反演
$${\eta }_i=\sum _{j⩾i}\left(\genfrac{}{}{0ex}{}{j}{i}\right)f_j$$

因为 $\mathrm{H}(x+1)=F(x)$ 嘛。总复杂度 $\mathcal{O}(n\log n)$ 。

代码

#include <cstdio>
#include <algorithm> // Almighty XJX yyds!!!
#include <cstring>  // oracle: ZXY yydBUS!!!
#include <cctype> // Huge Egg Dog eats me!!!
using llong = long long;
# define rep(i,a,b) for(int i=(a); i<=(b); ++i)
# define drep(i,a,b) for(int i=(a); i>=(b); --i)
# define rep0(i,a,b) for(int i=(a); i!=(b); ++i)
inline int readint(){
	int a = 0, c = getchar(), f = 1;
	for(; !isdigit(c); c=getchar()) if(c == '-') f = -f;
	for(; isdigit(c); c=getchar()) a = a*10+(c^48);
	return a*f;
}

const int MOD = 998244353, LOGMOD = 30;
inline void modAddUp(int &x, const int &y){
	if((x += y) >= MOD) x -= MOD; // faster?
}
inline llong qkpow(llong b, int q){
	llong a = 1;
	for(; q; q>>=1,b=b*b%MOD) if(q&1) a = a*b%MOD;
	return a;
}
const int MAXN = 400000;
int g[LOGMOD], inv2[LOGMOD];
inline void prepareNtt(){
	int p = MOD-1, x = 0; inv2[0] = 1;
	for(inv2[1]=(MOD+1)>>1; !(p&1); p>>=1,++x)
		inv2[x+1] = int(llong(inv2[1])*inv2[x]%MOD);
	for(g[x]=int(qkpow(3,p)); x; --x)
		g[x-1] = int(llong(g[x])*g[x]%MOD);
}
void ntt(int a[], int n){
	for(int w=1<<n>>1,x=n; x; w>>=1,--x)
	for(int *p=a; p!=a+(1<<n); p+=(w<<1))
	for(int *i=p,*j=p+w,v=1; i!=p+w; ++i,++j,v=int(llong(v)*g[x]%MOD)){
		const llong t = llong((*i)+MOD-(*j))*v%MOD;
		modAddUp(*i,*j), *j = int(t);
	}
}
void dntt(int a[], int n){
	for(int w=1,x=1; x<=n; w<<=1,++x)
	for(int *p=a; p!=a+(1<<n); p+=(w<<1))
	for(int *i=p,*j=p+w,v=1; i!=p+w; ++i,++j,v=int(llong(v)*g[x]%MOD)){
		const int t = int(llong(*j)*v%MOD);
		modAddUp(*j=*i,MOD-t), modAddUp(*i,t);
	}
	std::reverse(a+1,a+(1<<n));
	for(int *i=a; i!=a+(1<<n); ++i)
		*i = int(llong(*i)*inv2[n]%MOD);
}
/// @return minimal k s.t. 2^k > n
inline int getNttLen(int n){
	return 32-__builtin_clz(n);
}
inline void array_mul(int a[], const int b[], int n){
	for(; n; --n,++a,++b) *a = int(llong(*a)*(*b)%MOD);
}

int a[MAXN], b[MAXN];
int main(){
	prepareNtt();
	int n = readint(); llong m; scanf("%lld",&m);
	for(int i=0,jc=1; i<=n; ++i,jc=int(llong(i)*jc%MOD)){
		a[i] = int(llong(readint())*jc%MOD);
		b[i] = int(qkpow(jc,MOD-2)); // inv
	}
	int len = getNttLen(n), zhishu = int(m%(MOD-1));
	std::reverse(b,b+n+1); // minus convolution
	ntt(a,len+1), ntt(b,len+1), array_mul(a,b,2<<len), dntt(a,len+1);
	for(int i=0,jc=1; i<=n; ++i,jc=int(llong(i)*jc%MOD)){
		llong inv = qkpow(i+1,MOD-2); // multiplier
		a[i] = int(qkpow(inv,zhishu)*a[i+n]%MOD);
		b[i] = int(qkpow(jc,MOD-2)); // reset
		if(i&1) b[i] = MOD-b[i]; // negate
	}
	memset(a+n+1, 0, ((2<<len)-n-1)<<2);
	memset(b+n+1, 0, ((2<<len)-n-1)<<2);
	std::reverse(b, b+n+1); // minus convolution
	ntt(a,len+1), ntt(b,len+1), array_mul(a,b,2<<len), dntt(a,len+1);
	for(int i=0,jc=1; i<=n; ++i,jc=int(llong(i)*jc%MOD))
		printf("%lld ",a[i+n]*qkpow(jc,MOD-2)%MOD);
	putchar('\n');
	return 0;
}