Matrix快速幂 模板

讲解：http://www.cnblogs.com/SYCstudio/p/7211050.html

给定n*n的矩阵A，求A^k

https://www.luogu.org/problem/show?pid=3390

#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<algorithm>
using namespace std;

#define ll long long

const int maxN=201;
const ll Mod=1000000007;
const ll inf=2147483647;

int n;

class Matrix
{
public:
    ll M[maxN][maxN];
    Matrix(int x)
    {
        for (int i=0;i<n;i++)
            for (int j=0;j<n;j++)
                M[i][j]=x;
    }
    Matrix(ll Arr[maxN][maxN])
    {
        for (int i=0;i<n;i++)
            for (int j=0;j<n;j++)
                M[i][j]=Arr[i][j];
    }
    void print()
    {
        for (int i=0;i<n;i++)
        {
            for (int j=0;j<n;j++)
                cout<<M[i][j]<<' ';
            cout<<endl;
        }
    }
};

Matrix operator * (Matrix A,Matrix B)
{
    Matrix Ans(0);
    for (int i=0;i<n;i++)
        for (int j=0;j<n;j++)
            for (int k=0;k<n;k++)
                Ans.M[i][j]=(Ans.M[i][j]+A.M[i][k]*B.M[k][j]%Mod)%Mod;
    return Ans;
}

ll Arr[maxN][maxN];

ll read();
void Pow(ll Po);

int main()
{
    n=read();
    ll Po=read();
    for (int i=0;i<n;i++)
        for (int j=0;j<n;j++)
            Arr[i][j]=read();
    Pow(Po-1);//注意，这里为什么要-1呢，因为我们知道a^1是a，对于矩阵来说就是A^1是A,所以在传进去的时候先-1（相当于已经进行了一次操作），而若Po==1，则在Pow(Po-1)中不会执行循环，此时也正好是矩阵A（仔细揣摩一下）
    return 0;
}

ll read()
{
    ll x=0;
    ll k=1;
    char ch=getchar();
    while (((ch>'9')||(ch<'0'))&&(ch!='-'))
        ch=getchar();
    if (ch=='-')
    {
        k=-1;
        ch=getchar();
    }
    while ((ch>='0')&&(ch<='9'))
    {
        x=x*10+ch-48;
        ch=getchar();
    }
    return x*k;
}

void Pow(ll P)
{
    Matrix A(Arr);
    Matrix B(Arr);
    while (P!=0)
    {
        if (P&1)
            A=A*B;
        B=B*B;
        P=P>>1;
    }
    A.print();
    return;
}

例题 求斐波那契数列的第N项%MOD

#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<algorithm>
using namespace std;

#define ll long long//注意用长整形，因为有可能会爆int

const int Mod=1000000007;
const int inf=2147483647;

class Matrix//定义矩阵
{
public:
    ll M[2][2];
    Matrix()
    {
        memset(M,0,sizeof(M));
    }
    Matrix(int Arr[2][2])//定义两个方便的矩阵初始化
    {
        for (int i=0;i<2;i++)
            for (int j=0;j<2;j++)
                M[i][j]=Arr[i][j];
    }
};

Matrix operator * (Matrix A,Matrix B)//重载乘号操作
{
    Matrix Ans;
    for (int i=0;i<2;i++)
        for (int j=0;j<2;j++)
            for (int k=0;k<2;k++)
                Ans.M[i][j]=(Ans.M[i][j]+A.M[i][k]*B.M[k][j]%Mod)%Mod;
    return Ans;
}

ll n;

int main()
{
    cin>>n;
    if (n<=2)
    {
        cout<<1<<endl;
        return 0;
    }
    n=n-2;
    int a[2][2]={{1,1},{0,0}};//初始矩阵
    int b[2][2]={{1,1},{1,0}};//即上文的T
    Matrix A(a);
    Matrix B(b);
    while (n!=0)//快速幂
    {
        if (n&1)
            A=A*B;
        B=B*B;
        n=n>>1;
    }
    cout<<A.M[0][0]<<endl;
    return 0;
}

