高斯消元模板

入门去看：英雄哪里出来和诚叙(两位大佬)。

理解后的部分修改的两个板子：整型和浮点型。

//整型版
#include <string.h>
#include <cstdio>
#define ll long long
using namespace std;
struct GAUSS {
const static int mod = 1e9+7;
const static int maxn = 1005;
int r, c, free_x, cnt; //cnt表示x取值的值域
ll matrix[maxn][maxn];
ll x[maxn];
ll ans;
void init(int _r, int _c, int _cnt)
{
r = _r, c = _c, cnt = _cnt;
memset(matrix, 0, sizeof matrix);
}
//构造增广矩阵(系数, b)
void Link(int r, int c, ll val)
{
matrix[r][c] = val;
}
ll _abs(ll x)
{
return x < 0? -x: x;
}
ll qpow(int base, int n)
{
ll ans = 1;
for(int i = 1; i <= n; ++i) ans *= base;
return ans;
}
void ex_gcd(ll a, ll b, ll &d, ll &x, ll &y)
{
if(b == 0)
{
x = 1, y = 0, d = a;
return;
}
ex_gcd(b, a%b, d, y, x);
y -= a/b*x;
}
void swap_row(int a, int b)
{
for(int i = 0; i <= c; ++i)
{
ll tmp = matrix[a][i];
matrix[a][i] = matrix[b][i];
matrix[b][i] = tmp;
}
}
void swap_col(int a, int b)
{
for(int i = 0; i < r; ++i)
{
ll tmp = matrix[i][a];
matrix[i][a] = matrix[i][b];
matrix[i][b] = tmp;
}
}
bool gauss()
{
int row, col, maxrow, i, j;
for(row = 0, col = 0; row < r && col < c; ++row)
{
maxrow = row;
for(i = row+1; i < r; ++i)
{
if(_abs(matrix[i][col]) > _abs(matrix[maxrow][col])) maxrow = i;
}
if(matrix[maxrow][col] == 0)
{
++col;
--row;
continue;
}
if(maxrow != row) swap_row(row, maxrow);
for(i = row+1; i < r; ++i)
{
if(matrix[i][col])
{
ll tmp = matrix[i][col]; //matrix[i][col]必须暂时保存，在首次计算时就会被改变
for(j = col; j <= c; ++j)
{
matrix[i][j] = (matrix[i][j]*matrix[row][col]-matrix[row][j]*tmp)%mod;
if(matrix[i][j] < 0) matrix[i][j] += mod;
}
}
}
++col;
}
for(i = row; i < r; ++i)
if(matrix[i][c]) return false; //无解
free_x = c-row;
//转化为上三角
for(i = 0; i < r && i < c; ++i)
{
if(matrix[i][i] == 0)
{
for(j = i+1; j < c; ++j)//虽然此处条件是判断是否<c，但实际上要么一定能找到系数矩阵中一列，要么连带b向量中全找不到
if(matrix[i][j])
{
swap_col(i, j);
break;
}
}
}
ans = qpow(cnt, free_x);
return true;
}
//如果存在唯一解，下面计算每个未知数的值
void calc_x()
{
for(int i = c-
4000
1; i >= 0; --i)
{
ll tmp = matrix[i][c];
for(int j = i+1; j < c; ++j)
tmp = (tmp-matrix[i][j]*x[j]+mod)%mod;
ll D, X, Y;
ex_gcd(matrix[i][i], mod, D, X, Y);
X = (X%mod+mod)%mod;
x[i] = X*tmp%mod;
}
}
} AC;
void solve()
{
int r, c, cnt;
scanf("%d %d %d", &r, &c, &cnt);
AC.init(r, c, cnt);

for(int i = 0; i < r; ++i)
for(int j = 0; j <= c; ++j) //输入增广矩阵(系数, b)
{
int x;
scanf("%lld", &x);
AC.Link(i, j, x);
}

if(!AC.gauss()) printf("No Answer.\n");
{
printf("%lld\n", AC.ans);
if(AC.ans == 1)
{
AC.calc_x();
for(int i = 0; i < c; ++i)
printf("%lld ", AC.x[i]);
printf("\n");
}
}
}
int main()
{
solve();
return 0;
}

//浮点数版
#include <string.h>
#include <cstdio>
#include <cmath>
#define LL long long
using namespace std;
struct GAUSS {
const static int maxn = 105;
const static double eps = 1e-12;
int r, c, free_x;
double matrix[maxn][maxn];
double x[maxn];
void init(int _r, int _c)
{
r = _r, c = _c;
//浮点数也可以这么清零
memset(matrix, 0, sizeof matrix);
}
//构造增广矩阵(系数, b)
void Link(int r, int c, double val)
{
matrix[r][c] = val;
}
inline bool sgn(double x) //判断绝对值是否小于极小数, = 0代表为极小数
{
return (x > eps) - (x < -eps);
}
void swap_row(int a, int b)
{
for(int i = 0; i <= c; ++i)
{
double tmp = matrix[a][i];
matrix[a][i] = matrix[b][i];
matrix[b][i] = tmp;
}
}
void swap_col(int a, int b)
{
for(int i = 0; i < r; ++i)
{
double tmp = matrix[i][a];
matrix[i][a] = matrix[i][b];
matrix[i][b] = tmp;
}
}
bool gauss()
{
int row, col, maxrow, i, j;
for(row = 0, col = 0; row < r && col < c; ++row)
{
maxrow = row;
for(i = row+1; i < r; ++i) //先进行此处判断，因为绝对值小于极小数的数都被认为是0
{
if(fabs(matrix[i][col]) > fabs(matrix[maxrow][col])) maxrow = i;
}
if(sgn(matrix[maxrow][col]) == 0)
{
++col;
--row;
continue;
}
if(maxrow != row) swap_row(row, maxrow);
for(i = row+1; i < r; ++i)
{
if(sgn(matrix[i][col]) != 0)
{
double tmp = matrix[i][col];
for(j = col; j <= c; ++j)
{
matrix[i][j] -= tmp/matrix[row][col]*matrix[row][j];
}
}
}
++col;
}
for(i = row; i < r; ++i)
if(sgn(matrix[row][c]) != 0) return false; //无解
free_x = c-row;
for(i = 0; i < r && i < c; ++i)
{
if(sgn(matrix[i][i]) == 0)
{
for(j = i+1; j < c; ++j)
if(sgn(matrix[i][j]) != 0)
{
swap_col(i, j);
break;
}
}
}
return true;
}
void calc_x()
{
for(int i = c-1; i >= 0; --i)
{
double tmp = matrix[i][c];
for(int j = i+1; j < c; ++j)
tmp = tmp-matrix[i][j]*x[j];
x[i] = tmp/matrix[i][i];
}
}
} AC;
void solve()
{
int r, c;
scanf("%d %d", &r, &c);
AC.init(r, c);
for(int i = 0; i < r; ++i)
for(int j = 0; j < c; ++j)
{
double x;
scanf("%lf", &x);
AC.Link(i, j, x);
}
AC.gauss();
AC.calc_x();
for(int i = 0; i < c-1; ++i)
printf("%.2f ", AC.x[i]+AC.eps); //避免-0.0的输出
printf("%.2f\n", AC.x[c-1]+AC.eps);
}
int main()
{
solve();
return 0;
}

附：-0.0的产生是因为被用来表示一个可以四舍五入为零的负数，或者是一个从负方向上趋近于零的数。所以浮点数运算时最终结果需要加一个eps避免-0.0。

Go on.