POJ 3233 Matrix Power Series(矩陣+二分)

題目大意：求由矩陣 A構成的矩陣 S = A + A^2 + A^3 + … + A^k。k的取值範圍是：10^9資料很大，應該二分。

對於一個k來說，s(k) = (1+A^(k/2)) *( A+A^2+……+A^(k/2))。如果k為奇數的話需要加上A^(k/2 + 1)。

所以二分求和，複雜度就降下來了，當然還得用到矩陣快速冪。

Description

Input

Output

Sample Input

2 2 4
0 1
1 1

Sample Output

1 2
2 3

#include <algorithm>
#include <iostream>
#include <stdlib.h>
#include <string.h>
#include <iomanip>
#include <stdio.h>
#include <string>
#include <queue>
#include <cmath>
#include <stack>
#include <map>
#include <set>
#define eps 1e-10
///#define M 1000100
#define LL __int64
///#define LL long long
///#define INF 0x7ffffff
#define INF 0x3f3f3f3f
#define PI 3.1415926535898
#define zero(x) ((fabs(x)<eps)?0:x)

///#define mod 9973

int mod;
const int maxn = 2010;

using namespace std;

struct matrix
{
    int f[40][40];
};

matrix mul(matrix a, matrix b, int n)///矩陣乘法
{
    matrix c;
    memset(c.f, 0, sizeof(c.f));
    for(int i = 0; i < n; i++)
    {
        for(int j = 0; j < n; j++)
        {
            for(int k = 0; k < n; k++) c.f[i][j] += a.f[i][k]*b.f[k][j];
            c.f[i][j] %= mod;
        }
    }
    return c;
}

matrix pow_mod(matrix a, int b, int n)///矩陣快速冪
{
    matrix s;
    memset(s.f, 0 , sizeof(s.f));
    for(int i = 0; i < n; i++) s.f[i][i] = 1;
    while(b)
    {
        if(b&1) s = mul(s, a, n);
        a = mul(a, a, n);
        b >>= 1;
    }
    return s;
}

matrix Add(matrix a,matrix b, int n)  ///矩陣加法
{
    matrix c;
    for(int i = 0; i < n; i++)
    {
        for(int j = 0; j < n; j++)
        {
            c.f[i][j] = a.f[i][j]+b.f[i][j];
            c.f[i][j] %= mod;
        }
    }
    return c;
}

matrix Matrix_Sum(matrix a, int k, int n)
{
    if(k == 1) return a;
    matrix dx,dy;
    dx = Matrix_Sum(a, k/2, n);///二分，遞迴；
    if(k&1)
    {
        dy = pow_mod(a, k/2+1, n);
        dx = Add(dx, mul(dx, dy, n), n);
        dx = Add(dy,dx, n);
    }
    else 
    {
        dy = pow_mod(a, k/2, n);
        dx = Add(dx,mul(dx, dy, n), n);
    }
    return dx;
}

int main()
{
    int n, k;
    while(cin >>n>>k>>mod)
    {
        matrix c;
        for(int i = 0; i < n; i++)
        {
            for(int j = 0; j < n; j++)
            {
                scanf("%d",&c.f[i][j]);
                c.f[i][j] %= mod;
            }
        }
        matrix d = Matrix_Sum(c, k, n);
        for(int i = 0; i < n; i++)
        {
            for(int j = 0; j < n-1; j++) cout<<d.f[i][j]<<" ";
            cout<<d.f[i][n-1]<<endl;
        }
    }
    return 0;
}