SRM620 D1-Hard PerfectSquare - こうかいメモ

線形代数やらないといけないと思いつつやってない。

問題

n*n行列が与えられる。そのなかからいくつかの要素を以下の制約を満たすように選ぶ。

各行から選ばれた要素の数が奇数個となる
各列から選ばれた要素の数も奇数個となる
すべての選ばれた要素の積が平方数となる

このような条件を満たす要素の選び方の数 (MOD 1000000007) を求めよ。

n <= 20
行列の要素 <= 1,000,000,000

解法

01連立方程式を立てて階数を求める。

式の数は出現する素因数（多くても900個程度）と各行・各列の数の総和で、変数の数はn^2個。

各素因数に対応する式では、その数を因数として奇数個持つ要素に対応する係数を1にして、その右辺を0にする。

各行・列に対応する式では、その行・列に含まれる要素に対応する係数を1にして、その右辺を1にする。

あとは解があるかどうかを確かめた後に階数から値の定まらない変数の数を求めて、2^(その数) を求めればOK。

コード

typedef unsigned long long ull;
const int MOD = 1000000007;

template <typename T>
vector<int> eratosthenes(T &sieve){
	const int n = sieve.size();
	for(int i = 0; i < n; ++i){ sieve[i] = 1; }
	sieve[0] = sieve[1] = 0;
	vector<int> primes;
	for(int i = 2; i < n; ++i){
		if(!sieve[i]){ continue; }
		primes.push_back(i);
		for(int j = i + i; j < n; j += i){ sieve[j] = 0; }
	}
	return primes;
}

vector<int> eratosthenes(int n){
	vector<bool> sieve(n);
	return eratosthenes(sieve);
}

class BinaryMatrix {
private:
	int N, M;
	vector< vector<ull> > data;
public:
	BinaryMatrix(int N, int M) :
		N(N), M(M), data(N, vector<ull>((M + 63) / 64))
	{ }
	int operator()(int i, int j) const {
		return (data[i][j / 64] >> (j % 64)) & 1;
	}
	void set(int i, int j, int x){
		const ull t = x ? 1 : 0;
		data[i][j / 64] &= ~(1ull << (j % 64));
		data[i][j / 64] |= t << (j % 64);
	}
	void swap_rows(int a, int b){
		if(a == b){ return ; }
		data[a].swap(data[b]);
	}
	int rank() const {
		BinaryMatrix A = *this;
		int r = 0;
		for(int i = 0; r < N && i < M; ++i){
			int nz = -1;
			for(int j = r; nz < 0 && j < N; ++j){
				if(A(j, i)){ nz = j; }
			}
			if(nz < 0){ continue; }
			A.swap_rows(nz, r);
			for(int j = r + 1; j < N; ++j){
				if(A(j, i)){
					const vector<ull> &r0 = A.data[r];
					vector<ull> &r1 = A.data[j];
					for(int k = i / 64; k * 64 < M; ++k){
						r1[k] ^= r0[k];
					}
				}
			}
			++r;
		}
		return r;
	}
};

class PerfectSquare {
public:
	int ways(vector <int> x){
		const int n2 = x.size();
		int n = 0;
		while(n * n < n2){ ++n; }
		const vector<int> primes = eratosthenes(35000);
		vector<int> factors;
		for(int i = 0; i < n2; ++i){
			int t = x[i];
			for(int j = 0; j < primes.size(); ++j){
				const int p = primes[j];
				int c = 0;
				while(t % p == 0){
					t /= p;
					++c;
				}
				if(c % 2 == 1){ factors.push_back(p); }
			}
			if(t > 1){ factors.push_back(t); }
		}
		sort(factors.begin(), factors.end());
		factors.erase(unique(factors.begin(), factors.end()), factors.end());
		BinaryMatrix a(factors.size() + n + n, n2);
		BinaryMatrix b(factors.size() + n + n, n2 + 1);
		for(int i = 0; i < n2; ++i){
			int t = x[i];
			for(int j = 0; j < primes.size(); ++j){
				const int p = primes[j];
				int c = 0;
				while(t % p == 0){
					t /= p;
					++c;
				}
				if(c % 2 == 1){
					const int k = lower_bound(
						factors.begin(), factors.end(), p) - factors.begin();
					a.set(k, i, 1);
					b.set(k, i, 1);
				}
			}
			if(t > 1){
				const int k = lower_bound(
					factors.begin(), factors.end(), t) - factors.begin();
				a.set(k, i, 1);
				b.set(k, i, 1);
			}
		}
		for(int i = 0; i < n; ++i){
			b.set(factors.size() + i, n2, 1);
			for(int j = 0; j < n; ++j){
				a.set(factors.size() + i, i * n + j, 1);
				b.set(factors.size() + i, i * n + j, 1);
			}
			b.set(factors.size() + i + n, n2, 1);
			for(int j = 0; j < n; ++j){
				a.set(factors.size() + n + j, i * n + j, 1);
				b.set(factors.size() + n + j, i * n + j, 1);
			}
		}
		const int ra = a.rank();
		const int rb = b.rank();
		if(ra != rb){ return 0; }
		const int shift = n2 - ra;
		int answer = 1;
		for(int i = 0; i < shift; ++i){
			answer = (answer + answer) % MOD;
		}
		return answer;
	}
};