例题 求广义斐波那契数列的第N项%MOD

#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<algorithm>
using namespace std;

#define ll long long

const int inf=2147483647;

ll n,Mod,p,q,A1,A2;

class Matrix//定义矩阵结构体
{
public:
    ll M[2][2];
    Matrix()
    {
        memset(M,0,sizeof(M));
    }
    Matrix(ll Arr[2][2])
    {
        for (int i=0;i<2;i++)
            for (int j=0;j<2;j++)
                M[i][j]=Arr[i][j];
    }
};

Matrix operator * (Matrix A,Matrix B)//重载乘法操作
{
    Matrix Ans;
    for (int i=0;i<2;i++)
        for (int j=0;j<2;j++)
            for (int k=0;k<2;k++)
                Ans.M[i][j]=(Ans.M[i][j]+A.M[i][k]*B.M[k][j]%Mod)%Mod;
    return Ans;
}

int main()
{
    cin>>p>>q>>A1>>A2>>n>>Mod;
    if (n==1)//1和2的情况特殊处理（虽然我也不知道是否有特殊点）
    {
        cout<<A1<<endl;
        return 0;
    }
    if (n==2)
    {
        cout<<A2<<endl;
        return 0;
    }
    n=n-2;
    ll a[2][2]={{A2,A1},{0,0}};//初始矩阵
    ll b[2][2]={{p,1},{q,0}};
    Matrix A(a);
    Matrix B(b);
    while (n!=0)//快速幂
    {
        if (n&1)
            A=A*B;
        B=B*B;
        n=n>>1;
    }
    cout<<A.M[0][0]<<endl;
    return 0;
}

小 X 在上完生物课后对细胞的分裂产生了浓厚的兴趣。于是他决定做实验并

观察细胞分裂的规律。

他选取了一种特别的细胞，每天每个该细胞可以分裂出 x − 1 个新的细胞。

小 X 决定第 i 天向培养皿中加入 i 个细胞（在实验开始前培养皿中无细胞）。

现在他想知道第 n 天培养皿中总共会有多少个细胞。

由于细胞总数可能很多，你只要告诉他总数对 w 取模的值即可。

#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<algorithm>
#include<cmath>
using namespace std;

#define ll long long

ll n,X,Mod;

class Matrix//定义一个矩阵结构体
{
public:
    ll M[3][3];
    Matrix()//初始化0
    {
        for (int i=0;i<3;i++)
            for (int j=0;j<3;j++)
                M[i][j]=0;
    }
    Matrix(ll Arr[3][3])//用数组来初始化
    {
        for (int i=0;i<3;i++)
            for (int j=0;j<3;j++)
                M[i][j]=Arr[i][j];
    }
};

Matrix operator * (Matrix A,Matrix B)//重载乘法运算符
{
    Matrix Ans;
    for (int i=0;i<3;i++)
        for (int j=0;j<3;j++)
            for (int k=0;k<3;k++)
                Ans.M[i][j]=(Ans.M[i][j]+A.M[i][k]*B.M[k][j]%Mod)%Mod;
    return Ans;
}

const int inf=2147483647;

ll Pow();

int main()//无比简单的主函数
{
    cin>>n>>X>>Mod;
    cout<<Pow()<<endl;
    return 0;
}

ll Pow()//矩阵快速幂
{
    ll a[3][3]={{0,1,1},{0,0,0},{0,0,0}};//初始状态
    ll b[3][3]={{X,0,0},{1,1,0},{0,1,1}};//用来使状态发生转移的矩阵，即文中提到的T
    Matrix A(a);
    Matrix B(b);
    while (n!=0)
    {
        if (n&1)
        {
            A=A*B;//注意这里的乘法顺序，矩阵乘法不满足交换律
        }
        B=B*B;
        n=n>>1;
    }
    return A.M[0][0];//根据我们的定义，最后的值就在A.M[0][0]
}