题解 UVA11149 【Power of Matrix】

倍增+矩乘(O(T*n^3logk))

然而因为看到很多人的写法都是 (O(T*n^3log^2k)) 的做法才发的？

先说一下做法吧，注意到矩阵乘法具有分配律，大概就是：

((A+B)*C = A*C +B*C)

至于这里的加法是定义在两个矩阵大小相同的时候每个相应的位置上的元素和。

至于上述(log^2k)的做法简单提一下算了，大致思路就是将(k)分成两瓣((k/2,k/2+k/2))，然后对于左(k/2)的部分递归求，(k/2+k/2)则利用已经求出的(k/2)来计算，即(s[k] = s[k/2]*c[k/2]+s[k/2])然而递归(logk)层，每层都需要(logk)地去求一遍快速幂。

然后为了方便倍增的方法的讲解，还是先讲一下复杂度 (O(T*n^3sqrt k))的做法吧

具体就是把原式变成：

[S(k) = A^1+...+A^{sqrt k} + A^{sqrt k} * (A^1 ...+A^{sqrt k}) +A^{2*sqrt k} * (...)+... ]

然而这样做复杂度是 (O(T*n^3sqrt k))

我们发现这个思想貌似也可以倍增的去处理欸欸欸？

比如：

[S(k) = A^1 + A^1 * (A^1) + A^2 * (A^1 +A^2) +A^4 * (A^1+A^2+...+A^4)+ ]

然后我们就可以处理出两个数组，表示

[c[i] = A^{2^i},s[i]= (A^1+A^2+...+A^{2^i}) ]

然后我们可以发现其的处理规律：

[c[i] = c[i-1]*c[i-1],s[i]=s[i-1]+s[i-1]*c[i-1] ]

然后处理到最后一截还剩一点点其不为(2^i)

大概就是：(k-2^i)这一部分。

我们可以类似于倍增求 (LCA) 的部分，反向枚举，从某一次幂枚举到低次幂。

先处理出一个当前处理到那个位置 (now)，然后发现：比如现在想求(S^{1024+1023})（S表示和）

我们现在已经处理出来了(A^0,A^1,A^2,A^4...A^{1024},)和(S^0,S^1,S^2,S^4...S^{1024})

然后现在想要求出 (S^{1024+1023})

发现(S^{1024})已经处理出来了，然后我们可以考虑枚举(2)的某次幂，从高到底枚举，比如先枚举10，发现不行，然后枚举9，发现(now(1024)+2^9<=1024+1023)

然后就令(ans = ans(S^{1024}) + S^{512}*now(A^{1024}))

然后：(now = now * A^{512})

然后大概长这样：

drep( i, p, 0 ) { //orz表示当前的c矩阵，len表示当前的位置，now表示2^i次幂，
		if( len + now <= k ) 
        	ans = Add( ans, Mul( orz, s[i] ) ), len += now, orz = Mul( orz, c[i] ) ; 
		now /= 2 ;
	}

类似的处理下去，不难发现，剩余部分处理的复杂度也是(O(logk))

总复杂度(O(T*n^3logk))

详见代码：

#include<bits/stdc++.h>
using namespace std;
int read() {
    char cc = getchar(); int cn = 0, flus = 1;
    while(cc < '0' || cc > '9') {  if( cc == '-' ) flus = -flus;  cc = getchar();  }
    while(cc >= '0' && cc <= '9')  cn = cn * 10 + cc - '0', cc = getchar();
    return cn * flus;
}
const int mod = 10 ; 
const int N = 43 ; 
#define rep( i, s, t ) for( register int i = s; i <= t; ++ i )
#define drep( i, t, s ) for( register int i = t; i >= s; -- i )
int n, k ;
struct Mat {
    int m[N][N] ; 
    void init1() {
        rep( i, 1, n ) rep( j, 1, n ) m[i][j] = 0 ;
    }
} c[60], s[60]; 
Mat Mul( Mat a, Mat b ) {
    Mat ans ; ans.init1() ; 
    rep( i, 1, n ) rep( j, 1, n ) rep( k, 1, n ) 
        ans.m[i][j] += a.m[i][k] * b.m[k][j], ans.m[i][j] %= mod ;  
    return ans ; 
}
Mat Add( Mat a, Mat b ) {
    Mat ans ; ans.init1() ; 
    rep( i, 1, n ) rep( j, 1, n ) 
    ans.m[i][j] = a.m[i][j] + b.m[i][j], ans.m[i][j] %= mod ; 
    return ans ; 
}
void print( Mat a ) {
    rep( i, 1, n ) {
        rep( j, 1, n ) {
            printf("%lld", a.m[i][j] % mod ) ;
            if( j != n ) printf(" ");
        }
        puts("");
    }
    puts("");
}
void solve() {
    int now = 2, p = 1, len ; Mat orz, ans ;
    while( now <= k ) c[p] = Mul( c[p - 1], c[p - 1] ), s[p] = Add( s[p - 1], Mul( s[p - 1], c[p - 1] ) ), now *= 2, ++ p ;
    
    now /= 2, len = now, -- p, orz = c[p], ans = s[p] ;
    
    drep( i, p, 0 ) {
        if( len + now <= k ) ans = Add( ans, Mul( orz, s[i] ) ), len += now, orz = Mul( orz, c[i] ) ; 
        now /= 2;
    }
    print(ans) ; 
}
signed main()
{
    while(n = read()) {
        k = read();
        rep( i, 1, n ) rep( j, 1, n ) s[0].m[i][j] = c[0].m[i][j] = read() % mod ; 
        solve() ;
    }
    k = read() ; 
    return 0;
}

貌似(20ms quad QwQ)