[2016-02-04][矩阵快速幂]

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

#define MOD 10000

const int mat_size = 2;

struct matrix{

        long long a[mat_size][mat_size];

};

matrix matrixE;//单位矩阵

void inimatrixE(){

        for(int i = 0;i < mat_size;++i)

                for(int j = 0;j < mat_size;++j)

                        matrixE.a[i][j] = i == j ? 1 : 0;

}

void mat_mult(matrix & a,matrix & b,matrix & c,int m,int n,int s,long long mod){

        //a == m*n      b == n*s        c == m*s

        memset(c.a,0,sizeof(c.a));

        for(int i = 0 ;i < m ; ++i)

                for(int k = 0;k < n ; ++k)

                        for(int j = 0 ; j < s ; ++j){

                                c.a[i][j] = (c.a[i][j] + a.a[i][k]*b.a[k][j]) % mod;

//如果结果不会超过longlong范围,那么取模运算可以放在第二个for内(第3个for 外面)

                        }

}

void quick_matrix_pow(matrix & a,int p,matrix & res){

        //inimatrixE();

        memset(res.a,0,sizeof(res.a));

        matrix tmp = a,tmpres;

        res = matrixE;

        while(p >= 1){

                if(p & 1){

                        mat_mult(tmp,res,tmpres,mat_size,mat_size,mat_size,MOD);

                        res = tmpres;

                }

                p >>= 1;

                mat_mult(tmp,tmp,tmpres,mat_size,mat_size,mat_size,MOD);

                tmp = tmpres;

        }

}