网课学习笔记——矩阵行列式与高斯消元。
1行列式与基本概念
1.1一些概念与性质
1.1.1概念
- 逆序对( \(\tau\) )：设 A 为一个有 n 个数字的有序集 (n>1)，其中所有数字各不相同。
如果存在正整数 i, j 使得 1 ≤ i < j ≤ n 而且 A[i] > A[j]，则 <A[i], A[j]> 这个有序对称为 A 的一个逆序对。 - 排列：一般地，从n个不同元素中取出m（m≤n）个元素，按照一定的顺序排成一列，叫做从n个元素中取出m个元素的一个排列。特别地，当m=n时，这个排列被称作全排列，这个全排列被称作n阶排列。
如无特殊说明，一下排列都指全排列。
- 对换：交换排列的任意两个数，得到另一组排列。这个操作叫做对换。
- 如果一个排列的逆序对数为奇数，则称这个排列为奇排列，如果为偶数，则称这个排列为偶排列。
1.1.2性质
-
对换改变排列的奇偶性。
证明：设排列\(a_i(1{\le}i{\le}n)\)，易得如果交换两个相邻得数\(a_i,a_{i-1}\)则满足上述性质，因为除了\(a_i,a_{i-1}\)的之间逆序对数改变\(a_i\)左边的数和\(a_{i+1}\) 右边的数对这两个数的逆序对数不改变。如果交换的是两个不相邻的数，即\(a_1,a_2,...a_{i-1},a_i,a_{i+1},...a_{j-1},a_{j},a_{j+1}...a_n\)现在要交换\(a_i\)和\(a_j\)，那么可以看做是先把\(a_i\)做j-i次相邻的对换得到\(a_1,a_2,...a_{i-1},a_{i+1},...a_{j-1},a_{j},a_i,a_{j+1}...a_n\)，注意\(a_i\)的位置变换，在这个基础上，在进行j-i-1次相邻对换，得到\(a_1,a_2,...a_{i-1},a_{j},a_{i+1},...a_{j-1},a_i,a_{j+1}...a_n\)。这样，总共进行了2i-2j-1次对换，因为i、j为整数，那么2i-2j-1肯定为奇数，奇数次改变排列的奇偶性，则对换肯定改变排列的奇偶性。
-
在全部n阶排列中，奇、偶排列各占一半。
证明：对于一个奇排列，交换\(a_1,a_2\)一定能够得到一个对应的偶排列，而对于一个偶排列，可以通过相同的方式得到一个奇排列，所以奇排列与偶排列是一一对应的。故得证。 -
任意一个排列可以经过n次对换变为自然排列（n阶自然排列为\(1,2,3,4,5...n\)），且n与原排列奇偶性相同。
证明：首先前半句话肯定是对的，否则所有基于对换的排序算法都是错的。考虑每一次对换都改变排列的奇偶性，并且自然排列是一个偶排列（逆序对数为0），所以n一定与原排列奇偶性相同，这是显然的。
1.2行列式定义、概念和定理
1.2.1定义
N阶行列式由\(N^2\)个数\(a_{ij}(1{\le}i,j{\le},n)\)组成。
1.2.2行列式的值
求和公式:
这里sgn函数为如果该排列为偶，则返回1，否则则返回-1。即\(sgn(j1j2j3{\cdots}jn)=(-1)^{\tau(j_1,j_2,...j_n)}\)
j是一个1至n的排列
上面这个求和公式用普通话翻译过来就是，枚举j的所有全排列，用a去乘，把所得结果相加，实质是从每一行，每一列中都取出一个数来，相乘，把枚举所得到的所有结果相加，就是最终答案。
1.2.3定理
1.行列互换，值不变。
证明，个人认为这个是显然的。行列互换的意思是
变成
本质是以为枚举全排列还是那几个数。
2.用一个数乘行列式的某行等于用次数乘此行列式。即
注意这里可以不是第2行，可以把k乘到任意行上去。
证明：
根据求和公式:
在两边同时乘上k可以得到
利用乘法分配律，我们把k乘进去那么做加法的每一项都将会有且只有一个k，这与行列式
进行求和的答案是一样的，因为求和的本质是在每一行每一列取一个数，这样取的话，求和每一项都有且只有一个k。
3.如果行列式中某一行等于两个行列式之和，则这个行列式等于两个行列式之和，这两个行列式分别以这两组数为改行，其余各行相同。什么意思？我们画图来理解一下：
有行列式
（这里还是以第2行为例，其实可以是任意行），如果有\(b_{2j}+c_{2j}=a_{2j}(1{\le}j{\le}n)\) 那么就有
证明：把这三个行列式的求和公式写一下，会发现对右边可以用乘法分配律，这里不再赘述。
4.对换行列式两行位置，行列式反号。
证明：设第i行我们取的数为\(j_i\)那么我们假设随着行列式两行交换，所有我们枚举的全排列对应的\(j_i,j_k\) 这两个数也互换，易知，交换后的全排列没有重复，但所有的全排列奇偶性全部取反，由行列式计算公式可以得到，行列式反号。
5.如果行列式中有两行成比例，则行列式等于0。
证明：设这个行列式第i行与第j行成比例，则有\(a_{ik}=u*{a_{jk}} (1{\le}k{\le}n)\)，那么根据性质2，可以把这个u提出来，提出来后，这个行列式有相等的两行。设这个行列式的值为A，我们交换这两个相等的两行，根据定理4，交换后行列式的值为-A，但是又因为交换前后两个行列式完全等同，所以有A=-A，由此可得A=0。
6.把一行乘以某数加到另一行，行列式的值不变。
加完后的这个行列式，根据定理3，可以把它看作原行列式，和另一行列式A的和，由定理5得A=0，所以得证。
这是行列式的重要定理，下面我们用这些定理讲一讲高斯消元。
1.3高斯消元
本质上高斯消元是是我们计算行列式的一种方法。如果我们用计算公式来求行列式的话，很明显该算法复杂度可以达到\(O(n!*n)\)(n!枚举全排列，而n来计算）
如果这个行列式的从左上到右下的对角线一下全部为0，即
这个行列式的值如何计算？注意，任何有0参与的乘法结果都将是0，所以在选数计算的过程中不用考虑选0的情况，只需要考虑不选0的情况，又根据每行每列只取一个数的原则，我们只有一种方案可以选，即对角线，所以这个行列式的值就应该等于\(\prod\limits_{i=1}^na_{ii}\)。所以高斯消元所做的工作就是把一个行列式转化成上述行列式的形式。
如何转化？用行列式定义6。
考虑用\(a_{11}\)去消第一列，我们要做的就是先算个比值，即\(\frac{a_{21}}{a_{11}}\)用第二行的每一个数，都减去第一行同一列的数乘上这个比值，就可以在把\(s_{21}\)消成0的同时，保证行列式的值不变。由此用\(a_{11}\)去消第一列，再用\(a_{22}\)去消第二列，最终可以得到我们想要的行列式，再把对角线相乘即可得到行列式的值。
1.3.1代码1
#define dd double
#define N 101
dd z[N][N];int n;
inline dd guass(){
dd x=1.0;
for(int a=1;a<=n;a++){
for(int b=a+1;b<=n;b++){
dd d=z[b][a]/z[a][a];
for(int c=1;c<=n;c++) z[b][c]-=k*z[a][c];
}
}
for(int i=1;i<=n;i++) x*=z[i][i];
return x;
}
通过这个代码我们其实可以知道这个算法的时间复杂度为\(o(n^3)\)
但其实这个算法还有点小问题，就是z[a][a]可能原本就是个0，除以一个0这个程序就无法运行。我们利用定理4，把一个此位置非0的一行对换上来，行列式值取反。同时，如果这一列全都是0，那m这个行列式的值就肯定是0了，所以我们在兑换前加一步处理0，同时判断一列全为0的情况。即：(众所周知，实数中的0位1e-8）
#define dd double
#define N 101
dd a[N][N];int n;
inline dd guass(){
dd x=1;
for(int a=1;a<=n;a++){
for(int b=a;b<=n;b++)
if(fabs(z[b][a])>1e-8){//
if(b==a) break;
x=-x;//
for(int c=1;c<=n;c++) swap(z[a][c],z[b][c]);
break;//
}
if(fabs(z[a][a])<=1e-8) return 0.0;//
for(int b=a+1;b<=n;b++){
dd k=z[b][a]/z[a][a];
for(int c=1;c<=n;c++) z[b][c]-=k*z[a][c];
}
}
for(int i=1;i<=n;i++) x*=z[i][i];
return x;
}
时间复杂度为\(o(n^3)\)
这其中加了“//”的都是容易忘的，尤其是fabs这个实数取绝对值函数，一定不要忘。（还有定理四中行列式值取反也不要忘）
但是如果你要拿这个模板去过洛谷的高斯消元是过不去的，其原因在于，这个的精度不准的太严重了，上面这个代码，在竞赛中也是不会用的。试想一下，如果一个全部是整数的行列式，其值也一定是整数，但是我们用这个算法算出来的值却一定是一个实数。怎样优化？回答这个问题之前，先回答下面这个问题。
1.3.2优化
考虑一个实数1.12341214235341...,对他进行两个操作，乘100000000，和除以100000000，都存入一个double里，那个结果精度更高？
答案可能会让你大吃一惊————除以100000000的精度会更高。
原因：
double能存的下的数是一定的，也就是说，double的空间是一定的，对于乘100000000来说，他的空间有一部分给了整数部分，而对于除以10000000来说，整数部分就是个0，它大部分的空间都给了小数部分，所以除以100000000的精度会更高。也就是说，除以的这个数，绝对值越大，对精度的贡献也就越高。
对于上面的程序，我们的精度消耗在哪里了呢？显然:dd k=z[b][a]/z[a][a];这行代码消耗的精度更多。我们总希望z[a][a]的绝对值越大。所以我们就想办法把绝对值越大的弄到上面来，于是我们就从a开始，从上往下扫一遍让z[a][a]这个地方的绝对值最大，同时如果取了最大绝对之后，这个位置上的值仍然为0，那么说明这一列都是0。这种算法叫做主元高斯消元法。代码：
#define dd double
#define N 101
dd z[N][N];int n;
inline dd guass(){
dd x=1;
for(int a=1;a<=n;a++){
for(int b=a+1;b<=n;b++)
if(fabs(z[b][a])>fabs(z[a][a]){//
x=-x;//
for(int c=1;c<=n;c++) swap(z[a][c],z[b][c]);
}
if(fabs(z[a][a])<=1e-8) return 0.0;//
for(int b=a+1;b<=n;b++){
dd k=z[b][a]/z[a][a];
for(int c=1;c<=n;c++) z[b][c]-=k*z[a][c];
}
}
for(int i=1;i<=n;i++) x*=z[i][i];
return x;
}
这个代码可以过掉洛谷上的高斯消元模板。
但实际上，这个代码并没有解决掉精度问题，为了一劳永逸的解决精度问题，我们引出辗转相消高斯消元。
1.3.3辗转相消高斯消元
我们参考一下辗转相除求最大公因数的思路。通过gcd(a,b)=gcd(b,a%b)后，知道a能整除b，也就是a%b等于0，如果我们类似这个过程，对\(a_{ii}\)和\(a_{i+1i}\)进行辗转相消，以达到使\(a_{i+1i}\)等于0的目的，因为参与运算的都是整数，所以就不存在精度问题。代码：
int z[N][N];int n;
inline int guass(){
int x=1;
for(int a=1;a<=n;a++)
for(int b=a+1;b<=n;b++)
while(z[b][a]!=0){
int k=z[a][a]/z[b][a];
for(int c=1;c<=n;c++) z[a][c]-=z[b][c]*k;
for(int c=1;c<=n;c++) swap(z[a][c],z[b][c]);
x=-x;
}
for(int i=1;i<=n;i++) x*=z[i][i];
return x;
}
都是重灾区。。。一定注意辗转相消完全模拟了辗转相除的过程，一定是用第b行取消第a行，while下前两行相当于取模，后一行交换，在此过程中注意保证行列式的值不变。
这里有一个while，那么这个算法的复杂度是多少呢？事实上，因为int的计算要比double快很多，尤其是除法，所以这个算法消耗的时间和前两个差不多。实际上多了一个进行n个数取最大公约数的过程。bfor和while是对n个数辗转相除的过程，复杂度\(O(n+log_n)\)，afor和cfor是\(O(n^2)\) 所以总体的复杂度为\(O(n^2*(n+log_n))=O(n^3+n^2log_n)=O(n^3)\)
1.3.4例题
https://www.luogu.com.cn/problem/P3389
别看是个蓝题，其实很简单，像我这种蒟蒻都会。
对于一个n阶行列式来说，这个数据多了一列，但是没关系，我们就带着这一列去计算，我们最终想要的行列式应该是这样的：
你会发现\(a_{nn+1}\)为数b，\(a_{nn}\)是未知数\(x_n\)的系数，由此我们可以算出\(x_n\)的值，再把它的值往回代，算出其它未知数的值即可得出答案。
代码：
#include<iostream>
#include<cstdio>
#include<cmath>
#include<algorithm>
#include<cstring>
#include<sstream>
#include<queue>
#include<map>
#include<vector>
#include<set>
#include<deque>
#include<cstdlib>
#include<ctime>
#define dd double
#define ll long long
#define ull unsigned long long
#define N 102
#define M number
using namespace std;
int n;
dd z[N][N],ans[N];
inline void print(){
for(int i=1;i<=n;i++){
for(int j=1;j<=n+1;j++) printf("%0.2lf ",z[i][j]);
printf("\n");
}
printf("\n");
}
inline bool guass(){
for(int a=1;a<=n;a++){
for(int b=a+1;b<=n;b++)
if(fabs(z[b][a])>fabs(z[a][a]))
for(int c=1;c<=n+1;c++) swap(z[a][c],z[b][c]);
if(fabs(z[a][a])<=1e-8) return 0;
for(int b=a+1;b<=n;b++){
dd k=z[b][a]/z[a][a];
for(int c=1;c<=n+1;c++) z[b][c]-=z[a][c]*k;
}
// print();
}
return 1;
}
int main(){
scanf("%d",&n);
for(int i=1;i<=n;i++)
for(int j=1;j<=n+1;j++)
scanf("%lf",&z[i][j]);
if(!guass()){
printf("No Solution");
return 0;
}
for(int i=n;i>=1;i--){//n-i ans i+1->n
dd x=0.0;
for(int j=i+1;j<=n;j++){
x+=ans[j]*z[i][j];
}
ans[i]=(z[i][n+1]-x)/z[i][i];
}
for(int i=1;i<=n;i++){
printf("%0.2lf\n",ans[i]);
}
return 0;
}
2矩阵
在数学中，矩阵（Matrix）是一个按照长方阵列排列的复数或实数集合。（n*m）
aaaa还有矩阵，真是不想写。。。
2.1特殊矩阵
1.零矩阵:所有元素均为0
2.对角矩阵:除了主对角线，其余元素均为0（n*n）
3.单位矩阵:主对角线元素均为1的对角矩阵。
2.2普通计算与关系
相等：每个位置的元素对应相等，大小相等。
加法：大小相同的两个矩阵，对应位置相加得到的新矩阵。
减法：由加法的逆运算可得到。
数量乘法：一个数乘一个矩阵等于这个数乘上这个矩阵的每一个元素所形成的的新矩阵。
2.3矩阵乘法：
若有矩阵A、B、C，由AB=C，则A的大小为nm，B的大小为mk，C的大小就为nk，也就是说，两个能相乘的矩阵的大小关系是有关联的，这里的这个例子可以很明显的表示出来。
如何计算呢，这里给出计算规则：\({C_{ij}}={\sum}A_{ik}*B_{kj}{(1{\le}k{\le}m)}\)
很难用语言表达出来，\(C_{ij}\)这个位置等于A矩阵中第i行与B中第j列对应元素相乘再相加得到。
2.4矩阵乘法用途：
2.4.1矩阵加速递推。
例题https://www.luogu.com.cn/problem/P1939
原理：我们首先要构造出一个矩阵。1.确定矩阵大小。2.设未知数列方程。3.根据递推式求解方程。4.矩阵快速幂求解答案
这是做这个题的基本思路。
代码：
#include<iostream>
#include<cstdio>
#include<cmath>
#include<algorithm>
#include<cstring>
#include<sstream>
#include<queue>
#include<map>
#include<vector>
#include<set>
#include<deque>
#include<cstdlib>
#include<ctime>
#define dd double
#define ll long long
#define ull unsigned long long
#define N 4
#define M number
using namespace std;
const ll mod=1e9+7;
int t;
struct matrix{
int n,m;
ll a[N][N];
inline matrix(){
n=m=0;
memset(a,0,sizeof(a));
}
inline void dwh(){
m=n;
for(int i=1;i<=n;i++) a[i][i]=1;
}
inline void out(){
for(int i=1;i<=n;i++){
for(int j=1;j<=m;j++) printf("%d ",a[i][j]);
printf("\n");
}
}
};
inline matrix operator * (const matrix &a,const matrix &b){
matrix c;
c.n=a.n;c.m=a.m;
for(int i=1;i<=c.n;i++)
for(int j=1;j<=c.m;j++)
for(int k=1;k<=a.m;k++){
c.a[i][j]+=(a.a[i][k]*b.a[k][j])%mod;
c.a[i][j]%=mod;
}
return c;
}
inline void ksm(matrix &a,ll b){
matrix c;
c.n=a.n;
c.dwh();
while(b){
if(b&1) c=c*a;
// c.out();
a=a*a;
b>>=1;
}
a=c;
}
int main(){
scanf("%d",&t);
while(t--){
matrix a;a.n=1;a.m=3;a.a[1][1]=a.a[1][2]=a.a[1][3]=1;
matrix b;b.n=b.m=3;b.a[1][1]=b.a[1][2]=b.a[2][3]=b.a[3][1]=1;
ll n;scanf("%lld",&n);
ksm(b,n-1);
b=b*a;
printf("%d\n",b.a[1][1]);
}
return 0;
}
相同思路的题还有：https://www.luogu.com.cn/problem/P1962
代码：
#include<iostream>
#include<cstdio>
#include<cmath>
#include<algorithm>
#include<cstring>
#include<sstream>
#include<queue>
#include<map>
#include<vector>
#include<set>
#include<deque>
#include<cstdlib>
#include<ctime>
#define dd double
#define ll long long
#define ull unsigned long long
#define N 3
#define M number
using namespace std;
const ll mod=1e9+7;
struct matrix{
int n,m;
ll a[N][N];
inline matrix(){
n=m=0;
memset(a,0,sizeof(a));
}
inline void dwh(const matrix b){
m=n=b.n;
for(int i=1;i<=n;i++) a[i][i]=1;
}
inline void out(){
for(int i=1;i<=n;i++){
for(int j=1;j<=m;j++) printf("%d ",a[i][j]);
printf("\n");
}
}
};
inline matrix operator * (const matrix &a,const matrix &b){
matrix c;
c.n=a.n;c.m=b.m;
for(int i=1;i<=c.n;i++)
for(int j=1;j<=c.m;j++)
for(int k=1;k<=a.m;k++){
c.a[i][j]+=(a.a[i][k]*b.a[k][j])%mod;
c.a[i][j]%=mod;
}
return c;
}
inline void ksm(matrix &a,ll b){
matrix c;c.dwh(a);
while(b){
if(b&1) c=c*a;
// c.out();
a=a*a;
b>>=1;
}
a=c;
}
int main(){
ll n;scanf("%lld",&n);
if(n==1||n==2){
printf("1");
return 0;
}
matrix a;
a.n=1;a.m=2;a.a[1][1]=a.a[1][2]=1;
matrix b;
b.n=b.m=2;
b.a[1][2]=b.a[2][1]=b.a[2][2]=1;
ksm(b,n-1);
// b.out();
b=a*b;
printf("%lld",b.a[1][1]);
return 0;
}
2.4.2矩阵加速DP
凡是符合\(f_{ij}={\sum}f_{i-1k}*M_{kj}\)形式的DP方程都可以用矩阵优化。
把\(f_i\)看做是一个1行m列的矩阵，则有f\(f_i=f{i-1}*M\) 故可以用矩阵优化。
