luoguP4783 [模板]矩阵求逆 线性代数
求\(n^2\)的矩阵的逆
翻了翻题解,看到了初等矩阵这个东西,突然想起来在看线代的时候看到过....
然后又温习了一遍线性代数的知识
不妨设\(PA = E\),其中\(P\)是一堆初等矩阵的积(必须同时是行变换)
由于\(PA = E, PE = P\),因此\(P(A, E) = (E, P)\)
所以我们只要对矩阵\((A, E)\)来做初等变换
由于我们只做行变换
因此,两个分块矩阵之间互相不干扰
所以当左侧的\(A\)变化为\(E\)时,右边的\(E\)自然变成了\(P\)
复杂度\(O(n^3)\)
#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>
using namespace std;
#define ri register int
#define rep(io, st, ed) for(ri io = st; io <= ed; io ++)
#define drep(io, ed, st) for(ri io = ed; io >= st; io --)
#define gc getchar
inline int read() {
int p = 0, w = 1; char c = gc();
while(c > '9' || c < '0') { if(c == '-') w = -1; c = gc(); }
while(c >= '0' && c <= '9') p = p * 10 + c - '0', c = gc();
return p * w;
}
const int sid = 405;
const int mod = 1e9 + 7;
inline void inc(int &a, int b) { a += b; if(a >= mod) a -= mod; }
inline void dec(int &a, int b) { a -= b; if(a < 0) a += mod; }
inline int mul(int a, int b) { return 1ll * a * b % mod; }
inline int inv(int a) {
int ret = 1, k = mod - 2;
for( ; k; k >>= 1, a = mul(a, a))
if(k & 1) ret = mul(ret, a);
return ret;
}
int n;
int A[sid][sid], B[sid][sid];
inline int Guass() {
rep(i, 1, n) {
int pos = i;
rep(j, i + 1, n) if(A[j][i]) pos = j;
if(!A[pos][i]) return 0;
swap(A[i], A[pos]); swap(B[i], B[pos]);
int IA = inv(A[i][i]);
rep(j, 1, n) {
if(i == j) continue;
int ia = mul(A[j][i], IA);
rep(k, 1, n) {
if(k >= i) dec(A[j][k], mul(ia, A[i][k]));
dec(B[j][k], mul(ia, B[i][k]));
}
}
}
rep(i, 1, n) {
int IA = inv(A[i][i]);
rep(j, 1, n) B[i][j] = mul(B[i][j], IA);
}
return 1;
}
int main() {
n = read();
rep(i, 1, n) rep(j, 1, n) A[i][j] = read();
rep(i, 1, n) B[i][i] = 1;
if(Guass()) {
rep(i, 1, n) {
rep(j, 1, n) printf("%d ", B[i][j]);
printf("\n");
}
}
else printf("No Solution\n");
return 0;
}
也许下次我们可以出一道求\(AP = B\)或者\(PA = B\)的\(P\)
相信能卡死一片人QAQ
喵喵喵?喵喵喵! 喵喵喵......