高斯消元模板（转）

参考资料：http://blog.sina.com.cn/s/blog_afe6bb330101a59d.html

模板：

const int maxn=30;
int a[maxn][maxn+1],x[maxn];//a是系数矩阵和增广矩阵，x是最后存放的解

// a[][maxn]中存放的是方程右面的值（bi）

int equ,var;//equ是系数阵的行数，var是系数矩阵的列数（变量的个数）

int free_num,ans=100000000;

void Debug(void)  //调试输出，看消元后的矩阵值，提交时，就不用了
{
    int i, j;
    for (i = 0; i < equ; i++)
    {
        for (j = 0; j < var + 1; j++)
        {
            cout << a[i][j] << " ";
        }
        cout << endl;
    }
    cout << endl;
}
inline int gcd(int a, int b) //最大公约数
{
    int t;
    while (b != 0)
    {
        t = b;
        b = a % b;
        a = t;
    }
    return a;
}
inline int lcm(int a, int b)  //最小公倍数
{
    return a  / gcd(a, b)*b;
}

void swap(int &a,int &b)
{
    int temp=a;    //交换2个数
    a=b;
    b=temp;
}

int Gauss()
{
    int k,col = 0;  //当前处理的列
    for(k = 0; k < equ && col < var; ++k,++col)
    {
        int max_r = k;
        for(int i = k+1; i < equ; ++i)
            if(a[i][col] > a[max_r][col])
                max_r = i;
        if(max_r != k)
        {
            for(int i = k; i < var + 1; ++i)
                swap(a[k][i],a[max_r][i]);
        }
        if(a[k][col] == 0)
        {
            k--;
            continue;
        }
        for(int i = k+1; i < equ; ++i)
        {
            if(a[i][col] != 0)
            {
                int LCM = lcm(a[i][col],a[k][col]);
                int ta = LCM/a[i][col], tb = LCM/a[k][col];
                if(a[i][col]*a[k][col] < 0)
                    tb = -tb;
                for(int j = col; j < var + 1; ++j)
                    a[i][j] = ( (a[i][j]*ta)%2 - (a[k][j]*tb)%2  +  2 ) % 2;    //a[i][j]只有0和1两种状态
            }
        }
    }
    //上述代码是消元的过程，行消元完成
    //解下来2行，判断是否无解
    //注意K的值，k代表系数矩阵值都为0的那些行的第1行
    for(int i = k; i < equ; ++i)
        if(a[i][col] != 0)
            return -1;   // 无解返回 -1
    //Debug();
    //唯一解或者无穷解,k<=var
    //var-k==0  唯一解；var-k>0 无穷多解，自由解的个数=var-k
    //能执行到这，说明肯定有解了，无非是1个和无穷解的问题。
    //下面这几行很重要，保证秩内每行主元非0，且按对角线顺序排列，就是检查列
    for(int i = 0; i <equ; i++)
        if(!a[i][i])
        {
            int j;
            for(j=i+1;j<=var;j++)
                if(a[i][j])
                    break;
            if(j == var)
                break;
            for(int k = 0; k < equ; ++k)
                swap(a[k][i],a[k][j]);
        }
    // ----处理保证对角线主元非0且顺序，检查列完成
    // free_num=k;
    if (var-k>0)
    {
     //无穷多解，先枚举解，然后用下面的回带代码进行回带；
     //这里省略了下面的回带的代码；不管唯一解和无穷解都可以回带，只不过无穷解//回带时，默认为最后几个自由变元=0而已。
    }
    if (var-k==0)//唯一解时
    {
        //下面是回带求解代码，当无穷多解时，最后几行为0的解默认为0；
        for(int i = k-1; i >= 0; --i) //从消完元矩阵的主对角线非0的最后1行，开始往//回带
        {
            int tmp = a[i][var] % 2;
            for(int j = i+1; j < var; ++j) //x[i]取决于x[i+1]--x[var]啊，所以后面的解对前面的解有影响。
                if(a[i][j] != 0)
                    tmp = ( tmp - (a[i][j]*x[j])%2  +  2 ) % 2;
            //if (a[i][i]==0) x[i]=tmp;   //最后的空行时，即无穷解得
            //else
            x[i] = (tmp/a[i][i]) % 2; //上面的正常解
        }
        //回带结束了
    }
}