zoukankan      html  css  js  c++  java
  • 【题解】Matrix BZOJ 4128 矩阵求逆 离散对数 大步小步算法

    传送门:http://www.lydsy.com/JudgeOnline/problem.php?id=4128

    大水题一道

    使用大步小步算法,把数字的运算换成矩阵的运算就好了

    矩阵求逆?这么基础的线代算法我也不想多说,还是自行百度吧

    需要注意的是矩阵没有交换律,所以在计算$Bcdot A^{-m}$的时候不要把顺序搞混

    代码:

      1 #include <cstring>
      2 #include <cstdio>
      3 #include <algorithm>
      4 #include <map>
      5 #include <cmath>
      6 
      7 using namespace std;
      8 const int MAXN = 70;
      9 
     10 int n, p;
     11 
     12 struct Matrix {
     13     int a[MAXN][MAXN];
     14     int *operator[]( int idx ) {
     15         return a[idx];
     16     }
     17     const int *operator[]( int idx ) const {
     18         return a[idx];
     19     }
     20     bool operator<( const Matrix &rhs ) const { // 用于map
     21         for( int i = 0; i < n; ++i )
     22             for( int j = 0; j < n; ++j )
     23                 if( a[i][j] < rhs[i][j] ) return true;
     24                 else if( a[i][j] > rhs[i][j] ) return false;
     25         return false;
     26     }
     27     void unit() { // 填成n阶单位矩阵
     28         memset( a, 0, sizeof(a) );
     29         for( int i = 0; i < n; ++i ) a[i][i] = 1;
     30     }
     31     void clear() { // 填成全0
     32         memset( a, 0, sizeof(a) );
     33     }
     34 };
     35 Matrix A, B; // 题目中给定的A,B矩阵
     36 
     37 int pow_mod( int a, int b ) { // 整数快速幂
     38     int ans = 1;
     39     while(b) {
     40         if( b&1 ) ans = ans * a % p;
     41         a = a * a % p;
     42         b >>= 1;
     43     }
     44     return ans;
     45 }
     46 int inv( int a ) { // 整数求逆
     47     return pow_mod(a,p-2);
     48 }
     49 void input_matrix( Matrix &A ) {
     50     for( int i = 0; i < n; ++i )
     51         for( int j = 0; j < n; ++j ) {
     52             scanf( "%d", &A[i][j] );
     53             A[i][j] %= p;
     54         }
     55 }
     56 void mul( const Matrix &A, const Matrix &B, Matrix &C ) { // 矩阵乘法,结果存入C
     57     Matrix tmp; tmp.clear();
     58     for( int i = 0; i < n; ++i )
     59         for( int j = 0; j < n; ++j )
     60             for( int k = 0; k < n; ++k )
     61                 tmp[i][j] = (tmp[i][j] + A[i][k] * B[k][j] % p) % p;
     62     C = tmp;
     63 }
     64 void gauss_jordan( int A[MAXN][MAXN<<1] ) { // 高斯约当消元求逆
     65     for( int i = 0; i < n; ++i ) {
     66         int r;
     67         for( r = i; r < n; ++r )
     68             if( A[r][i] ) break;
     69         if( r != i )
     70             for( int j = i; j < (n<<1); ++j )
     71                 swap( A[r][j], A[i][j] );
     72         int iv = inv( A[i][i] );
     73         for( int j = i; j < (n<<1); ++j )
     74             A[i][j] = A[i][j] * iv % p;
     75         for( int j = 0; j < n; ++j )
     76             if( j != i && A[j][i] ) {
     77                 int tmp = (p - A[j][i]) % p;
     78                 for( int k = i; k < (n<<1); ++k )
     79                     A[j][k] = (A[j][k] + A[i][k] * tmp % p) % p;
     80             }
     81     }
     82 }
     83 void inv( const Matrix &A, Matrix &B ) { // 矩阵求逆,结果存入B
     84     int tmp[MAXN][MAXN<<1] = {0};
     85     for( int i = 0; i < n; ++i )
     86         for( int j = 0; j < n; ++j )
     87             tmp[i][j] = A[i][j];
     88     for( int i = 0; i < n; ++i ) tmp[i][i+n] = 1;
     89     gauss_jordan(tmp);
     90     for( int i = 0; i < n; ++i )
     91         for( int j = 0; j < n; ++j )
     92             B[i][j] = tmp[i][j+n];
     93 }
     94 
     95 map<Matrix,int> mp;
     96 int bsgs( Matrix A, Matrix B ) { // 大步小步算法
     97     int m = sqrt(p)+1;
     98     Matrix E; E.unit();
     99     for( int i = 0; i < m; ++i ) {
    100         if( !mp.count(E) ) mp[E] = i;
    101         mul(E,A,E);
    102     }
    103     inv(E,E);
    104     for( int i = 0; i < p; i += m ) {
    105         if( mp.count(B) ) return i + mp[B];
    106         mul(B,E,B); // 注意矩阵乘法顺序
    107     }
    108     return -1;
    109 }
    110 
    111 int main() {
    112     scanf( "%d%d", &n, &p );
    113     input_matrix(A), input_matrix(B);
    114     printf( "%d
    ", bsgs(A,B) );
    115     return 0;
    116 }
  • 相关阅读:
    Encryption (hard) CodeForces
    cf 1163D Mysterious Code (字符串, dp)
    AC日记——大整数的因子 openjudge 1.6 13
    AC日记——计算2的N次方 openjudge 1.6 12
    Ac日记——大整数减法 openjudge 1.6 11
    AC日记——大整数加法 openjudge 1.6 10
    AC日记——组合数问题 落谷 P2822 noip2016day2T1
    AC日记——向量点积计算 openjudge 1.6 09
    AC日记——石头剪刀布 openjudge 1.6 08
    AC日记——有趣的跳跃 openjudge 1.6 07
  • 原文地址:https://www.cnblogs.com/mlystdcall/p/6666892.html
Copyright © 2011-2022 走看看