矩阵求逆和整数逆元比较类似,是通过高斯消元的思路实现的。
思路
求 $n$ 阶矩阵 $A$ 的逆矩阵,也就是 $A^{-1}$,$\frac{1}{A}$,求出一个 $n$ 阶矩阵 $B$ 使得 $AB=I$(其中 $I$ 表示单位矩阵)。
先回顾一下矩阵的三种初等行变换(具体内容详见「算法笔记」矩阵):
- 交换某两行。
- 将某一行的所有元素乘上 $k$(其中 $k\neq 0$)。
- 将某一行的所有元素乘上 $k$ 后加到另一行上去。
假设我们通过 $P_1,P_2,P_3,\dots,P_k$ 等初等矩阵分别左乘 $A$,即 $P_1P_2P_3\cdots P_kA=PA$,其中的 $P$ 为 $P_1$ 到 $P_k$ 的乘积。我们想让 $A$ 变为 $P$ 非常简单,只需要使用高斯消元先变成上三角矩阵,然后变成对角矩阵。
考虑如何求出这个 $P$,我们首先有 $PA=I$,$PI=P$,如果我们同时维护两个矩阵 $A$ 和 $B$ 并使得 $B$ 一开始等于 $I$,在对 $A$ 进行初等行变换变为 $I$ 的过程中也对 $B$ 进行相同的初等行变换。由于初等行变换等价于被对应的初等矩阵左乘,所以当 $A$ 变成 $P$ 后,$B$ 也就从 $I$ 变成了 $P$。
对于无解情况,如果在高斯消元的过程中发现 $A$ 无法变成 $I$ 就说明无解。
时间复杂度:$\mathcal O(n^3)$
代码
#include <cstdio>
#include <cmath>
#include <algorithm>
const int N = 305;
const int mod = 1e9 + 7;
int n, a[N][N << 1];
int pow(int x, int p) {
int ans = 1;
for (; p; p >>= 1, x = 1LL * x * x % mod) {
if (p & 1) ans = 1LL * ans * x % mod;
}
return ans;
}
bool Gauss(int n, int m) {
for (int i = 1; i <= n; i++) {
int p = i;
for (int k = i + 1; k <= n; k++) {
if (std::abs(a[k][i]) > std::abs(a[p][i])) p = k;
}
if (i != p) std::swap(a[i], a[p]);
if (!a[i][i]) return 0;
int d = pow(a[i][i], mod - 2);
for (int j = i; j <= m; j++) {
a[i][j] = 1LL * a[i][j] * d % mod;
}
for (int k = 1; k <= n; k++) {
if (i == k) continue;
int d = a[k][i];
for (int j = i; j <= m; j++) {
a[k][j] = (a[k][j] - 1LL * a[i][j] * d % mod + mod) % mod;
}
}
}
return 1;
}
int main() {
scanf("%d", &n);
int m = n + n;
for (int i = 1; i <= n; i++) {
for (int j = 1; j <= n; j++) {
scanf("%d", &a[i][j]);
}
a[i][n + i] = 1;
}
if (!Gauss(n, m)) {
puts("No Solution");
return 0;
}
for (int i = 1; i <= n; i++) {
for (int j = n + 1; j <= m; j++) {
printf("%d%c", a[i][j], " \n"[j == m]);
}
}
return 0;
}