cp-library

This documentation is automatically generated by online-judge-tools/verification-helper

View the Project on GitHub rniya/cp-library

:heavy_check_mark: Matrix on $\mathbb{F}_2$
(src/matrix/F2Matrix.hpp)

概要

std::bitset<> を利用して $\mathbb{F}_2$ 上での各種行列処理を高速化したライブラリ。 固定長の都合上、列数以上の値を引数 MAX_W として渡す必要がある。

下表ではワードサイズを $w = 64$ とする。

メンバ関数 効果 時間計算量
empty() 行列が空か否かを返す． $\mathrm{O}(1)$
size() 行数を返す． $\mathrm{O}(1)$
height() 行数を返す． $\mathrm{O}(1)$
width() 列数を返す． $\mathrm{O}(1)$
identity(N) $N \times N$ 単位行列を返す． $\mathrm{O} \left( \frac{N^2}{w} \right)$
加算 $N \times M$ 行列 $A$ に $B$ を加算する． $\mathrm{O} \left( \frac{N M}{w} \right)$
減算 $N \times M$ 行列 $A$ に $B$ を減算する． $\mathrm{O} \left( \frac{N M}{w} \right)$
乗算 $N \times M$ 行列 $A$ に $M \times L$ 行列 $B$ を減算する． $\mathrm{O} \left( \frac{N M L}{w} \right)$
pow(n) $N \times N$ 正方行列 $A$ を $n$ 乗した行列を返す． $\mathrm{O} \left( \frac{N^3 \log n}{w} \right)$
rank() $N \times M$ 行列 $A$ の rank を返す． $\mathrm{O} \left( \frac{N M^2}{w} \right)$
det() $N \times N$ 正方行列 $A$ の determinant を返す． $\mathrm{O} \left( \frac{N^3}{w} \right)$
inv() $N \times N$ 正方行列 $A$ の逆行列を返す（行列が正則である必要がある）． $\mathrm{O} \left( \frac{N^3}{w} \right)$
system_of_linear_equations(b) $N \times M$ 行列 $A$ と長さ $N$ のｔ縦ベクトル $b$ について $A x = b$ という線形方程式系を考える．この方程式系に解が存在しない場合は空配列を返し，存在する場合は解のうちの 1 つ及び解空間の基底ベクトルをこの順にまとめた配列を返す． $\mathrm{O} \left( \frac{N M^2}{w} \right)$

出題例

Verified with

Code

#pragma once
#include <bitset>
#include <cassert>
#include <utility>
#include <vector>

