[模板]矩阵的快速幂

 原理：

$left( {egin{array}{*{20}{c}}
{{F_{n + 2}}}\
{{F_{n + 1}}}
end{array}} ight) = left( {egin{array}{*{20}{c}}
1&1\
1&0
end{array}} ight)left( egin{array}{l}
{F_{{ m{n + }}1}}\
{F_n}
end{array} ight)$

 记这个矩阵为A，则有

$left( {egin{array}{*{20}{c}}
{{F_{n + 1}}}\
{{F_n}}
end{array}} ight) = {A^n}left( egin{array}{l}
{F_1}\
{F_0}
end{array} ight)$

注意这种初始化方式，matrix为只有一个二维数组的结构体

struct matrix
{
    long long m[4][4];
};
matrix ans = {  1,0,0,0,
                0,1,0,0,
                0,0,1,0,
                0,0,0,1,};
matrix base = { 1,1,1,0,
                1,0,0,0,
                0,1,0,0,
                1,0,0,1,};

第一种要注意，界限的确定，要保证符合所有矩阵。

#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<algorithm>
using namespace std;
typedef long long ll;
struct mat{
    ll m[2][2];
};
ll n=10;
ll mod=10000;
mat mul(mat &A,mat &B){
    mat C={0};
    for(int i=0;i<2;i++){
        for(int j=0;j<2;j++){
            for(int k=0;k<2;k++){
                C.m[i][j]=(C.m[i][j]+A.m[i][k]*B.m[k][j])%mod;
            }
        }
    }
    return C;
}
mat pow(mat A,ll n){
    mat B={0};
    B.m[0][0]=B.m[1][1]=1;
    while(n>0){
        if(n&1) B=mul(B,A);
        A=mul(A,A);
        n>>=1;
    }
    return B;
}
int main(){
    mat A={0};
    A.m[0][0]=A.m[1][0]=A.m[0][1]=1;
    A=pow(A, n);
    printf("%lld
",A.m[1][0]);
}

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef vector<ll> vec;
typedef vector<vec> mat;
ll n=10;
const int M=10000;
mat mul(mat &A,mat &B){
    mat C(A.size(),vec(B[0].size()));
    for(int i=0;i<(int)A.size();i++){
        for(int j=0;j<(int)B[0].size();j++){
            for(int k=0;k<(int)B.size();k++){
                C[i][j]=(C[i][j]+A[i][k]*B[k][j])%M;
            }
        }
    }
    return C;
} 

mat pow(mat A,ll n){
    mat B(A.size(),vec(A.size()));
    for(int i=0;i<(int)A.size();i++){
        B[i][i]=1;//单位矩阵 
    }
    while(n>0){
        if(n&1) B=mul(B,A);
        A=mul(A,A);
        n>>=1;
    }
    return B;
}
int main(){
    mat A(2,vec(2));
    A[0][0]=1;A[0][1]=1;
    A[1][0]=1;A[1][1]=0;
    A=pow(A,n);
    printf("%lld",A[1][0]);//为什么输出这个，按照递推公式得出 
     
}

http://acm.xidian.edu.cn/problem.php?id=1145

再介绍一题http://acm.hit.edu.cn/hoj/problem/view?id=2255

是哈工大的在线oj上的一个题目，一个类似于斐波那契的题目，题目会给出a,b,p,q,s,e, 其中f(0)=a,f(1)=b,当n>=2时 f(n)=P*f(n-1)+q*f((n-2) 求它组成的数列从第s项起一直加到第e项的和是多少，这题就不能推f(n)的矩阵乘法的递推式了，得推出他的前n项和s(n)的递推式,当然如果会推斐波那契的矩阵形式的递推形式的想必这个就不难了

推理过程如下：

$Fleft( N ight) = Sleft( N ight) - Sleft( {N - 1} ight);$

又有$Fleft( N ight) = P*Fleft( {N - 1} ight) + Q*Fleft( {N - 2} ight)$

所以 $Sleft( N ight) - Sleft( {N - 1} ight) = P*left( {Sleft( {N - 1} ight) - Sleft( {N - 2} ight)} ight) + Q*left( {Sleft( {N - 2} ight) - Sleft( {N - 3} ight)} ight)$

$Sleft( N ight) = left( {P + 1} ight)*Sleft( {N - 1} ight) + left( {Q - P} ight)*Sleft( {N - 2} ight) - Q*Sleft( {N - 3} ight)$

推到这里你应该知道怎么把它化成矩阵乘法的形式了吧？注意到右边有三个S(x)项所以矩阵递推式左边也要有三项写出来就是

化简一下

所以第s项到第e项的和就是S(E)-S(s-1) 注意是S(s-1)

解题方式与上面类似

#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<algorithm>
using namespace std;
typedef long long ll;
ll mod=1000000009;
struct mat{
    ll m[4][4];
};
ll l,r,res1,res2;
mat A={ 2,0,0,-1,
        1,0,0,0,
        0,1,0,0,
        0,0,1,0};
mat mul(mat &A,mat &B){
    mat C={0};
    for(int i=0;i<4;i++){
        for(int j=0;j<4;j++){
            for(int k=0;k<4;k++){
                C.m[i][j]=(C.m[i][j]+A.m[i][k]*B.m[k][j]%mod+mod)%mod;
            }
        }
    }
    return C;
}

mat mod_pow(mat A,ll n){
    mat B={0};
    B.m[0][0]=B.m[1][1]=B.m[2][2]=B.m[3][3]=1;
    while(n>0){
        if(n&1) B=mul(B,A);
        A=mul(A,A);
        n>>=1;
    }
    return B;
}

ll res(mat &a){
    return (a.m[3][0]*6+a.m[3][1]*3+a.m[3][2]*2+a.m[3][3])%mod;
}

int main(){
    while(scanf("%lld%lld",&l,&r)!=EOF){
        if(l==0) res1=0;
        else{
            mat l2=mod_pow(A,l-1);
            res1=res(l2);
        }
        if(r==0) res2=1;
        else{
        mat rr=mod_pow(A,r);
            res2=res(rr);
        }
        ll ans=(res2-res1+mod)%mod;
        printf("%lld
",ans);
    }
}