kyopro-lib-cpp

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

View the Project on GitHub packer-jp/kyopro-lib-cpp

:heavy_check_mark: test/judge.yosupo.jp/Determinant_of_Matrix.0.test.cpp

Depends on

Code

#define PROBLEM "https://judge.yosupo.jp/problem/matrix_det"
#include "../../math/matrix.hpp"
#include "../../math/modint.hpp"

int main() {
    using mint = modint998244353;
    struct mint_field {
        using V = mint;
        static V zero() { return 0; }
        static V one() { return 1; }
    };
    ll n;
    cin >> n;
    matrix<mint_field> a(n, n);
    for (ll i : rep(n)) {
        for (ll j : rep(n)) { cin >> a[i][j]; }
    }
    cout << a.det() << endl;
}
#line 1 "test/judge.yosupo.jp/Determinant_of_Matrix.0.test.cpp"
#define PROBLEM "https://judge.yosupo.jp/problem/matrix_det"
#line 2 "math/matrix.hpp"

#line 2 "template.hpp"

#include <bits/stdc++.h>
using namespace std;

#define all(a) begin(a), end(a)
#define rall(a) rbegin(a), rend(a)
#define uniq(a) (a).erase(unique(all(a)), (a).end())
#define t0 first
#define t1 second
using ll = long long;
using ull = unsigned long long;
using pll = pair<ll, ll>;
using vll = vector<ll>;
constexpr double pi = 3.14159265358979323846;
constexpr ll dy[9] = {0, 1, 0, -1, 1, 1, -1, -1, 0};
constexpr ll dx[9] = {1, 0, -1, 0, 1, -1, -1, 1, 0};
constexpr ll sign(ll a) { return (a > 0) - (a < 0); }
constexpr ll fdiv(ll a, ll b) { return a / b - ((a ^ b) < 0 && a % b); }
constexpr ll cdiv(ll a, ll b) { return -fdiv(-a, b); }
constexpr ll pw(ll n) { return 1ll << n; }
constexpr ll flg(ll n) { return 63 - __builtin_clzll(n); }
constexpr ll clg(ll n) { return n == 1 ? 0 : flg(n - 1) + 1; }
constexpr ll safemod(ll x, ll mod) { return (x % mod + mod) % mod; }
template <typename T> using priority_queue_rev = priority_queue<T, vector<T>, greater<T>>;
template <typename T> constexpr T sq(const T &a) { return a * a; }
template <typename T, typename U> constexpr bool chmax(T &a, const U &b) { return a < b ? a = b, true : false; }
template <typename T, typename U> constexpr bool chmin(T &a, const U &b) { return a > b ? a = b, true : false; }
template <typename T, typename U> ostream &operator<<(ostream &os, const pair<T, U> &a) {
    os << "(" << a.first << ", " << a.second << ")";
    return os;
}
template <typename T, typename U, typename V> ostream &operator<<(ostream &os, const tuple<T, U, V> &a) {
    os << "(" << get<0>(a) << ", " << get<1>(a) << ", " << get<2>(a) << ")";
    return os;
}
template <typename T> ostream &operator<<(ostream &os, const vector<T> &a) {
    os << "(";
    for (auto itr = a.begin(); itr != a.end(); ++itr) os << *itr << (next(itr) != a.end() ? ", " : "");
    os << ")";
    return os;
}
template <typename T> ostream &operator<<(ostream &os, const set<T> &a) {
    os << "(";
    for (auto itr = a.begin(); itr != a.end(); ++itr) os << *itr << (next(itr) != a.end() ? ", " : "");
    os << ")";
    return os;
}
template <typename T> ostream &operator<<(ostream &os, const multiset<T> &a) {
    os << "(";
    for (auto itr = a.begin(); itr != a.end(); ++itr) os << *itr << (next(itr) != a.end() ? ", " : "");
    os << ")";
    return os;
}
template <typename T, typename U> ostream &operator<<(ostream &os, const map<T, U> &a) {
    os << "(";
    for (auto itr = a.begin(); itr != a.end(); ++itr) os << *itr << (next(itr) != a.end() ? ", " : "");
    os << ")";
    return os;
}
#ifdef ONLINE_JUDGE
#define dump(...) (void(0))
#else
void debug() { cerr << endl; }
template <typename Head, typename... Tail> void debug(Head &&head, Tail &&... tail) {
    cerr << head;
    if (sizeof...(Tail)) cerr << ", ";
    debug(tail...);
}
#define dump(...) cerr << __LINE__ << ": " << #__VA_ARGS__ << " = ", debug(__VA_ARGS__)
#endif
struct rep {
    struct itr {
        ll v;
        itr(ll v) : v(v) {}
        void operator++() { ++v; }
        ll operator*() const { return v; }
        bool operator!=(itr i) const { return v < *i; }
    };
    ll l, r;
    rep(ll l, ll r) : l(l), r(r) {}
    rep(ll r) : rep(0, r) {}
    itr begin() const { return l; };
    itr end() const { return r; };
};
struct per {
    struct itr {
        ll v;
        itr(ll v) : v(v) {}
        void operator++() { --v; }
        ll operator*() const { return v; }
        bool operator!=(itr i) const { return v > *i; }
    };
    ll l, r;
    per(ll l, ll r) : l(l), r(r) {}
    per(ll r) : per(0, r) {}
    itr begin() const { return r - 1; };
    itr end() const { return l - 1; };
};
struct io_setup {
    static constexpr int PREC = 20;
    io_setup() {
        cout << fixed << setprecision(PREC);
        cerr << fixed << setprecision(PREC);
    };
} iOS;
#line 4 "math/matrix.hpp"

