题解 [CF923E] Perpetual Subtraction

传送门

等退役了，我要去把 CC 协议通读一遍然后尝试搬一个基金会版式

先口胡（暴力）做法：
发现很矩阵乘，那就矩阵乘，没法优化（摊手

然后康题解，Dr.x义x 说它是可以优化的
于是尝试入门一些线性代数
link，link and link
所有可以表示为给定向量线性组合的向量的集合叫做张成的空间（Span）
那么单个向量张成的空间为一条直线
一个矩阵左乘一个列向量相当于对这个向量做一个线性变换
如果这个向量在线性变换后并没有离开其变换前所张成的空间，那么这个向量是关于这个矩阵的特征向量
虽然并未离开其变换前所张成的空间，但这个向量可能受到了一些拉伸，其被乘上的那个标量被称作这个特征向量的特征值
注意并不是所有变换都有对应的特征向量

那么要如何对一个给定的矩阵 \(A\) 计算特征向量 \(\vec{v}\) 和其特征值 \(\lambda\) 呢？
有

\[\begin{aligned} A\vec{v}&=\lambda\vec{v}\\ A\vec{v}&=\lambda I\vec{v}\\ (A-\lambda I)\vec{v}&=\vec{0}\\ \end{aligned}\]

注意 \(\vec{v}\neq\vec{0}\)，所以 \((A-\lambda I)\) 不是满秩的
否则会有 \(\vec{0}(A-\lambda I)^{-1}=\vec{v}\neq \vec{0}\)，这显然是不可能的
那么现在限制变成了 \(\det(A-\lambda I)=0\)
发现未知数 \(\lambda\) 集中于主对角线上，那么现在要解一个一元 \(n\) 次方程，就变成了[HEOI2012] Akai 的数学作业……？
并不是，能这么出题的话求特征值应该都基于矩阵的特殊性质
比如说这个题的矩阵长这样：

进行一些初等行变换后变成这样：

那么 \(n\) 个可能的 \(\lambda\) 就是显然的了
否则大概只能手解或当成上面那题来做

接下来需要根据特征值来还原特征向量 \(\vec{v}\)
设

\[\vec{v}=\begin{bmatrix} x_1\\x_1\\ \cdots \\ x_n\end{bmatrix} \]

又有

\[(A-\lambda I)\vec{v}=\vec{0} \]

那么就可以对 \(x_1\cdots x_n\) 列方程，然后高斯消元解了
但是你需要解 \(n\) 个特征向量，所以复杂度 \(O(n^4)\) 的
不过要是复杂度支持 \(O(n^4)\) 那为啥不直接快速幂啊
除非特征向量与输入无关可以本地挂机跑？
哦空间开不下啊。那没事了

现在特征向量可以求了
让我们来康康矩阵的基变换
我们知道一般的矩阵是以 \((1, 0, 0, \cdots, 0), (0, 1, 0, \cdots, 0), \cdots, (0, 0, 0, \cdots, 1)\) 为基的
但现在我要是想让它以一组别的向量为基呢？
感觉上面这篇博讲得很清楚了
大意就是发现从目标基变换回当前基是容易的
发现从目标基变换回当前基可以通过乘一个矩阵 \(A\) 实现
然后这个 \(A\) 实际上就是当前基表示的目标基，肯定是线性无关的
那么 \(A\) 有逆，而 \(A^{-1}\) 就是当前基转目标基的变换矩阵

现在矩阵也可以随便换基底了
来康康矩阵对角化是什么人间疾苦
发现由特征向量构成的矩阵大概长这样

这样的矩阵是很容易做快速幂的，像这样

那么假如我们使用一组特征向量作为基（特征基）
考虑在这组新的基意义下，原矩阵对应的变换将仅包括拉伸而不包括旋转
也就是说，基变换后的矩阵将是一个对角矩阵，且对角元为对应的特征值
那么这个矩阵是容易快速幂的了
只要先做基变换，再快速幂，最后换回原来的基即可
基变换一般要靠特殊性质优化复杂度

回到这个题
现在特征值已经求出来了
发现特征值很有规律，猜测特征向量也很有规律
打表应该可以发现第 \(i\) 个特征向量的第 \(j\) 维是

\[(-1)^{i+j}\binom{i}{j} \]

再对逆矩阵打表，发现第 \(i\) 列第 \(j\) 行是

\[\binom{i}{j} \]

那么现在要优化矩阵乘向量的过程
令基变换矩阵为 \(F\)，输入的列向量为 \(\vec{p}\)
那么要算

\[F(F^{-1}D^mF)F^{-1}p \]

利用结合率从右到左算
将 \(Fp\) 和 \(F^{-1}p\) 的式子写开，发现是减法卷积
然后中间的部分是对角矩阵，直接乘就好了
最终复杂度 \(O(n\log n)\)

