BZOJ4128 Matrix 【BSGS】

BZOJ4128 Matrix

---

Description

给定矩阵A,B和模数p,求最小的x满足
A^x = B (mod p)

Input

第一行两个整数n和p，表示矩阵的阶和模数,接下来一个n * n的矩阵A.接下来一个n * n的矩阵B

Output

输出一个正整数，表示最小的可能的x，数据保证在p内有解

Sample Input

2 7
1 1
1 0
5 3
3 2

Sample Output

4

HINT

对于100%的数据，n <= 70,p <=19997，p为质数,0<= A_{ij},B_{ij}< p
保证A有逆

---

一看就是BSGS嘛
别人出题人都保证了矩阵有逆你还在犹豫什么呢

Guass消元求逆
大概就是拿一个单位矩阵出来，原矩阵该怎么消怎么消，单位矩阵跟着原矩阵一起操作就行了
感觉非常简单

然后我居然非常傻逼地忽略了题目中已给的模数自己定义了一个模数

阿弥陀佛

---

#include<bits/stdc++.h>
using namespace std;
#define N 80
#define P 20000
#define M 37373
int n,p;
struct Matrix{
    int t[80][80];
    Matrix(){memset(t,0,sizeof(t));}
};
struct Edge{
    int v,next,val;
    Edge(){}
    Edge(int v,int next,int val):v(v),next(next),val(val){}
}E[M+2];
int head[M+2],tot=0;
void add(int u,int v,int val){
    E[++tot]=Edge(v,head[u],val);
    head[u]=tot;
}
int query(int Hash_value){
    for(int i=head[Hash_value%M];i!=-1;i=E[i].next)
        if(E[i].val==Hash_value)return E[i].v;
    return -1;
}
int add(int a,int b){return (a+b)%p;}
int dec(int a,int b){return (a-b+p)%p;}
int mul(int a,int b){return 1ll*a*b%p;}
Matrix mul(Matrix a,Matrix b){
    Matrix c;
    for(int i=1;i<=n;i++)
        for(int j=1;j<=n;j++)
            for(int k=1;k<=n;k++)
                c.t[i][j]=add(c.t[i][j],mul(a.t[i][k],b.t[k][j]));
    return c;
}
int fast_pow(int a,int b){
    int res=1;
    while(b){
        if(b&1)res=mul(res,a);
        b>>=1;
        a=mul(a,a);
    }
    return res;
}
Matrix inv_Matrix(Matrix a){
    Matrix res;
    for(int i=1;i<=n;i++)res.t[i][i]=1;
    for(int i=1;i<=n;i++){
        int k=0;
        for(int j=i;j<=n;j++)if(a.t[j][i]){k=j;break;}
        for(int j=1;j<=n;j++){
            swap(a.t[i][j],a.t[k][j]);
            swap(res.t[i][j],res.t[k][j]);
        }
        int inv=fast_pow(a.t[i][i],p-2);
        for(int j=1;j<=n;j++){
            a.t[i][j]=mul(a.t[i][j],inv);
            res.t[i][j]=mul(res.t[i][j],inv);
        }
        for(int j=1;j<=n;j++){
            if(i==j)continue;
            int pre=a.t[j][i];
            for(int w=1;w<=n;w++){
                a.t[j][w]=dec(a.t[j][w],mul(pre,a.t[i][w]));
                res.t[j][w]=dec(res.t[j][w],mul(pre,res.t[i][w]));
            }
        }
    }
    return res;
}
int Hash(Matrix a){return a.t[n][n];}
int BSGS(Matrix A,Matrix B,int c){
    int siz=sqrt(c)+1;
    Matrix tmp;
    for(int i=1;i<=n;i++)tmp.t[i][i]=1;
    for(int i=0;i<siz;i++){
        int Hash_value=Hash(tmp);
        if(query(Hash_value)==-1)add(Hash_value%M,i,Hash_value);
        tmp=mul(tmp,A);
    }
    Matrix inv_matrix=inv_Matrix(tmp);
    for(int i=0;i<=siz;i++){
        int Hash_value=Hash(B);
        int tip=query(Hash_value);
        if(tip!=-1)return tip+i*siz;
        B=mul(B,inv_matrix);
    }
    return -1;
}
int main(){
    freopen("4128.in","r",stdin);
    memset(head,-1,sizeof(head));
    Matrix A,B;
    scanf("%d%d",&n,&p);
    for(int i=1;i<=n;i++)
        for(int j=1;j<=n;j++)scanf("%d",&A.t[i][j]);
    for(int i=1;i<=n;i++)
        for(int j=1;j<=n;j++)scanf("%d",&B.t[i][j]);
    printf("%d\n",BSGS(A,B,p));
    return 0;
}