template <int MAX_W> struct F2Matrix {
    int H, W;

    std::vector<std::bitset<MAX_W>> A;

    F2Matrix(int H, int W) : H(H), W(W), A(H) {
        assert(W <= MAX_W);
        for (int i = 0; i < H; i++) A[i].reset();
    }

    bool empty() const { return A.empty(); }

    int size() const { return H; }

    int height() const { return H; }

    int width() const { return W; }

    inline const std::bitset<MAX_W>& operator[](int i) const { return A[i]; }

    inline std::bitset<MAX_W>& operator[](int i) { return A[i]; }

    static F2Matrix identity(int n) {
        F2Matrix res(n, n);
        for (int i = 0; i < n; i++) res[i][i] = 1;
        return res;
    }

    F2Matrix& operator+=(const F2Matrix& B) {
        assert(H == B.height() and W == B.width());
        for (int i = 0; i < H; i++) (*this)[i] ^= B[i];
        return *this;
    }

    F2Matrix& operator-=(const F2Matrix& B) { return *this += B; }

    F2Matrix& operator*=(const F2Matrix& B) {
        assert(W == B.height());
        F2Matrix C(H, B.width());
        for (int i = 0; i < H; i++) {
            for (int j = 0; j < W; j++) {
                if (A[i][j]) {
                    C[i] ^= B[j];
                }
            }
        }
        std::swap(A, C.A);
        return *this;
    }

    F2Matrix operator+(const F2Matrix& B) const { return F2Matrix(*this) += B; }

    F2Matrix operator-(const F2Matrix& B) const { return F2Matrix(*this) -= B; }

    F2Matrix operator*(const F2Matrix& B) const { return F2Matrix(*this) *= B; }

    F2Matrix pow(long long n) const {
        assert(H == W);
        assert(0 <= n);
        F2Matrix x = *this, r = identity(H);
        while (n) {
            if (n & 1) r *= x;
            x *= x;
            n >>= 1;
        }
        return r;
    }

    int rank() const { return F2Matrix(*this).gauss_jordan(); }

    int det() const {
        assert(H == W);
        return rank() == H ? 1 : 0;
    }

    F2Matrix inv() const {
        assert(height() == width());
        int n = height();
        F2Matrix B(*this), C = identity(n);
        for (int j = 0; j < n; j++) {
            int pivot = -1;
            for (int i = j; i < n; i++) {
                if (B[i][j]) {
                    pivot = i;
                    break;
                }
            }
            assert(pivot != -1);
            if (pivot != j) {
                std::swap(B[pivot], B[j]);
                std::swap(C[pivot], C[j]);
            }
            for (int i = 0; i < n; i++) {
                if (i != j and B[i][j]) {
                    B[i] ^= B[j];
                    C[i] ^= C[j];
                }
            }
        }
        return C;
    }

    std::vector<std::bitset<MAX_W>> system_of_linear_equations(const std::vector<bool>& b) {
        assert(W + 1 <= MAX_W);
        assert(int(b.size()) == H);
        F2Matrix B(*this);
        for (int i = 0; i < H; i++) B[i][W] = b[i];
        int rank = B.gauss_jordan(W);
        for (int i = rank; i < H; i++) {
            if (B[i][W]) {
                return {};
            }
        }
        std::vector<std::bitset<MAX_W>> res(1);
        std::vector<int> pivot(W, -1);
        for (int i = 0, j = 0; i < rank; i++) {
            while (not B[i][j]) j++;
            res[0][j] = B[i][W];
            pivot[j] = i;
        }
        for (int j = 0; j < W; j++) {
            if (pivot[j] != -1) continue;
            std::bitset<MAX_W> x;
            x[j] = 1;
            for (int k = 0; k < j; k++) {
                if (pivot[k] != -1) {
                    x[k] = B[pivot[k]][j];
                }
            }
            res.emplace_back(x);
        }
        return res;
    }

  private:
    int gauss_jordan(int pivot_end = -1) {
        if (empty()) return 0;
        if (pivot_end == -1) pivot_end = W;
        assert(pivot_end <= MAX_W);
        int rank = 0;
        for (int j = 0; j < pivot_end; j++) {
            int pivot = -1;
            for (int i = rank; i < H; i++) {
                if ((*this)[i][j]) {
                    pivot = i;
                    break;
                }
            }
            if (pivot == -1) continue;
            if (pivot != rank) std::swap((*this)[pivot], (*this)[rank]);
            for (int i = 0; i < H; i++) {
                if (i != rank and (*this)[i][j]) {
                    (*this)[i] ^= (*this)[rank];
                }
            }
            rank++;
        }
        return rank;
    }
};
#line 2 "src/matrix/F2Matrix.hpp"
#include <bitset>
#include <cassert>
#include <utility>
#include <vector>

