【bzoj2854】高斯消元

求解N元一次方程组，N<=200，每个未知元前系数范围在[0,10^18]（常数不保证），保证答案为[0,10^9]间的整数且唯一解。

考虑大素数取模下做，整数除法用逆元即可，再用中国剩余定理合并，注意选的素数不能使方程多解(矩阵满秩)。

#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
#include <string>
#define Rep(i, x, y) for (int i = x; i <= y; i ++)
using namespace std;
typedef long long LL;
const int N = 205;
int n, mod[4] = {1000000007, 1000000021, 1000000607, 1000000009}, A[N][N];
LL a[N][N], p[N], X[4][N], x1 = -1, x2 = -1, G[4][N];
char s[N];
LL Pow(LL x, LL P) {
	LL z = 1, y = P - 2;
    while (y) {
        if (y % 2 != 0) (z *= x) %= P;
        (x *= x) %= P, y /= 2;
    }
    return z;
}
bool Work(int x) {
    int P = mod[x];
    Rep(i, 1, n) {
        Rep(j, 1, n) a[i][j] = A[i][j];
        p[i] = G[x][i] % P;
    }
    Rep(i, 1, n) {
        int k = i;
        while (k <= n && a[k][i] == 0) k ++;
        if (k > n) return 0;
        Rep(o, 1, n) swap(a[k][o], a[i][o]);
        Rep(j, i + 1, n) {
            LL z = (a[j][i] * Pow(a[i][i], P)) % P;
            Rep(o, 1, n) {
                (a[j][o] -= a[i][o] * z) %= P;
                if (a[j][o] < 0) a[j][o] += P;
            }
            (p[j] -= p[i] * z) %= P;
            if (p[j] < 0) p[j] += P;
        }
    }
    for (int i = n; i >= 1; i --) {
        Rep(j, i+1, n) {
            (p[i] -= a[i][j] * X[x][j]) %= P;
            if (p[i] < 0) p[i] += P;
        }
        X[x][i] = (p[i] * Pow(a[i][i], P)) % P;
    }
    return 1;
}
LL exgcx ( LL a, LL b, LL &x, LL &y )
{
	!b ? ( x = 1, y = 0 ) : ( exgcx ( b, a % b, y, x ), y -= x * ( a / b ) );
	return x;
}
LL Muit(LL x, LL y, LL md) {
	LL z = 0;
	while (y) {
		if (y % 2 == 1) (z += x) %= md;
		(x += x) %= md, y /= 2;
	}
	return z;
}
void Print(LL X1, LL p1, LL X2, LL p2) {
    LL x, y;
    LL ans = Muit(p2 * ((exgcx(p2%p1, p1, x, y) + p1) % p1) , X1, p1*p2) + Muit(p1 * ((exgcx(p1%p2, p2, x, y) + p2) % p2), X2, p1*p2);
    cout << ans % (p1*p2) <<endl;
}
int main()
{
    cin >> n;
    Rep(i, 1, n) {
        Rep(j, 1, n) scanf ("%d", &A[i][j]);
        scanf ("%s", s);
        Rep(j, 0, strlen(s) - 1) Rep(o, 0, 3) G[o][i] = (G[o][i] * 10 + s[j] - '0') % mod[o];
    }
    Rep(i, 0, 3) {
        if (!Work(i)) continue ;
        if (x1 == -1) x1 = i;
        else if (x2 == -1) x2 = i;
    }
    Rep(i, 1, n) {
        Print(X[x1][i], mod[x1], X[x2][i], mod[x2]);
    }
    return 0;
}

考虑大素数取模下做，整数除法用逆元即可，再用中国剩余定理合并，注意选的素数不能使方程多解(矩阵满秩)。