矩阵

矩阵

//从a[1][1]存到a[n][n],n等于矩阵的长宽
struct Matrix {
    ll mat[N][N];
    Matrix operator * (const Matrix& b)const {
        Matrix ans;
        for (int i = 1; i <= n; i++) {
            for (int j = 1; j <= n; j++) {
                ans.mat[i][j] = 0;
                for (int k = 1; k <= n; k++) {
                    if (mat[i][k] && b.mat[k][j]) {
                        ans.mat[i][j] += mat[i][k] * b.mat[k][j];
                        ans.mat[i][j] %= mod;
                    }
                }
            }
        }
        return ans;
    }
};
​
Matrix q_pow(Matrix a, int b) {
    Matrix ans;
    memset(ans.mat, 0, sizeof(ans));
    for (int i = 1; i <= n; i++) {
        ans.mat[i][i] = 1;
    }
    while (b) {
        if (b & 1)
            ans = ans * a;
        b >>= 1;
        a = a * a;
    }
    return ans;
}

//龟速乘
ll Wuguidechengfa(ll x, ll y)
{
    ll ans = 0;
    while (y)
    {
        if (y & 1) (ans += x) %= mod;
        (x += x) %= mod;
        y >>= 1;
    }
    return ans;
}

//封装

struct matrix {
    vector<vector<ll>> a;
    matrix(ll n, ll m, ll init) : a(vector<vector<ll>>(n, vector<ll>(m, init))) {}  // n rows, m columns
    matrix(ll n, ll m) : matrix(n, m, 0) {}
    matrix(ll n) : matrix(n, n) {}
    matrix() = default;
    matrix(const vector<vector<ll>>& a) : a(a) {}
    vector<ll>& operator[](ll n) { return a[n]; }
    const vector<ll>& operator[](ll n) const { return a[n]; }
    ll sum() const {
        ll ret = 0;
        for (auto& p : a)
            for (ll i : p) ret = (ret + i) % mod;
        return ret;
    }
    matrix operator+(const matrix& b) const {
        ll n, m;
        matrix ret(n = a.size(), m = a[0].size());
        for (ll i = 0; i < n; i++) {
            for (ll j = 0; j < m; j++) {
                ret[i][j] = (a[i][j] + b[i][j]) % mod;
            }
        }
        return move(ret);
    }
    matrix operator-(const matrix& b) const {
        ll n, m;
        matrix ret(n = a.size(), m = a[0].size());
        for (ll i = 0; i < n; i++) {
            for (ll j = 0; j < m; j++) {
                ret[i][j] = ((a[i][j] - b[i][j]) % mod + mod) % mod;
            }
        }
        return move(ret);
    }
    matrix operator*(const matrix& b) const {
        ll n = a.size(), m = a[0].size(), o = b[0].size();
        matrix ret(n, o);
        for (ll i = 0; i < n; i++) {
            for (ll j = 0; j < o; j++) {
                for (ll k = 0; k < m; k++) {
                    ret[i][j] = (ret[i][j] + a[i][k] * b[k][j] % mod) % mod;
                }
            }
        }
        return move(ret);
    }
    static matrix E(ll n) {
        matrix ret(n, n);
        for (ll i = 0; i < n; i++) ret[i][i] = 1;
        return move(ret);
    }
    void print() const {
        for (auto& p : a) {
            for (ll i : p) cerr << i << " ";
            cerr << endl;
        }
    }
};
matrix qpow(matrix a, ll b) {
    matrix c = matrix::E(a[0].size());
    for (; b; b >>= 1, a = a * a) {
        if (b & 1) c = c * a;
    }
    return (c);
}