template <typename S> struct matrix {
    using V = typename S::V;
    vector<vector<V>> val;
    matrix(int n, int m) : matrix(vector(n, vector(m, S::zero()))) {}
    matrix(const vector<vector<V>> &src) : val(src) {}
    vector<V> &operator[](int i) { return val[i]; }
    const vector<V> &operator[](int i) const { return val[i]; }
    int height() const { return val.size(); }
    int width() const { return val[0].size(); }
    static matrix id(int n) {
        matrix ret(n, n);
        for (int i : rep(n)) ret[i][i] = S::one();
        return ret;
    }
    void row_add(int i, int j, V a) {
        for (int k : rep(width())) { val[i][k] += val[j][k] * a; }
    }
    bool place_nonzero(int i, int j) {
        for (int k : rep(i, height())) {
            if (val[k][j] != S::zero()) {
                if (k > i) row_add(i, k, S::one());
                break;
            }
        }
        return val[i][j] != S::zero();
    }
    matrix upper_triangular() const {
        matrix ret(*this);
        for (int i = 0, j = 0; i < height() && j < width(); j++) {
            if (!ret.place_nonzero(i, j)) continue;
            for (int k : rep(i + 1, height())) ret.row_add(k, i, -ret[k][j] / ret[i][j]);
            i++;
        }
        return ret;
    }
    V det() const {
        V ret = S::one();
        matrix ut = upper_triangular();
        for (int i : rep(height())) ret *= ut[i][i];
        return ret;
    }
    matrix inv() const {
        matrix ex(height(), width() << 1);
        for (int i : rep(height())) {
            for (int j : rep(width())) { ex[i][j] = val[i][j]; }
        }
        for (int i : rep(height())) ex[i][width() + i] = S::one();
        matrix ut = ex.upper_triangular();
        for (int i : per(height())) {
            ut.row_add(i, i, S::one() / ut[i][i] - S::one());
            for (int j : rep(i)) ut.row_add(j, i, -ut[j][i] / ut[i][i]);
        }
        matrix ret(height(), width());
        for (int i : rep(height())) {
            for (int j : rep(width())) { ret[i][j] = ut[i][width() + j]; }
        }
        return ret;
    }
    matrix pow(ll k) const {
        matrix ret = matrix::id(height()), mul(*this);
        while (k) {
            if (k & 1) ret *= mul;
            mul *= mul;
            k >>= 1;
        }
        return ret;
    }
    matrix &operator+=(const matrix &a) {
        for (int i : rep(height())) {
            for (int j : rep(width())) { val[i][j] += a[i][j]; }
        }
        return *this;
    }
    matrix &operator-=(const matrix &a) {
        for (int i : rep(height())) {
            for (int j : rep(width())) { val[i][j] -= a[i][j]; }
        }
        return *this;
    }
    matrix &operator*=(const matrix &a) {
        matrix res(height(), a.width());
        for (int i : rep(height())) {
            for (int j : rep(a.width())) {
                for (int k : rep(width())) { res[i][j] += val[i][k] * a[k][j]; }
            }
        }
        val.swap(res.val);
        return *this;
    }
    matrix &operator/=(const matrix &a) { return *this *= a.inv(); }
    bool operator==(const matrix &a) const { return val == a.val; }
    bool operator!=(const matrix &a) const { return rel_ops::operator!=(*this, a); }
    matrix operator+() const { return *this; }
    matrix operator-() const { return matrix(height(), width()) -= *this; }
    matrix operator+(const matrix &a) const { return matrix(*this) += a; }
    matrix operator-(const matrix &a) const { return matrix(*this) -= a; }
    matrix operator*(const matrix &a) const { return matrix(*this) *= a; }
    matrix operator/(const matrix &a) const { return matrix(*this) /= a; }
};

