线性代数相关算法小记
高斯消元
高斯消元是对矩阵行化简的算法,可以化成阶梯型或者简化阶梯型。《线性代数及其应用》给出的步骤如下:
- 选取最左边的非零列;
- 在该列中任意选取一个非零元素,通过对换变换将该行移到最上面;
- 通过倍加变换将下面的行的该列元素全部变成 \(0\);
- 暂时不管该行(即第一行),将剩下的子矩阵代回第一步进行循环,直到没有非零列(矩阵内全为 \(0\))为止。此时该矩阵已经是阶梯型矩阵;
- 对每个非零行先用倍加变换将先导元素上面的元素全部变成 \(0\),再用倍乘变换将先导元素变成 \(1\)。此时该矩阵已经是简化阶梯型矩阵。
复杂度瓶颈在于这若干个倍加变换。设矩阵 \(A\) 为 \(m\times n\),一次任意初等行变换的复杂度显然都是 \(\mathrm O(n)\)。不管目标是阶梯型还是简化阶梯型,对每个主元位置的产生,都需要将 \(\mathrm O(m)\) 个行的该列位置清零,所以处理一个主元位置的复杂度是 \(\mathrm O(mn)\)。然后一共有 \(\operatorname{rank}A\) 个主元,所以总复杂度是 \(\mathrm O(mn\operatorname{rank}A)\),有 \(\operatorname{rank}A\leq\min(m,n)\)(但可能有特殊情况使得矩阵的秩很小)。如果认为 \(m,n\) 同阶,复杂度可以粗略的认为是 \(\mathrm O\!\left(n^3\right)\)。
但是书上是为了便于说明,才刻意搞的先化成阶梯型,然后再化成简化阶梯型。事实上,如果你最终目标是简化阶梯型的话,完全可以把第 5 步合到第 3 步里一起搞。
高斯消元还有个比较重要的关于模数的问题,分三类讨论一下:
- 无模数,那肯定是要在实数上搞了,也没什么好说的。只是你看第五步,先把该列其它元素搞成 \(0\),再把自己搞成 \(1\),比反过来的精度误差要小
(其实也没小到哪去的说)。upd 2022.2.22:有一种误差更小的做法:只要每次找到主元后,找到当前主元所在列绝对值最大的一个元素与该行交换,效果极好。这个 trick 在《线性代数及其应用》的数值计算注解里有提到,然而我数值计算注解都没看( - 模数为质数。这是比较爽的情况,不用考虑精度误差,而且原来的线性运算在带模之后仍然保持一切优良性质。有的时候有元素模成零了,但实际上不为零,于是可能出现的情况是在不模的时候矩阵的各列线性无关性、各列生成 \(\R^m\) 性要比模的时候好,但模的时候它就是没有这些性质,不需要考虑别的。倍乘和倍加变换需要求逆元,看上去复杂度多个 log,但其实只需要求主元的逆元,于是可以把 log 给平行化掉,相当于复杂度不变。
- 任意模数。这就比较坑爹了,它破坏了一些线性代数性质:倍乘和倍加变换可能没有逆操作,导致行等价关系不对称,因某些数没有逆元(你可能要说了,模数为质数的情况下,\(0\) 不也没有逆元么?但实数情况下 \(0\) 也没有逆元……所以并没有损失性质)。但前辈们依然勉强找到了化成阶梯型的方法:倍加变换的根本目的是把另一行的某列的位置变成 \(0\),避免求逆元的情况下可以类似欧几里得算法的辗转相除,终究能把目标位置消成 \(0\),并且只使用对换和倍加变换。这样看似复杂度又多了个 log,但是每个主元的总辗转次数其实是 log 的,所以这个 log 又可以平行化掉。至于要化成简化阶梯型,那是真的做不到了。第一,往上消零的话用上述辗转相除不可取,因为不难发现辗转相除的一个必要条件是对其它列都没有要求,而往上消零还需要保持上面行的先导元素位置不变;第二,倍乘变换是真的没办法用奇技淫巧替代了……
mol ban tea:
- 解线性方程组,实数运算,如果最后一列有主元位置就无解,否则如果有自由变量就无限解,否则唯一解,code。
- 矩阵求逆,只需要放个 \(I_n\) 在右边,然后对左边行化简即可,质数模数不慌,code。
- 行列式求值,任意模数,但是求行列式只需要化成阶梯型然后把主对角线乘起来,不慌,注意对换操作只有交换不同行才需要取反,另外这题很卡常,倍加变换如果倍数是 \(0\) 就
return
这个剪枝加上去就能过了,code。
线性基
给定一些数,考虑所有子集的异或和的集合。异或的话可以把数按二进制位拆成向量,那么异或就是 \(\bmod 2\) 意义下的加法。所以这个集合是线性空间,线性基就是这个集合的一组基。考虑求这组基。
我们知道,将一个矩阵 \(A\) 行化简,非零行是 \(\operatorname{Row}A\) 的一组基,主元列的列号对应的原来的向量是 \(\operatorname{Col}A\) 的一组基。所以我们求 \(A\) 的列的基有两种方法:一种是看成 \(\operatorname{Row}A^{\mathrm T}\),一种是看成 \(\operatorname{Col}A\)。我们采用前者求线性基,因为求出来有良好的性质(阶梯型 / 简化阶梯型的性质呗),况且本来就不需要基包含于原集合。
线性基一般采用动态高斯消元的方法,动态插入数,按道理插入一次是 \(\mathrm O\!\left(\log^2v\right)\)(其中 \(v\) 是值域)的,但是行变换可以直接异或(相当于 bitmask
优化),所以插入一次是 \(\mathrm O(\log v)\) 的。我们也选择不用 vector
来维护,而对每个主元位置维护对应的行,如果暂时是空的就是 \(0\),这样好写很多,不用 insert
之类,还可以直接提取想要的主元位置对应行。而且我们并不需要消成简化阶梯型,普通阶梯型就可以了。折腾一番之后代码变成了这个鬼样子(是不是炒鸡好写!):
void insert(int x){
for(int i=logv;~i;i--)if(x>>i&1)
if(b[i])x^=b[i];
else{b[i]=x;break;}
}
这种方法构造出来的线性基的特殊性质就是阶梯型,贪心起来特爽。线性基大小显然是 \(\mathrm O(\log v)\) 吧。
能解决的一些问题(都特别简单啦,都是线代基础):
- 数集线性相关性:线性基大小是否等于原数集。
- 生成空间大小:那显然就是 \(2^{|\mathcal B|}\)。
- 任选子集,异或和异或上 \(x\) 的最大 / 小值:从高往低按位贪心,跟据阶梯型性。
- 异或和等于 \(x\) 的子集数:按高斯消元解方程模拟。
- 异或和不超过 \(x\) 的子集数:按数位 DP 统计方式。
矩阵树定理
用来计算无向图生成树个数或者有向图的内 / 外向树个数。
下面允许图有重边,但不允许有自环:
- 无向图:设 \(D\) 为无向图的度数矩阵(对角线是度数,其它位置是 \(0\)),\(T\) 为邻接矩阵,基尔霍夫矩阵为 \(K=D-T\),则生成树数量为 \(K\) 的任意一个 \((i,i)\) 余子式(都相等)。
- 有向图:设 \(iD\) 为入度矩阵,\(oD\) 为出度矩阵,\(T\) 依然是邻接矩阵,基尔霍夫矩阵 \(iK=iD-T,oK=oD-T\),则以 \(r\) 为根的外向树个数为 \(iK\) 的 \((r,r)\) 余子式,内向树个数为 \(oK\) 的 \((r,r)\) 余子式。
证明的话就略了吧。网上找到的主要有两种证法,一种是基于关联矩阵和柯西-比内公式的证法,较能反应矩阵树定理的本质,但只能证无向图情况;另一种是对边数归纳,有向无向都能证,但不能反应本质。
如果图中有自环要手动去掉,因为有没有自环对生成树数量毫无影响,而矩阵树算法算上自环的话明显会受到干扰。如果要求有向图以所有点为根的内 / 外向树的个数和,那只能四方地算,因为 \(iK,oK\)(包括无向图的 \(K\))都不满秩,无法套用 \(A^*=|A|A^{-1}\) 来算伴随矩阵(upd 2022.8.18:是能三方的,见 2022 多校第十场 1010。而且这题是我做的)。有一个利用矩阵树定理反过来证明不满秩的方法:考虑再加一个点不连边,此时图不连通,生成树数量显然为 \(0\),用矩阵树定理计算的话去掉新加点所在行列得到原矩阵,行列式就为 \(0\) 了。
扩展:它可以带权,此时度数矩阵、邻接矩阵里装的都是边权和,这样矩阵树定理算出来的是所有生成树的边权积的边权和(仔细想一下,重边也恰可以用这个扩展来理解)。边权甚至可以为多项式(GF)或者其它奇怪的东西(需要定义加法、乘法,并且有某些良好的性质)放到行列式里去消元,至于为什么能这样以后再深究。
这样一来我们有办法计算所有生成树的边权和的和(其实可以枚举边缩边跑矩阵树,但那样复杂度最高是五方的):注意到 \(\prod\limits_{e\in T}(w_ex+1)\bmod x^2\) 的一次项就是生成树 \(T\) 的边权和。那么可以矩阵树套一次多项式,加减乘都是显然的,除法的话就是 \(\bmod x^2\) 意义下求逆,稍微解个方程知道 \(\dfrac{ax+b}{cx+d}\equiv\dfrac{ad-cb}{d^2}x+\dfrac bd\),如果 \(c\neq 0,d=0\) 就挂了,那怎么办呢?解决方法是:如果该列只有 \(ax\) 形状的多项式,那直接可以消,否则一定能找到 \(ax+b(b\neq 0)\)(但大多数情况下这种情况不会出现,概率是不是太小了)。
BEST 定理
实现基于矩阵树,但理论不基于矩阵树的一个定理,用来求有向图欧拉回路数量。
定理内容:连通欧拉图中从 \(s\) 出发的欧拉回路数量为 \(T\deg_s\prod\limits_{i\in V}(\deg_i-1)!\),其中 \(T\) 是以 \(s\) 为根的内 / 外向树个数(由于是欧拉图,所以所有点入度等于出度,易得两者相等)。
BEST 定理的证明
对一棵内向树,给每个点的出边集指定一个顺序,满足对非根节点它的树边出边是最后一个,构造路径:从 \(s\) 出发每到一个点沿着顺序中最前面还没走的边走下去,直到走不下去了(当前点的出边全部被走过了)停止。我们只需要证明对每个 (内向树, 顺序) 二元组构造出来的路径是欧拉回路,以及每个起点为 \(s\) 的欧拉回路有恰好一个二元组映射到它(即 (内向树, 顺序) 二元组与起点为 \(s\) 的欧拉回路形成双射)。
为证明前者,只要证不会在非根节点走不下去,以及在根节点走不下去时已经访问完所有边。
- 对非根节点,每次刚到达它的时候,显然它的已访问入边比已访问出边多 \(1\)。而由于原图是欧拉图,所以入度等于出度,此时必然还有出边没访问,不可能走不下去。事实上这个结论可以脱离「内向树」而存在,就是说在欧拉图上随便走(每条边只能访问一次)都不可能在非起点节点停下。
- 如果在根节点停下时还存在边未访问,任何边都一定是某个点的出边(废话),那么该点的出边顺序中位于这条边后面的边也没被访问。而树边一定在顺序的最后,所以必定存在某条树边 \((x,fa_x)\) 未被访问。那么 \(fa_x\) 的入边没有被访问完,跟据入度等于出度,它一定有出边没访问完,所以 \(\left(fa_x,fa_{fa_x}\right)\) 也一定没被访问。如此一直向上蔓延,知道根节点有入边没被访问,跟据入度等于出度,它一定有出边还没访问,就不会停下来,矛盾!
接下来证明后者。那显然对某个起点为 \(s\) 的欧拉回路,最多只可能有一棵内向树映射到它,只要证取出每个非起点节点最后访问的出边一定构成内向树即可。一个必要条件是基图构成无向有根树,而跟据每个非根节点恰有一条出边易证无向有根树恰好存在一种合法的定向方式,并且得到内向树,所以也是充分条件。\(n\) 个点 \(n-1\) 条无向边,只要证无环即可知它们构成树。如果有环的话,跟据每个非根节点恰有一条出边可知无向环还原成有向边之后,一定是有向环,而且根节点不会在环上(因为根没有出边)。取这个有向环上在欧拉回路中第一条被访问的边 \((x,y)\),那么环上其他边还没被访问,此时 \(x\) 的树边出边已经被访问,说明 \(x\) 所有出边都已经被访问,跟据入度等于出度知道 \(x\) 的所有入边也已被访问,而显然存在环上的一条未访问的 \(x\) 的入边,矛盾!
BEST 定理得证!
如果不限定起点的话,从 \(s\) 出发的欧拉回路每 \(\deg_s\) 形成一个循环同构,所以数量是 \(T\prod\limits_{i\in V}(\deg_i-1)!\)(无论选取哪个 \(s\) 都得到这个结果,恰体现了合理性)。如果不算循环同构的话,每 \(m\) 个欧拉回路构成一个循环同构,总数是 \(Tm\prod\limits_{i\in V}(\deg_i-1)!=T\sum\limits_{i\in V}\deg_i\prod\limits_{i\in V}(\deg_i-1)!\),很合理对否。
考虑实现的时候需要注意的事情:
- 需要判断是否欧拉图,不是欧拉图的话使用 BEST 会出错,应该直接输出
0
跑路。 - 自环问题,这个在矩阵树的时候就说过了,要手动去自环消元,但是欧拉回路数量计算公式 \(\deg\) 要算上自环。
- 孤立点问题,出现孤立点的时候导致图不连通,导致 \(T=0\),但是是否欧拉图只关心边的连通性。此时应该把孤立点直接删除,一个比较方便的操作是直接把基尔霍夫的 \(K_{i,i}\) 由 \(0\) 改成 \(1\),跟据行列式的按行 / 列展开易证正确性。但是此时如果指定起点是孤立点的时候,答案便是 \(0\) 了;但是又有一个特例(烦吧!):图上没有边时,答案是 \(1\)。
lgv 引理
在一个 DAG 上,有一个起点集合 \(A\) 和终点集合 \(B\),满足 \(|A|=|B|\)。一条路径 \(P\) 上边权的积是 \(P\) 的权值 \(w(P)\)。设 \(e(x,y)=\sum\limits_{P:x\to y}w(p)\),即 \(x\to y\) 所有路径的权值和。设 \(A\to B\) 的一个不相交路径方案为对于任意的 \(1\sim |A|\) 排列 \(p\) 的 \(\{P_i:A_i\to B_{p_i}\}\) 这样一个两两不相交的路径集合,它的权值为 \(w(P)=\prod\limits_{i=1}^{|A|}w(P_i)\)。
lgv 引理:设 \(|A|\) 阶方阵 \(M\) 满足 \(M_{i,j}=e(A_i,B_j)\),则有
即 \(A\to B\) 的所有不相交路径方案的权值代数和。
证明:将行列式按第一定义展开,它显然等于 \(A\to B\) 所有路径方案(不保证相交)的权值代数和。那为什么它又等于不相交的呢?我们只需要在偶相交路径方案集合和奇相交路径方案集合之间建立一个双射,并且每组映射满足原像和像权值相等,这样所有相交路径方案就可以抵消掉了。构造映射的方案很简单:对任意相交路径方案,找到第一个交点(这里第一个指的是以相交甲方为第一关键字、相交乙方为第二关键字、交点离两个起点的距离为第三关键字的第一个),将交点后的两条路径互换身份(就像 P1007 一样),这样排列中交换了元素,奇偶性就变化了。易证双射性。
如果想计数的话,只需要把边权都设成 \(1\) 即可。与矩阵树定理类似,边权可以是多项式或其它类型,放到行列式里跑。这样一来把边权积变成边权和显然也是可以做的了。甚至你发现 \(e\) 的计算和行列式计算是分开的,所以甚至可以解决「路径权值是边权积,路径方案权值是路径权值和」这样的问题。
\(e\) 函数一般很好算,毕竟是 DAG,复杂度瓶颈一般在高斯消元。lgv 引理的题有一大半都是只有排列 \(p=[1,2,\cdots,|A|]\) 的不相交路径方案权值非零(包括某些网格图的题),它们是用求代数和的方法来求数值和,虽然有点杀鸡用牛刀的感觉,但并没有更好的方法。洛谷上的模板题感觉并不模板,就当习题了。