数学---矩阵快速幂

hdu6470，斐波那契，矩阵快速幂，递推

解决哪些问题

求解一些递推公式的第n项的时候，通过递推公式构造转移矩阵，并用矩阵快速幂可以快速得到第n项的值。特别对于n很大的时候不能循环迭代，只能用矩阵快速幂解决。
比如斐波那契数列
0 1 1 2 3 5 8...
如果要求第n项，n比较小的话，可以直接循环迭代出来。
如果n比较大，第一亿个是多少，用循环就太慢了，要用矩阵快速幂
根据矩阵乘法性质，构造转移矩阵

那么

这样问题就变成了求转移矩阵的幂。这个转移矩阵一定构造成方阵。

矩阵快速幂

#include <stdio.h>
#include <string.h>
#define MATRIX_SIZE 2

const int mod = 1e9+7;
struct Matrix //构造一个方阵
{
    int data[MATRIX_SIZE][MATRIX_SIZE];
    Matrix(){
        memset(data, 0, sizeof(data));
        for(int i=0; i<MATRIX_SIZE; i++){
            data[i][i]=1; // 初始化为单位矩阵
        }
    }
    Matrix operator + (Matrix  o)const{ //(a+b)%mod
        Matrix re; 
        for(int i=0; i<MATRIX_SIZE; i++){
            for(int j=0; j<MATRIX_SIZE; j++){
                re.data[i][j] = (this->data[i][j] + o.data[i][j]) % mod;
            }
        }
        return re;
    }
    Matrix operator * (Matrix  o)const{ //(a*b)%mod
        Matrix re;
        for(int i=0; i<MATRIX_SIZE; i++){
            for(int j=0; j<MATRIX_SIZE; j++){
                re.data[i][j] = 0;
                for(int k=0; k<MATRIX_SIZE; k++)
                    re.data[i][j] = (re.data[i][j] + ((this->data[i][k] * o.data[k][j]) % mod)) % mod;
            }
        }
        return re;
    }
    Matrix operator ^ (int n)const{ // (a^n)%mod
        Matrix re,base;
        base = *this;
        while(n){
            if(n&1)
                re = re * base;
            n>>=1;
            base = base * base;
        }
        return re;
    }
    Matrix Psum(int n)const{ //(a+a^2+a^3.....+a^n)%mod
        Matrix a,ans,pre;
        int m;
        a = *this;
        if(n==1) return a;
        m = n/2;
        pre = a.Psum(m); // a^[1,n/2] 相加
        ans = pre + (pre * (a^m)); // ans = [1,n/2] + a^(n/2)*[1,n/2]
        if(n&1)
            ans = ans + (a^n);    //n为奇数时候a^n会漏掉，补上
        return ans;
    }
    void out(){
        for(int i=0; i<MATRIX_SIZE; i++){
            for(int j=0; j<MATRIX_SIZE; j++){
                printf("%d ",data[i][j]);
            }printf("\n");
        }
    }
};

int main(){
    Matrix ma,swa;
    ma.data[0][0]=1;
    ma.data[0][1]=1;
    ma.data[1][0]=1;
    ma.data[1][1]=0;
    int n;
    ma.out();
    while(scanf("%d" , &n) , n>=1){
        swa = ma^(n-1);
        //swa.out();
        printf("is %d\n", swa.data[0][0]*1 + swa.data[0][1]*0 );
    }
    return 0;
}

构造矩阵例题

hdu6470
http://acm.hdu.edu.cn/showproblem.php?pid=6470

Farmer John有n头奶牛.
某天奶牛想要数一数有多少头奶牛,以一种特殊的方式:
第一头奶牛为1号,第二头奶牛为2号,第三头奶牛之后,假如当前奶牛是第n头,那么他的编号就是2倍的第n-2头奶牛的编号加上第n-1头奶牛的编号再加上自己当前的n的三次方为自己的编号.
现在Farmer John想知道,第n头奶牛的编号是多少,估计答案会很大,你只要输出答案对于123456789取模.

根据题意可以得到递推公式：

由于的存在，这个递推公式并不是线性的，所以要想办法把展开：

那么转移矩阵为：

取：

所以：

AC代码：

#include <stdio.h>
#include <string.h>
#include <iostream>
#define MATRIX_SIZE 6
using namespace std;
const int mod = 123456789;
struct Matrix //构造一个方阵
{
    long long data[MATRIX_SIZE][MATRIX_SIZE];
    Matrix(){
        memset(data, 0, sizeof(data));
        for(int i=0; i<MATRIX_SIZE; i++){
            data[i][i]=1; // 初始化为单位矩阵
        }
    }
    Matrix operator * (Matrix  o)const{ //(a*b)%mod
        Matrix re;
        for(int i=0; i<MATRIX_SIZE; i++){
            for(int j=0; j<MATRIX_SIZE; j++){
                re.data[i][j] = 0;
                for(int k=0; k<MATRIX_SIZE; k++)
                    re.data[i][j] = (re.data[i][j] + ((this->data[i][k] * o.data[k][j]) % mod)) % mod;
            }
        }
        return re;
    }
    Matrix operator ^ (long long n)const{ // (a^n)%mod
        Matrix re,base;
        base = *this;
        while(n){
            if(n&1)
                re = re * base;
            n>>=1;
            base = base * base;
        }
        return re;
    }
};
Matrix ma,swa;
int main(){
    int mm[6][6]={
        {0,1,0,0,0,0},
        {2,1,1,3,3,1},
        {0,0,1,3,3,1},
        {0,0,0,1,2,1},
        {0,0,0,0,1,1},
        {0,0,0,0,0,1}
    };
    for(int i=0;i<6;i++)
        for(int j=0;j<6;j++)
            ma.data[i][j]=mm[i][j];
    long long t,a;
    cin>>t;
    while(t--){
        cin>>a;
        swa = ma^(a-2);
        long long ac;
        ac= swa.data[1][0]*1;
        ac+=swa.data[1][1]*2;
        ac+=swa.data[1][2]*8;
        ac+=swa.data[1][3]*4;
        ac+=swa.data[1][4]*2;
        ac+=swa.data[1][5]*1;
        cout<<ac%mod<<endl;
    }
    return 0;
}