点击查看代码

#include <bits/stdc++.h>
using namespace std;
#define INF 0x3f3f3f3f
#define N 600010
#define ll long long
//#define int long long

char buf[1<<21], *p1=buf, *p2=buf;
#define getchar() (p1==p2&&(p2=(p1=buf)+fread(buf, 1, 1<<21, stdin)), p1==p2?EOF:*p1++)
inline ll read() {
	ll ans=0, f=1; char c=getchar();
	while (!isdigit(c)) {if (c=='-') f=-f; c=getchar();}
	while (isdigit(c)) {ans=(ans<<3)+(ans<<1)+(c^48); c=getchar();}
	return ans*f;
}

int n;
int rev[N], bln, bct;
const ll mod=998244353, phi=mod-1, rt=3;
ll fac[N], inv[N], f[N], g[N], tem[N], m;
inline ll qpow(ll a, ll b) {ll ans=1; for (; b; a=a*a%mod,b>>=1) if (b&1) ans=ans*a%mod; return ans;}

void ntt(ll* a, int len, int op) {
	for (int i=0; i<len; ++i) if (i<rev[i]) swap(a[i], a[rev[i]]);
	ll w, wn, t;
	for (int i=1; i<len; i<<=1) {
		wn=qpow(rt, (op*phi/(i<<1)+phi)%phi);
		for (int j=0,step=i<<1; j<len; j+=step) {
			w=1;
			for (int k=j; k<j+i; ++k,w=w*wn%mod) {
				t=w*a[k+i]%mod;
				a[k+i]=(a[k]-t)%mod;
				a[k]=(a[k]+t)%mod;
			}
		}
	}
	if (op==-1) {
		ll inv=qpow(len, mod-2);
		for (int i=0; i<len; ++i) a[i]=a[i]*inv%mod;
	}
}

signed main()
{
	n=read(); m=read();
	for (int i=0; i<=n; ++i) f[i]=read();
	fac[0]=fac[1]=1; inv[0]=inv[1]=1;
	for (bln=1; bln<=2*n; bln<<=1,++bct);
	for (int i=0; i<bln; ++i) rev[i]=(rev[i>>1]>>1)|((i&1)<<(bct-1));
	for (int i=2; i<bln; ++i) fac[i]=fac[i-1]*i%mod;
	for (int i=2; i<bln; ++i) inv[i]=(mod-mod/i)*inv[mod%i]%mod;
	for (int i=2; i<bln; ++i) inv[i]=inv[i-1]*inv[i]%mod;
	for (int i=0; i<bln; ++i) f[i]=f[i]*fac[i]%mod;
	for (int i=0; i<=n; ++i) g[i]=inv[i];
	reverse(g+1, g+bln);
	ntt(f, bln, 1); ntt(g, bln, 1);
	for (int i=0; i<bln; ++i) f[i]=f[i]*g[i]%mod;
	ntt(f, bln, -1);
	// memset(tem, 0, sizeof(tem));
	// for (int i=0; i<bln; ++i)
	// 	for (int j=0; j<bln; ++j) if (i>=j)
	// 		tem[i-j]=(tem[i-j]+f[i]*g[j])%mod;
	// for (int i=0; i<bln; ++i) f[i]=tem[i];
	// cout<<"f: "; for (int i=0; i<bln; ++i) cout<<(f[i]%mod+mod)%mod<<' '; cout<<endl;
	for (int i=0; i<=n; ++i) f[i]=f[i]*inv[i]%mod;
	for (int i=n+1; i<bln; ++i) f[i]=0;
	for (int i=0; i<=n; ++i) f[i]=f[i]*qpow(qpow(i+1, mod-2), m)%mod;
	for (int i=0; i<=n; ++i) f[i]=(i&1?-1:1)*fac[i]*f[i]%mod;
	ntt(f, bln, 1);
	for (int i=0; i<bln; ++i) f[i]=f[i]*g[i]%mod;
	ntt(f, bln, -1);
	// memset(tem, 0, sizeof(tem));
	// for (int i=0; i<bln; ++i)
	// 	for (int j=0; j<bln; ++j) if (i>=j)
	// 		tem[i-j]=(tem[i-j]+f[i]*g[j])%mod;
	// for (int i=0; i<bln; ++i) f[i]=tem[i];
	for (int i=0; i<=n; ++i) f[i]=(i&1?-1:1)*inv[i]*f[i]%mod;
	// cout<<"f: "; for (int i=0; i<bln; ++i) cout<<(f[i]%mod+mod)%mod<<' '; cout<<endl;
	for (int i=0; i<=n; ++i) printf("%lld%c", (f[i]%mod+mod)%mod, " \n"[i==n]);

	return 0;
}