template <int MAX_W> struct F2Matrix {
    int H, W;

    std::vector<std::bitset<MAX_W>> A;

    F2Matrix(int H, int W) : H(H), W(W), A(H) {
        assert(W <= MAX_W);
        for (int i = 0; i < H; i++) A[i].reset();
    }

    bool empty() const { return A.empty(); }

    int size() const { return H; }

    int height() const { return H; }

    int width() const { return W; }

    inline const std::bitset<MAX_W>& operator[](int i) const { return A[i]; }

    inline std::bitset<MAX_W>& operator[](int i) { return A[i]; }

    static F2Matrix identity(int n) {
        F2Matrix res(n, n);
        for (int i = 0; i < n; i++) res[i][i] = 1;
        return res;
    }

    F2Matrix& operator+=(const F2Matrix& B) {
        assert(H == B.height() and W == B.width());
        for (int i = 0; i < H; i++) (*this)[i] ^= B[i];
        return *this;
    }

    F2Matrix& operator-=(const F2Matrix& B) { return *this += B; }

    F2Matrix& operator*=(const F2Matrix& B) {
        assert(W == B.height());
        F2Matrix C(H, B.width());
        for (int i = 0; i < H; i++) {
            for (int j = 0; j < W; j++) {
                if (A[i][j]) {
                    C[i] ^= B[j];
                }
            }
        }
        std::swap(A, C.A);
        return *this;
    }

    F2Matrix operator+(const F2Matrix& B) const { return F2Matrix(*this) += B; }

    F2Matrix operator-(const F2Matrix& B) const { return F2Matrix(*this) -= B; }

    F2Matrix operator*(const F2Matrix& B) const { return F2Matrix(*this) *= B; }

    F2Matrix pow(long long n) const {
        assert(H == W);
        assert(0 <= n);
        F2Matrix x = *this, r = identity(H);
        while (n) {
            if (n & 1) r *= x;
            x *= x;
            n >>= 1;
        }
        return r;
    }

    int rank() const { return F2Matrix(*this).gauss_jordan(); }

    int det() const {
        assert(H == W);
        return rank() == H ? 1 : 0;
    }

    F2Matrix inv() const {
        assert(height() == width());
        int n = height();
        F2Matrix B(*this), C = identity(n);
        for (int j = 0; j < n; j++) {
            int pivot = -1;
            for (int i = j; i < n; i++) {
                if (B[i][j]) {
                    pivot = i;
                    break;
                }
            }
            assert(pivot != -1);
            if (pivot != j) {
                std::swap(B[pivot], B[j]);
                std::swap(C[pivot], C[j]);
            }
            for (int i = 0; i < n; i++) {
                if (i != j and B[i][j]) {
                    B[i] ^= B[j];
                    C[i] ^= C[j];
                }
            }
        }
        return C;
    }

    std::vector<std::bitset<MAX_W>> system_of_linear_equations(const std::vector<bool>& b) {
        assert(W + 1 <= MAX_W);
        assert(int(b.size()) == H);
        F2Matrix B(*this);
        for (int i = 0; i < H; i++) B[i][W] = b[i];
        int rank = B.gauss_jordan(W);
        for (int i = rank; i < H; i++) {
            if (B[i][W]) {
                return {};
            }
        }
        std::vector<std::bitset<MAX_W>> res(1);
        std::vector<int> pivot(W, -1);
        for (int i = 0, j = 0; i < rank; i++) {
            while (not B[i][j]) j++;
            res[0][j] = B[i][W];
            pivot[j] = i;
        }
        for (int j = 0; j < W; j++) {
            if (pivot[j] != -1) continue;
            std::bitset<MAX_W> x;
            x[j] = 1;
            for (int k = 0; k < j; k++) {
                if (pivot[k] != -1) {
                    x[k] = B[pivot[k]][j];
                }
            }
            res.emplace_back(x);
        }
        return res;
    }

  private:
    int gauss_jordan(int pivot_end = -1) {
        if (empty()) return 0;
        if (pivot_end == -1) pivot_end = W;
        assert(pivot_end <= MAX_W);
        int rank = 0;
        for (int j = 0; j < pivot_end; j++) {
            int pivot = -1;
            for (int i = rank; i < H; i++) {
                if ((*this)[i][j]) {
                    pivot = i;
                    break;
                }
            }
            if (pivot == -1) continue;
            if (pivot != rank) std::swap((*this)[pivot], (*this)[rank]);
            for (int i = 0; i < H; i++) {
                if (i != rank and (*this)[i][j]) {
                    (*this)[i] ^= (*this)[rank];
                }
            }
            rank++;
        }
        return rank;
    }
};
Back to top page