高斯消元
高斯消元法（Gauss-Jordan elimination）是求解线性方程组的经典算法，它在当代数学中有着重要的地位和价值，是线性代数课程教学的重要组成部分。
高斯消元法除了用于线性方程组求解外，还可以用于行列式计算、求矩阵的逆，以及其他计算机和工程方面。
夏建明等人之前提出了应用图形处理器 (GPU) 加速求解线性方程组的高斯消元法，所提出的算法与基于 CPU 的算法相比较取得更快的运算速度。二是提出各种变异高斯消元法以满足特定工作的需要。
高斯消元可以用于求解多元线性方程组，形如：
实现
将各项系数转化为矩阵。
每列都进行消元后，每列中只有一行该列有值，其余行该列全是零。
什么叫消元？例如：
现在用它来消去其余行的 x_1 系数，构造：
然后对于任意的：
都可以做对应的减法，得到：
(a_{i,1} - a_{i,1} times 1)x_1 (a_{i,2} - a_{i,1} times dfrac{a_{1,2}}{a_{1,1}})x_2 ldots = b_i - a_{i,1} times dfrac{b_1}{a_{1,1}}
这样消去后，该行 x_1 系数为 0.
使用每行的方程都消去一元，则每行最终系数都为：
0x_1 0x_2 ldots kx_i ldots 0x_n = b'
即可解出对应的 x_i.
在具体实现时，无须将用于消去的行的对应系数设为 1，只需要计算出消去相应所需的倍数即可。
a_{1,1}x_1 ldots
a_{2,1}x_1 ldots
只需要令 t = dfrac{a_{2,1}}{a_{1,1}}，而减去一式的 t 倍即可。
- 此外，若某列的系数全为 0，即该列无限制，说明方程组无唯一解。
- 消去 x_i 时，可以选择 a_{j,i} 最大的，即该列对应系数最大的行，作为用于消去的行，可以提高精度。
时间复杂度 O(n^3).
代码语言：javascript复制#include <algorithm>
#include <cmath>
#include <cstdio>
const int maxn = 110;
int n;
double a[maxn][maxn];
const double eps = 1e-8;
int main()
{
scanf("%d", &n);
for (int i = 1; i <= n; i) for (int j = 1; j <= n 1; j) scanf("%lf", &a[i][j]);
for (int i = 1; i <= n; i) // 选择列
{
//寻找系数最大行
int maxx = i;
for (int j = i 1; j <= n; j) if (fabs(a[j][i]) > fabs(a[maxx][i])) maxx = j;
//将系数最大行交换至第 i 行
for (int j = 1; j <= n 1; j) std::swap(a[i][j], a[maxx][j]);
//系数全零
if (fabs(a[i][i]) < eps) goto fail;
for (int j = 1; j <= n; j)
{
if (j == i) continue;
//消元
double t = a[j][i] / a[i][i];
for (int k = i 1; k <= n 1; k) a[j][k] -= a[i][k] * t;
}
}
//解一元一次方程
for (int i = 1; i <= n; i) printf("%.2fn", a[i][n 1] / a[i][i]);
return 0;
fail:
puts("No Solution");
return 0;
}矩阵求逆
高斯消元法可以用于矩阵求逆。
什么叫矩阵的逆？对于矩阵 A,B，满足 A times B = B times A = I_n，其中 I_n 为 n 阶单位矩阵，则 B 为 A 的逆矩阵。
形象地说：
有：
A times B = begin{bmatrix} 1 & 0 & cdots & 0 cdots & cdots & cdots & cdots 0 & 0 & cdots & 1 end{bmatrix}
则 B 为 A 的逆矩阵。
矩阵求逆的做法：
- 将 A 与 I 放在同一个矩阵中
- 对 A 进行消元，将 A 化为单位矩阵
- 此时原单位矩阵转化为 A 的逆矩阵
可以发现，高斯消元后，原矩阵化为一个对角矩阵，即只有 a_{i,i} 上有值。
而每行除以 a_{i,i} 后，即可化作单位矩阵。
而 A 通过一系列的线性变化变成了 I，即单位矩阵，每次变化都可以写成矩阵的形式。
可以这样表达：A times P_1 times P_2 times ldots P_m = I
即 A times prod limits _{i=1}^m P_i = I
根据逆矩阵定义，A times B = I，所求即 prod limits_{i=1}^m P_i.
而将单位矩阵放在 A 的右边，每次变化都会对单位矩阵同样变化。
A times B = I
I times B = B
即可求出 A 的逆矩阵 B.
代码语言：javascript复制#include <algorithm>
#include <cstdio>
const int maxn = 810;
const int p = 1e9 7;
#define ll long long
ll a[maxn][maxn];
int n;
void exgcd(ll a,ll b,ll &x,ll &y)
{
if (b == 0) return (void)(x = 1, y = 0);
exgcd(b, a % b, y, x), y -= (a / b) * x;
}
ll inv(ll x)
{
ll k1,k2;
exgcd(x,p,k1,k2);
return ((k1 % p) p) % p;
}
signed main()
{
scanf("%d", &n);
for (int i = 1; i <= n; i) for (int j = 1; j <= n; j) scanf("%lld", &a[i][j]);
for (int i = 1; i <= n; i) a[i][i n] = 1;
for (int i = 1; i <= n; i)
{
int maxx = i;
for (int j = i 1; j <= n; j) if (a[j][i] > a[maxx][i]) maxx = j;
if (i != maxx) std::swap(a[i], a[maxx]);
if (a[i][i] == 0) goto fail;
ll tmp = inv(a[i][i]);//求逆元
for (int j = 1; j <= n; j)
{
if (i == j) continue;
ll frac = (a[j][i] * tmp) % p; // a[j][i] / a[i][i]
for (int k = i; k <= 2 * n; k)
a[j][k] = (((a[j][k] - frac * a[i][k]) % p) p) % p;//消元
}
for (int j = 1; j <= 2 * n; j) a[i][j] = (a[i][j] * tmp) % p;//处理为 1
}
for (int i = 1; i <= n; puts(""), i) for (int j = 1; j <= n; j) printf("%lld ", a[i][j n]);
return 0;
fail:
puts("No Solution");
return 0;
}
