[高斯消元及理论]线性方程组整数/浮点数，模线性方程组，异或方程组模板

理论

高斯消元法，是线性代数规划中的一个算法，可用来为线性方程组求解。但其算法十分复杂，不常用于加减消元法，求出矩阵的秩，以及求出可逆方阵的逆矩阵

用系数矩阵表示n元一次方程，如
$$\{\begin{array}{l}2x+3y+4z=20\\ x-2y+3z=6\\ 6x+4y-3z=5\end{array}\Rightarrow \left(\begin{array}{cccc}2& 3& 4& 20\\ 1& -2& 3& 6\\ 6& 4& -3& 5\end{array}\right)$$
系数矩阵的运算法则：
• 1.将某一行同时乘以/除以一个非0常数。
• 2.某一行加上或减去另一行的若干倍

---

演示过程

$$\left(\begin{array}{cccc}2& 3& 4& 20\\ 1& -2& 3& 6\\ 6& 4& -3& 5\end{array}\right){\Rightarrow }^{r_1\ast \frac{1}{2}}\left(\begin{array}{cccc}1& \frac{3}{2}& 2& 10\\ 1& -2& 3& 6\\ 6& 4& -3& 5\end{array}\right){\Rightarrow }^{r_2-r_1,r_3-6r_1}\left(\begin{array}{cccc}1& \frac{3}{2}& 2& 10\\ 0& -\frac{7}{2}& 1& -4\\ 0& -5& -15& -55\end{array}\right)$$
$${\Rightarrow }^{r_2\ast (-\frac{2}{7})}\left(\begin{array}{cccc}1& \frac{3}{2}& 2& 10\\ 0& 1& -\frac{2}{7}& \frac{8}{7}\\ 0& -5& -15& -55\end{array}\right){\Rightarrow }^{r_1-\frac{3}{2}{r}_2,r_3+5r_2}\left(\begin{array}{cccc}1& 0& \frac{17}{7}& \frac{58}{7}\\ 0& 1& -\frac{2}{7}& \frac{8}{7}\\ 0& 0& -\frac{115}{7}& -\frac{345}{7}\end{array}\right)$$$${\Rightarrow }^{r_3\ast (-\frac{7}{115})}\left(\begin{array}{cccc}1& 0& \frac{17}{7}& \frac{58}{7}\\ 0& 1& -\frac{2}{7}& \frac{8}{7}\\ 0& 0& 1& 3\end{array}\right){\Rightarrow }^{r_1-\frac{17}{7}{r}_3,r_2+\frac{2}{7}{r}_3}\left(\begin{array}{cccc}1& 0& 0& 1\\ 0& 1& 0& 2\\ 0& 0& 1& 3\end{array}\right)$$
即$$\{\begin{array}{l}{x}_1=1\\ x_2=2\\ x_3=3\end{array}$$

---

算法流程：
• 从1到n枚举i，第i步定未知数i为主元，从上到下找到第一个未知数i不为0的方程定为主方程（如不存在则跳过这一步）
• 通过乘除系数使得主方程未知数i前面的系数为1 • 对于其他方程，若未知数i前面的系数为k，那么该方程减去主方程的k倍，使得该方程未知数i系数为0
• 时间复杂度O(n^3)

---

一般来讲，n个未知数n个方程在方程之间线性不相关的情况下，恰好有一组解
如果算法结束后，存在一行全为0，那么被消去的未知数可以任意取值，称为自由元
如果有一行前n个数为0，最后一个数不为0，那么方程无解
当方程组数和未知数数不相等时，也可通过上面方法判定解数

提一句：多解
解的个数等于自由变量的取值范围的幂
若自由变量的取值范围为$p$ ,自由变量有$m$个，则解的个数为$p^m$

线性方程组整数类型解

int a[maxn][maxn];//增广矩阵 
int ans[maxn];
bool Free[maxn];//标记是否为自由变元 

int GCD ( int a, int b ) {
	if ( ! b )
		return a;
	return GCD ( b, a % b );
}
int LCM ( int a, int b ) {
	return a / GCD ( a, b ) * b;
}
int Fabs ( int x ) {
	if ( x < 0 )
		return -x;
	return x;
}

int Gauss ( int equ, int var ) {
	for ( int i = 0;i <= var;i ++ ) {
		ans[i] = 0;
		Free[i] = 1;
	}
	int row, col, MaxRow;
	col = 1;
	for ( row = 1;row <= equ && col < var;row ++, col ++ ) {
		MaxRow = row;
		for ( int i = row + 1;i <= equ;i ++ )//寻找当前列绝对值最大的行 
			if ( Fabs ( a[i][col] ) > Fabs ( a[MaxRow][col] ) )
				MaxRow = i;
		if ( MaxRow != row ) {//与第row行交换 
			for ( int i = row;i <= var;i ++ )
				swap ( a[row][i], a[MaxRow][i] );
		}
		if ( ! a[row][col] ) {//说明col列第row行及以下都是0，就处理当前行的下一列 
			row --;
			continue;
		}
		for ( int i = row + 1;i <= equ;i ++ ) {//枚举要删去各自col列的变元的系数的行 
			if ( a[i][col] ) {
				int lcm = LCM ( Fabs ( a[i][col] ), Fabs ( a[row][col] ) );
				int T1 = lcm / Fabs ( a[i][col] );
				int T2 = lcm / Fabs ( a[row][col] );
				if ( a[i][col] * a[row][col] < 0 )
					T2 = -T2;
				for ( int j = col;j <= var;j ++ )
					a[i][j] = a[i][j] * T1 - a[row][j] * T2;
			}
		}
	}
	//无解：化简的增广矩阵中存在(0,0,...,a)这种类型 
	for ( int i = row;i <= equ;i ++ )
		if ( a[i][col] )
			return -1;
	int temp;
	//无穷解，矩阵中出现(0,0,...0) 
	if ( row < var ) {
		for ( int i = row - 1;i > 0;i -- ) {
			// 第i行一定不会是(0, 0, ..., 0)的情况，因为这样的行是在第k行到第equ行.
            // 同样，第i行一定不会是(0, 0, ..., a), a != 0的情况，这样的无解的
			int free_num = 0, idx;
			//用于判断该行中的不确定的变元的个数，如果超过1个，则无法求解，它们仍然为不确定的变元
			for ( int j = 1;j < var;j ++ )	
				if ( a[i][j] && Free[j] ) {
					free_num ++;
					idx = j;
				}
			if ( free_num > 1 ) //无法求解出确定的变元 
				continue;
			temp = a[i][var];//说明只有一个不确定的变元，那么就可以解出该变元
			for ( int j = 1;j < var;j ++ ) {
				if ( a[i][j] && j != idx )
					temp -= a[i][j] * ans[j];
			}
			ans[idx] = temp / a[i][idx];
			Free[idx] = 0;
		}
		return var - row;//自由变元的个数 
	}
	//唯一解的情况，是一个严格的上三角形
	//1 0 ... 0
	//0 1 ... 0
	//0 0 1 . 0 
	for ( int i = var - 1;i > 0;i -- ) {
		temp = a[i][var];
		for ( int j = i + 1;j < var;j ++ )
			if ( a[i][j] )
				temp -= a[i][j] * ans[j];//--因为x[i]存的是temp/a[i][i]的值，即是a[i][i]=1时x[i]对应的值
//		if ( temp % a[i][i] )//可以用来判断有浮点数解，无整数解，但这个是整数解模板 
//			return -2;
		ans[i] = temp / a[i][i];
	}
	return 0;
}

线性方程组浮点类型解

#define eps 1e-6
double a[maxn][maxn];
double ans[maxn];
bool Free[maxn];

int GCD ( int a, int b ) {
	if ( ! b )
		return a;
	return GCD ( b, a % b );
}
int LCM ( int a, int b ) {
	return a / GCD ( a, b ) * b;
}

int Gauss ( int equ, int var ) {
	for ( int i = 0;i <= var;i ++ ) {
		ans[i] = 0;
		Free[i] = 1;
	}
	int row, col, MaxRow;
	col = 0;
	for ( row = 0;row < equ && col < var;row ++, col ++ ) {
		MaxRow = row;
		for ( int i = row + 1;i < equ;i ++ ) 
			if ( fabs ( a[i][col] ) > fabs ( a[MaxRow][col] ) )
				MaxRow = i;
		if ( MaxRow != row ) {
			for ( int i = row;i <= var;i ++ )
				swap ( a[row][i], a[MaxRow][i] );
		}
		if ( fabs ( a[row][col] ) < eps ) {
			row --;
			continue;
		}
		for ( int i = row + 1;i < equ;i ++ ) {
			if ( fabs ( a[i][col] ) > eps ) {
				double temp = a[i][col] / a[row][col];
				for ( int j = col;j <= var;j ++ )
					a[i][j] -= a[row][j] * temp;
				a[i][col] = 0;
			}
		}
	}
	for ( int i = row;i < equ;i ++ )
		if ( fabs ( a[i][col] ) > eps )
			return -1;
	double temp;
	if ( row < var ) {
		for ( int i = row - 1;i >= 0;i -- ) {
			int free_num = 0, idx;
			for ( int j = 0;j < var;j ++ )	
				if ( a[i][j] && Free[j] ) {
					free_num ++;
					idx = j;
				}
			if ( free_num > 1 )
				continue;
			temp = a[i][var];
			for ( int j = 0;j < var;j ++ ) {
				if ( a[i][j] && j != idx )
					temp -= a[i][j] * ans[j];
			}
			ans[idx] = temp / a[i][idx];
			Free[idx] = 0;
		}
		return var - row;
	}
	for ( int i = var - 1;i >= 0;i -- ) {
		temp = a[i][var];
		for ( int j = i + 1;j < var;j ++ )
			if ( a[i][j] )
				temp -= a[i][j] * ans[j];
		ans[i] = temp / a[i][i];
	}
	return 0;
}

模线性方程组

int gcd( int x, int y ) {
	if( ! y ) return x;
	else return gcd( y, x % y );
}

int inv( int x, int m ) {
	if( x == 0 ) return 0;
	if( x == 1 ) return 1;
	return inv( m % x, m ) * ( m - m / x ) % m;
}

int Gauss( int equ, int var ) {
	int row, col = 1;
	for( row = 1;row <= equ && col < var;row ++, col ++ ) {
		int MaxRow = row;
		for( int i = row + 1;i <= equ;i ++ )
			if( Fabs( a[i][col] ) > Fabs( a[MaxRow][col] ) )
				MaxRow = i;
		if( row != MaxRow ) {
			for( int i = row;i <= var;i ++ )
				swap( a[row][i], a[MaxRow][i] );
		}
		if( ! a[row][col] ) { row --; continue; }
		for( int i = row + 1;i <= equ;i ++ ) {
			if( a[i][col] ) {
				int d = gcd( Fabs( a[row][col] ), Fabs( a[i][col] ) );
				int lcm = Fabs( a[row][col] ) / d * Fabs( a[i][col] );
				int T1 = lcm / Fabs( a[i][col] );
				int T2 = lcm / Fabs( a[row][col] );
				if( a[i][col] * a[row][col] < 0 ) T2 = -T2;
				for( int j = col;j <= var;j ++ ) {
					a[i][j] = a[i][j] * T1 - a[row][j] * T2;
					a[i][j] = ( a[i][j] % mod + mod ) % mod;
				}
			}
		}
	}
	for( int i = row;i <= equ;i ++ )
		if( a[i][col] ) return -1;
	if( row < var ) return var - row;
	for( int i = var - 1;i;i -- ) {
		int temp = a[i][var];
		for( int j = i + 1;j < var;j ++ ) {
			if( a[i][j] ) {
				temp -= a[i][j] * x[j];
				temp = ( temp % mod + mod ) % mod;
			}
		}
		x[i] = temp * inv( a[i][i], mod ) % mod;
	}
	return 0;
}

异或方程组

int cnt;
int a[maxn][maxn];
int ans[maxn];
bool Free[maxn];

int GCD ( int a, int b ) {
	if ( ! b )
		return a;
	return GCD ( b, a % b );
}
int LCM ( int a, int b ) {
	return a / GCD ( a, b ) * b;
}

int Gauss ( int equ, int var ) {
	int row, col, MaxRow;
	col = 0;
	for ( row = 0;row < equ && col < var;row ++, col ++ ) {
		MaxRow = row;
		for ( int i = row + 1;i < equ;i ++ ) 
			if ( fabs ( a[i][col] ) > fabs ( a[MaxRow][col] ) )
				MaxRow = i;
		if ( MaxRow != row ) {
			for ( int i = row;i <= var;i ++ )
				swap ( a[row][i], a[MaxRow][i] );
		}
		if ( a[row][col] == 0 ) {
			row --;
			Free[++ cnt] = col;
			continue;
		}
		for ( int i = row + 1;i < equ;i ++ ) {
			if ( a[i][col] ) {
				for ( int j = col;j <= var;j ++ )
					a[i][j] ^= a[row][j];
			}
		}
	}
	for ( int i = row;i < equ;i ++ )
		if ( a[i][col] )
			return -1;
	if ( row < var ) 
		return var - row;
	for ( int i = var - 1;i >= 0;i -- ) {
		ans[i] = a[i][var];
		for ( int j = i + 1;j < var;j ++ )
			if ( a[i][j] )
				ans[i] ^= ( a[i][j] && ans[j] );
	}
	return 0;
}

高斯约旦消元

约旦消元

void gauss_jordan ( int equ, int var ) {
	int row, col = 0;
	for ( row = 0;row < equ && col < var;row ++, col ++ ) {
		int MaxRow = row;
		for ( int i = row + 1;i < equ;i ++ )
			if ( fabs ( a[i][col] ) > fabs ( a[MaxRow][col] ) )
				MaxRow = i;
		if ( MaxRow != row ) {
			for ( int i = 0;i <= var;i ++ )
				swap ( a[row][i], a[MaxRow][i] );
		}
		if ( fabs ( a[row][col] ) < eps ) {
			row --;
			continue;
		}
		for ( int i = 0;i < equ;i ++ )
			if ( i != row )
				for ( int j = var;j >= row;j -- )
					a[i][j] -= a[i][col] / a[row][col] * a[row][j];
	}
}

无解

for ( int i = 0;i < equ;i ++ ) {
		bool flag = 1;
		for ( int j = i;j < var;j ++ )
			if ( fabs ( a[i][j] ) > eps )
				flag = 0;
		if ( flag && fabs ( a[i][var] ) > eps )
			return ! printf ( "-1" );
	}

无穷解

	for ( int i = 0;i < equ;i ++ ) {
		int free_num = 0;
		for ( int j = i;j < var;j ++ )
			if ( fabs ( a[i][j] ) > eps )
				free_num ++;
		if ( free_num > 1 )
			return ! printf ( "0" );
	}

唯一解

	for ( int i = 0;i < equ;i ++ )
		ans[i] = a[i][var] / a[i][i];
	for ( int i = 0;i < equ;i ++ )
		if ( fabs ( ans[i] ) < eps )
			printf ( "x%d=0\n", i + 1 );
		else
			printf ( "x%d=%.2f\n", i + 1, ans[i] );

温馨提示，不保证code的绝对正确，反正我是过了
错了别打我，至少别打脸（跑ε=ε=ε=(_￣▽￣)