HDU 5755 Gambler Bo(高斯消元解同余方程)

http://acm.hdu.edu.cn/showproblem.php?pid=5755
一眼高斯消元的水题啊，只是数据是900，900^3的高斯消元可能过不去，所以比赛的时候怂了没有写，早知道就试一试了毕竟也就10分钟写完了。

---

高斯消元求解同余方程，模3的，然后每一位对自己的贡献是2，对四周是1，然后900的数据范围感觉是过不去的，因为高斯消元是$O(n^3)$，but，各位注意到一个条件没有，900个元素中，只有5个不是0啊，所以每次操作的时候，其实只有900 * (895 + 5 * 900)的运算量，如果是0的话直接continue了啊。还是太年轻了啊。

---

代码：

#include <map>
#include <set>
#include <stack>
#include <queue>
#include <cmath>
#include <string>
#include <vector>
#include <cstdio>
#include <cctype>
#include <cstring>
#include <sstream>
#include <cstdlib>
#include <iostream>
#include <algorithm>
#pragma comment(linker,"/STACK:102400000,102400000")

using namespace std;
#define   MAX           905
#define   MAXN          1000005
#define   maxnode       15
#define   sigma_size    30
#define   lson          l,m,rt<<1
#define   rson          m+1,r,rt<<1|1
#define   lrt           rt<<1
#define   rrt           rt<<1|1
#define   middle        int m=(r+l)>>1
#define   LL            long long
#define   ull           unsigned long long
#define   mem(x,v)      memset(x,v,sizeof(x))
#define   lowbit(x)     (x&-x)
#define   pii           pair<int,int>
#define   bits(a)       __builtin_popcount(a)
#define   mk            make_pair
#define   limit         10000

//const int    prime = 999983;
const int    INF   = 0x3f3f3f3f;
const LL     INFF  = 0x3f3f;
const double pi    = acos(-1.0);
const double inf   = 1e18;
const double eps   = 1e-4;
//const LL    mod    = 1e9+7;
const ull    mx    = 133333331;

/*****************************************************/
inline void RI(int &x) {
      char c;
      while((c=getchar())<'0' || c>'9');
      x=c-'0';
      while((c=getchar())>='0' && c<='9') x=(x<<3)+(x<<1)+c-'0';
 }
/*****************************************************/

const int mod = 3;
int exgcd(int a,int b,int &x,int &y){
    if(!b){x = 1; y = 0; return a;}
    else{
        int r = exgcd(b,a%b,y,x);
        y -= x * (a/b);
        return r;
    }
}
int lcm(int a,int b){
    int x = 0, y =0;
    return a / exgcd(a,b,x,y) * b;
}
int A[MAX][MAX],free_x[MAX],x[MAX];
void Gauss(int n,int m){
    int r,c;
    for(r=0,c=0;r<n && c<m;c++){
        int maxr = r;
        for(int i=r+1;i<n;i++) if(abs(A[i][c]) > abs(A[maxr][c])) maxr = i;
        if(maxr != r) for(int i=c;i<=m;i++) swap(A[r][i],A[maxr][i]);
        if(!A[r][c]) continue;
        for(int i=r+1;i<n;i++) if(A[i][c]){
            ///这里要保证运算都是整数，所以要求最小公倍数
            ///范围不大可以用int，int运算更快
            int d = lcm(A[i][c],A[r][c]);
            int t1 = d / A[i][c], t2 = d / A[r][c];
            for(int j=c;j<=m;j++)
                A[i][j] = ((A[i][j] * t1 - A[r][j] * t2) % mod + mod) % mod;
        }
        r++;
    }
    for(int i=r;i<n;i++) if(A[i][m]) return ;
    ///这里保证是没有自由变元的情况下。
    ///有自由变元的时候，不一定是x[i] 对应 A[i][m]应该找那一行最开始的一列不为0的那个
    for(int i=r-1;i>=0;i--){
        x[i] = A[i][m];
        for(int j=i+1;j<m;j++){
            x[i] = ((x[i] - A[i][j] * x[j]) % mod + mod) % mod;
        }
        int x1 = 0,y1 = 0;
        ///这里是用exgcd求逆元，也可以用费马小定理求，如果mod是素数
        int d = exgcd(A[i][i],mod,x1,y1);
        x1 = ((x1 % mod) + mod) % mod;
        x[i] = x[i] * x1 % mod;
    }
}
void Gauss_init(){
    mem(A,0); mem(free_x,0); mem(x,0);
}
int dx[4]={0,1,0,-1};
int dy[4]={1,0,-1,0};
int a[34][34];
int main(){
    //freopen("in.txt","r",stdin);
    int t;
    cin>>t;
    while(t--){
        int n,m;
        cin>>n>>m;
        for(int i=0;i<n;i++){
            for(int j=0;j<m;j++) scanf("%d",&a[i][j]);
        }
        Gauss_init();
        for(int i=0;i<n*m;i++){
            int r=i/m;
            int c=i%m;
            A[i][i]=2;
            for(int j=0;j<4;j++){
                int x=r+dx[j];
                int y=c+dy[j];
                if(x>=0&&y>=0&&x<n&&y<m){
                    A[i][x*m+y]=1;
                }
            }
            A[i][n*m]=(3-a[r][c])%3;
        }
        Gauss(n*m,n*m);
        vector<int> v;
        for(int i=0;i<n*m;i++){
            while(x[i]){
                v.push_back(i);
                x[i]--;
            }
        }
        printf("%d\n",v.size());
        for(int i=0;i<v.size();i++){
            int r=v[i]/m;
            int c=v[i]%m;
            r++,c++;
            printf("%d %d\n",r,c);
        }
    }
    return 0;
}