firiexp's Library

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

View the Project on GitHub firiexp/library

:heavy_check_mark: 逆行列(Inverse Matrix)
(math/inverse_matrix.cpp)

説明

modint 行列の逆行列を計算する。 Gauss-Jordan 消去を使い、計算量は $O(N^3)$。

できること

使い方

Avector<vector<mint>> で渡す。 返り値が空でなければ A^{-1} が入っている。

Depends on

Verified with

Code

#include "../util/modint.cpp"

vector<vector<mint>> inverse_matrix(const vector<vector<mint>> &A) {
    int n = A.size();
    vector<vector<mint>> B(n, vector<mint>(2 * n));
    for (int i = 0; i < n; ++i) {
        for (int j = 0; j < n; ++j) B[i][j] = A[i][j];
        B[i][n + i] = 1;
    }
    int rank = 0;
    for (int col = 0; col < n; ++col) {
        int pivot = -1;
        for (int row = rank; row < n; ++row) {
            if (B[row][col].val) {
                pivot = row;
                break;
            }
        }
        if (pivot == -1) return {};
        swap(B[pivot], B[rank]);
        mint inv = B[rank][col].inv();
        for (int j = 0; j < 2 * n; ++j) B[rank][j] *= inv;
        for (int row = 0; row < n; ++row) {
            if (row != rank && B[row][col].val) {
                mint coeff = B[row][col];
                for (int j = 0; j < 2 * n; ++j) {
                    B[row][j] -= B[rank][j] * coeff;
                }
            }
        }
        ++rank;
    }
    vector<vector<mint>> res(n, vector<mint>(n));
    for (int i = 0; i < n; ++i) {
        for (int j = 0; j < n; ++j) res[i][j] = B[i][n + j];
    }
    return res;
}

/**
 * @brief 逆行列(Inverse Matrix)
 */
#line 1 "util/modint.cpp"



template <uint Mod>
struct modint {
    uint val;
public:
    static modint raw(int v) { modint x; x.val = v; return x; }
    static constexpr uint get_mod() { return Mod; }
    static constexpr uint M() { return Mod; }
    modint() : val(0) {}
    template <class T>
    modint(T v) { ll x = (ll)(v % (ll)(Mod)); if (x < 0) x += Mod; val = uint(x); }
    modint(bool v) { val = ((unsigned int)(v) % Mod); }
    uint &value() noexcept { return val; }
    const uint &value() const noexcept { return val; }
    modint& operator++() { val++; if (val == Mod) val = 0; return *this; }
    modint& operator--() { if (val == 0) val = Mod; val--; return *this; }
    modint operator++(int) { modint result = *this; ++*this; return result; }
    modint operator--(int) { modint result = *this; --*this; return result; }
    modint& operator+=(const modint& b) { val += b.val; if (val >= Mod) val -= Mod; return *this; }
    modint& operator-=(const modint& b) { val -= b.val; if (val >= Mod) val += Mod; return *this; }
    modint& operator*=(const modint& b) { ull z = val; z *= b.val; val = (uint)(z % Mod); return *this; }
    modint& operator/=(const modint& b) { return *this = *this * b.inv(); }
    modint operator+() const { return *this; }
    modint operator-() const { return modint() - *this; }
    modint pow(long long n) const { modint x = *this, r = 1; while (n) { if (n & 1) r *= x; x *= x; n >>= 1; } return r; }
    modint inv() const { return pow(Mod - 2); }
    friend modint operator+(const modint& a, const modint& b) { return modint(a) += b; }
    friend modint operator-(const modint& a, const modint& b) { return modint(a) -= b; }
    friend modint operator*(const modint& a, const modint& b) { return modint(a) *= b; }
    friend modint operator/(const modint& a, const modint& b) { return modint(a) /= b; }
    friend bool operator==(const modint& a, const modint& b) { return a.val == b.val; }
    friend bool operator!=(const modint& a, const modint& b) { return a.val != b.val; }
};
using mint = modint<MOD>;
#define FIRIEXP_LIBRARY_MINT_ALIAS_DEFINED

/**
 * @brief modint(固定MOD)
 */


#line 2 "math/inverse_matrix.cpp"

vector<vector<mint>> inverse_matrix(const vector<vector<mint>> &A) {
    int n = A.size();
    vector<vector<mint>> B(n, vector<mint>(2 * n));
    for (int i = 0; i < n; ++i) {
        for (int j = 0; j < n; ++j) B[i][j] = A[i][j];
        B[i][n + i] = 1;
    }
    int rank = 0;
    for (int col = 0; col < n; ++col) {
        int pivot = -1;
        for (int row = rank; row < n; ++row) {
            if (B[row][col].val) {
                pivot = row;
                break;
            }
        }
        if (pivot == -1) return {};
        swap(B[pivot], B[rank]);
        mint inv = B[rank][col].inv();
        for (int j = 0; j < 2 * n; ++j) B[rank][j] *= inv;
        for (int row = 0; row < n; ++row) {
            if (row != rank && B[row][col].val) {
                mint coeff = B[row][col];
                for (int j = 0; j < 2 * n; ++j) {
                    B[row][j] -= B[rank][j] * coeff;
                }
            }
        }
        ++rank;
    }
    vector<vector<mint>> res(n, vector<mint>(n));
    for (int i = 0; i < n; ++i) {
        for (int j = 0; j < n; ++j) res[i][j] = B[i][n + j];
    }
    return res;
}

/**
 * @brief 逆行列(Inverse Matrix)
 */
Back to top page