[BZOJ 4128]Matrix

最新推荐文章于 2019-03-30 12:07:00 发布

原创最新推荐文章于 2019-03-30 12:07:00 发布 · 641 阅读

0 ·

本内容遵循CC 4.0 BY-SA版权协议

数学--高斯消元。线性基专栏收录该内容

5 篇文章

订阅专栏

本文介绍了一种使用裸BSGS（Baby-Step Giant-Step）算法进行矩阵求逆的方法，通过定义特定的矩阵运算符，实现了矩阵的乘法，并利用了BSGS算法来寻找逆元。

裸BSGS矩阵求逆。。

#include <bits/stdc++.h>

using namespace std;
typedef long long ll;

int n, md;
#define maxn 72
#define mod 13331
struct Matrix{
    int a[maxn][maxn];
    void read(){
        for(int i = 0; i < n; i ++)
            for(int j = 0; j < n; j ++)
                scanf("%d", &a[i][j]);
    }
    void print(){
        for(int i = 0; i < n; i ++){
            for(int j = 0; j < n; j ++)
                printf("%d ", a[i][j]);
            printf("\n");
        }
    }

    unsigned long long hash(){
        unsigned long long ret = 0;
        for(int i = 0; i < n; i ++)
            for(int j = 0; j < n; j ++)
                ret = ret * mod + a[i][j];
        return ret;
    }

    void clear(){memset(a, 0, sizeof a);}
    void set(){clear(); for(int i = 0; i < n; i ++)a[i][i] = 1;}
}A, B, mat;

Matrix operator*(const Matrix& a, const Matrix& b){
    Matrix c; c.clear();
    for(int i = 0; i < n; i ++)
        for(int j = 0; j < n; j ++)
            for(int k = 0; k < n; k ++)
                c.a[i][j] = (c.a[i][j] + (long long)a.a[i][k] * b.a[k][j]) % md;
    return c;
}

map<unsigned long long, int> M;

ll power_mod(ll a, ll b){
    ll ret = 1;
    while(b > 0){
        if(b & 1)ret = ret * a % md;
        b >>= 1;
        a = a * a % md;
    }return ret;
}

Matrix getinv(Matrix a){
    Matrix ret; ret.set();
    for(int i = 0; i < n; i ++){
        for(int j = i; j < n; j ++){
            if(a.a[j][i]){
                for(int k = 0; k < n; k ++){
                    swap(a.a[i][k], a.a[j][k]);
                    swap(ret.a[i][k], ret.a[j][k]);
                }

                ll inv = power_mod(a.a[i][i], md-2);
                for(int k = 0; k < n; k ++){
                    a.a[i][k] = a.a[i][k] * inv % md;
                    ret.a[i][k] = ret.a[i][k] * inv % md;
                }
                break;
            }
        }
        //if(!a.a[i][i])continue;
        for(int j = 0; j < n; j ++){
            if(i == j)continue;
            ll tmp = a.a[j][i];
            for(int k = 0; k < n; k ++){
                a.a[j][k] = (a.a[j][k] - a.a[i][k] * tmp % md + md) % md; 
                ret.a[j][k] = (ret.a[j][k] - ret.a[i][k] * tmp % md + md) % md;
            }
        }
    }
    return ret;
}

int main(){
    scanf("%d%d", &n, &md);
    A.read(), B.read();
    int q = sqrt(md) + 1; 
    mat.set();
    unsigned long long now = B.hash(), tmp = mat.hash();
    if(tmp == now)return puts("0"), 0;

    for(int i = 1; i <= q; i ++){
        mat = mat * A;
        tmp = mat.hash();
        if(tmp == now){
            printf("%d\n", i);
            return 0;
        }
        M[tmp] = i;
    }

    Matrix inv = getinv(mat);

    for(int i = 1; i <= q; i ++){
        B = B * inv;
        now = B.hash();
        if(M.count(now)){
            printf("%d\n", i * q + M[now]);
            return 0;
        }
    }

    return 0;
}