struct double_field {
    using V = double;
    static V zero() { return 0.0; }
    static V one() { return 1.0; }
};

template <> bool matrix<double_field>::place_nonzero(int i, int j) {
    static constexpr double EPS = 1e-12;
    for (int k : rep(i + 1, height())) {
        if (abs(val[k][j]) > abs(val[i][j])) {
            swap(val[i], val[k]);
            row_add(i, i, -2.0);
        }
    }
    return abs(val[i][j]) > EPS;
}
#line 2 "math/modint.hpp"

#line 4 "math/modint.hpp"

template <typename M> struct modint {
    ll val;
    modint(ll val = 0) : val(val >= 0 ? val % M::mod : (M::mod - (-val) % M::mod) % M::mod) {}
    static ll mod() { return M::mod; }
    modint inv() const {
        ll a = val, b = M::mod, u = 1, v = 0, t;
        while (b > 0) {
            t = a / b;
            swap(a -= t * b, b);
            swap(u -= t * v, v);
        }
        return u;
    }
    modint pow(ll k) const {
        modint ret = 1, mul = val;
        while (k) {
            if (k & 1) ret *= mul;
            mul *= mul;
            k >>= 1;
        }
        return ret;
    }
    modint &operator+=(const modint &a) {
        if ((val += a.val) >= M::mod) val -= M::mod;
        return *this;
    }
    modint &operator-=(const modint &a) {
        if ((val += M::mod - a.val) >= M::mod) val -= M::mod;
        return *this;
    }
    modint &operator*=(const modint &a) {
        (val *= a.val) %= M::mod;
        return *this;
    }
    modint &operator/=(const modint &a) { return *this *= a.inv(); }
    modint operator+() const { return *this; }
    modint operator-() const { return modint(-val); }
    friend bool operator==(const modint &a, const modint &b) { return a.val == b.val; }
    friend bool operator!=(const modint &a, const modint &b) { return rel_ops::operator!=(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 modint operator/(const modint &a, const modint &b) { return modint(a) /= b; }
    friend istream &operator>>(istream &is, modint &a) {
        ll val;
        is >> val;
        a = modint(val);
        return is;
    }
    friend ostream &operator<<(ostream &os, const modint &a) { return os << a.val; }
};

struct _998244353 {
    constexpr static ll mod = 998244353;
};
struct _1000000007 {
    constexpr static ll mod = 1000000007;
};
using modint998244353 = modint<_998244353>;
using modint1000000007 = modint<_1000000007>;

struct arbitrary {
    static ll mod;
};
ll arbitrary::mod;
#line 4 "test/judge.yosupo.jp/Determinant_of_Matrix.0.test.cpp"

int main() {
    using mint = modint998244353;
    struct mint_field {
        using V = mint;
        static V zero() { return 0; }
        static V one() { return 1; }
    };
    ll n;
    cin >> n;
    matrix<mint_field> a(n, n);
    for (ll i : rep(n)) {
        for (ll j : rep(n)) { cin >> a[i][j]; }
    }
    cout << a.det() << endl;
}
Back to top page