• 矩阵求逆·学习笔记 $ imes$ [$LuoguP4783$]矩阵求逆


    哦?今天在(luogu)上fa♂现了矩阵求逆的板子……于是就切了切。

    那么我们考虑一个矩阵(A),它的逆矩阵记作(A^{-1}),其中对于矩阵这个群来讲,会有(A cdot A^{-1} = I) 其中(I)表示单位矩阵,主对角线均为(1)

    那么我们对于矩阵(A:)

    (egin{bmatrix} a_{1,1} & a_{1,2} & a_{1,3}\ a_{2,1} & a_{2,2} & a_{2,3} \ a_{3,1} & a_{3,2} & a_{3,3} end{bmatrix} quad)

    我们定义他的单位增广矩阵长底下这个样:

    (egin{bmatrix} a_{1,1} & a_{1,2} & a_{1,3} & | & 1 & 0 & 0\ a_{2,1} & a_{2,2} & a_{2,3} & | & 0 & 1 & 0 \ a_{3,1} & a_{3,2} & a_{3,3} & | & 0 & 0 & 1 end{bmatrix} quad)

    最终这个矩阵应该是(n imes 2m)的。还有,中间那三条线是连起来的!连起来的!但是我并不知道怎么用( m{LaTex})去搞这东西,于是(GG)

    然后呢,我们对这个大矩阵进行暴力消元,直到发现左边的原矩阵变成(I)为止,右边的自然就是(A^{-1}) 那么我们思考它的正确性,以下是我想的根本不严谨只能意会的证明:

    我们思考对于多个基于矩阵初等变换的变换,我们不妨称其为“变换集”。那我们思考对于原矩阵而言,从初态变成(I)的时候,经过的初等变换集虽然不唯一,但是一定组成一个变换集。那么如果右边的单位矩阵的变换集相同的话(意会证明来了),那么也就是说经历了左半个矩阵从原矩阵变成单位矩阵的全过程,那么原来的单位矩阵自然就变为了逆矩阵

    好的,以上都是我不知道多长时间之前学的了(233)。啥?证明太水?没办法啊博主太弱啊……

    那么思考代码:高消

    (hhh)矩阵求逆·完

    好吧,其实有个地方值得一说,就是我们朴素的高消都是消成上三角即可,包括求矩阵的秩,或者消元,或者求行列式。但是此处我们要求消成单位矩阵……嗯,这就很(GG).但是我从(hjc)那里听来一个很牛逼的思路:我们消到最后一定是可以把对角线消成(1),这是很简单的,但,同时如果我们再用(A_{i,i})去消同一列的其他元素,那么就一定可以将其消成(0),因为此时窝萌的(A_{i,i})(1)就可以为所欲为了!那么也就是先正着消一遍,再反着消,用(1)元素去消其他的,然后就(over)了。那么大概效果图如下:

    此处我们用(F**k by A_{i,i})表示被(1)元素(A_{i,i})消成(0),那么

    (egin{bmatrix} 1 & F**k by A_{2,4} & F**k by A_{3,3} & F**k by A_{4,4} \ 0 & 1 & F**k by A_{3,3} &F**k by A_{4,4} \ 0 & 0& 1 & F**k by A_{4,4} \ 0 & 0 & 0 & 1end{bmatrix})

    (233)此处的(F**k)是一个你不可能见过的英文单词!不可能!!

    上代码吧,此处的是(LuoguP4783)

    #include <cstdio>
    #include <iostream>
    #define MAXN 401
    #define ll long long
    #define mod 1000000007
    
    using namespace std ; 
    ll NN, Inv ; bool flag ;
    ll N, A[MAXN][MAXN << 1], i, j, k ;
    
    inline ll qr(){
    	ll k = 0, f = 1 ; char c = getchar() ;
    	while (!isdigit(c)) {if (c == '-') f = -1 ; c = getchar() ;}
    	while (isdigit(c)) k = (k << 3) + (k << 1) + c - 48, c = getchar() ;
    	return k * f ; 
    }
    inline ll expow(ll P, ll Q){
    	ll res = 1 ;
    	while(Q){
    		if (Q & 1) res = res * P % mod ; 
    		P = P * P % mod ;
    		Q >>= 1 ;
    	}  return res ; 
    }
    inline void Gauss(){
    	for (i = 1 ; i <= N ; ++ i){
    		if (A[i][i] == 0) {flag = 1 ; return ;}
    		Inv = expow(A[i][i], mod - 2) , A[i][i] = 1 ;
    		for (j = i + 1 ; j <= NN ; ++ j) A[i][j] = A[i][j] * Inv % mod ;
    		for (j = 1 ; j <= N ; ++ j)
    			if (i == j) continue ;
    			else {
    				Inv = A[j][i], A[j][i] = 0 ;
    				for (k = i + 1 ; k <= NN ; ++ k){
    					A[j][k] = (A[j][k] - Inv * A[i][k]) % mod ;
    					if (A[j][k] < 0) A[j][k] += mod ;
    				}
    			}
    	}
    }
    int main(){
    	cin >> N ; NN = N << 1 ; 
    	for (i = 1; i <= N ; A[i][i + N] = 1, ++ i) 
    		for (j = 1; j <= N ; ++ j) A[i][j] = qr() ;
    	Gauss() ;
    	if (flag) {
    		cout << "No Solution" << endl ;
    		return 0 ;
    	}
    	for (i = 1; i <= N ; ++ i){
    		for(j = N + 1 ; j <= NN ; ++ j)
    			printf("%lld ", A[i][j]) ;
    		putchar('
    ') ;
    	}
    	return 0 ;
    } 
    

    (Upd)出锅了出锅了哈

    上面的代码忘记(swap)了,所以很(GG)……

    // luogu-judger-enable-o2
    #include <cstdio>
    #include <iostream>
    #define MAXN 401
    #define ll long long
    #define mod 1000000007
    
    using namespace std ; 
    ll NN, Inv ; bool flag ;
    ll N, A[MAXN][MAXN << 1], i, j, k ;
    
    inline ll qr(){
        ll k = 0, f = 1 ; char c = getchar() ;
        while (!isdigit(c)) {if (c == '-') f = -1 ; c = getchar() ;}
        while (isdigit(c)) k = (k << 3) + (k << 1) + c - 48, c = getchar() ;
        return k * f ; 
    }
    inline ll expow(ll P, ll Q){
        ll res = 1 ;
        while(Q){
            if (Q & 1) res = res * P % mod ; 
            P = P * P % mod ;
            Q >>= 1 ;
        }  return res ; 
    }
    inline void Gauss(){
        for (i = 1 ; i <= N ; ++ i){
            for(j = i ; j <= N ; ++ j)
            if(A[j][i]){
                for(k = 1 ; k <= NN ; ++ k) swap(A[i][k], A[j][k]);  break ;
            }
            if (A[i][i] == 0) {flag = 1 ; return ;}
            Inv = expow(A[i][i], mod - 2) , A[i][i] = 1 ;
            for (j = i + 1 ; j <= NN ; ++ j) A[i][j] = A[i][j] * Inv % mod ;
            for (j = 1 ; j <= N ; ++ j)
                if (i == j) continue ;
                else {
                    Inv = A[j][i], A[j][i] = 0 ;
                    for (k = i + 1 ; k <= NN ; ++ k)
                        A[j][k] = ((A[j][k] - Inv * A[i][k]) % mod + mod) % mod ;
                }
        }
    }
    int main(){
        cin >> N ; NN = N << 1 ;
        for (i = 1; i <= N ; ++ i) 
            for (j = 1; j <= N ; ++ j)
                A[i][j] = qr() ;
        for (i = 1; i <= N ; ++ i) A[i][i + N] = 1 ;
        Gauss() ;
        if (flag) {
            cout << "No Solution" << endl ;
            return 0 ;
        }
        for (i = 1; i <= N ; ++ i){
            for(j = N + 1 ; j <= NN ; ++ j)
                printf("%lld ", A[i][j]) ;
            putchar('
    ') ;
        }
    } 
    
  • 相关阅读:
    网址收藏
    Linux创建swap文件
    vim命令大全
    char * 和字符数组
    JSR 203终于要出来啦
    对象关系技术的探讨
    最近编码更流畅了
    孤独终止的地方,就是广场开始的地方......
    不要奢望.NET能够跨平台
    实现了HTTP多线程下载
  • 原文地址:https://www.cnblogs.com/pks-t/p/9508081.html
Copyright © 2020-2023  润新知