矩阵快速幂

一、矩阵

1 运算

1.1加减法

要求：两矩阵行列数完全一样。

1.2乘除法

要求：第一个矩阵的列数等于第二个矩阵的行数。
一个a行b列乘b行c列 = a行c列
乘法不满足交换律！
但是满足结合律：$a\times b\times c=a\times (b\times c)$

1.3 数乘

一个矩阵乘以一个数，等于矩阵中的每个数乘以这个数。

1.4 转置

原先的第一行变成第一列，以此类推。

2 单位元

在一个运算系统中，如果某个元素x和其他元素y做运算，结果依旧是y，那么x就称为这个运算系统中的“单位元”。

2.1 单位矩阵

其主对角线上的元素都是1，其余都是0。
例：
$$\left[\begin{array}{ccc}1& 0& 0\\ 0& 1& 0\\ 0& 0& 1\end{array}\right]$$

二、封装

使用重载运算符和结构体封装会使代码更简洁，同时易于更改。以下是封装代码；

#include<bits/stdc++.h>
#define ll long long
#define bug printf("---OK---")
#define pa printf("A: ")
#define pr printf("\n")
#define pi acos(-1.0)
using namespace std;
struct Mat{
	ll a[105][105];
	ll r,c;
	Mat(ll _r=0,ll _c=0){
		r=_r,c=_c;
		memset(a,0,sizeof a);
		if(c==0){
			c=r;//构造一个方阵
		}
	}
	void unit(){//赋单位矩阵
        memset(a,0,sizeof a);
        for(int i=1;i<=r;i++){
        	a[i][i]=1;
		}
    }
	Mat operator+(const Mat t)const{//加法
		Mat ans(r,c);
		for(int i=1;i<=r;i++){
			for(int j=1;j<=c;j++){
				ans.a[i][j]=a[i][j]+t.a[i][j];
			}
		}
		return ans;
	}
	Mat operator-(const Mat t)const{//减法
		Mat ans(r,c);
		for(int i=1;i<=r;i++){
			for(int j=1;j<=c;j++){
				ans.a[i][j]=a[i][j]-t.a[i][j];
			}
		}
		return ans;
	}
	Mat operator*(const Mat t)const{//乘法
		Mat ans(r,t.r);
		ll n=r,m=t.c;
		for(int i=1;i<=n;i++){
			for(int j=1;j<=m;j++){
				for(int k=1;k<=c;k++){
					ans.a[i][j]+=a[i][k]*t.a[j][k];
				}
			}
		}
		return ans;
	}
	Mat operator%(const ll t)const{//取余
		Mat ans=*this;
		for(int i=1;i<=r;i++){
			for(int j=1;j<=c;j++){
				ans.a[i][j]%=t;
			}
		}
		return ans;
	}
	Mat qpow(ll b,ll p){
        Mat ans(r,c),a= *this;
        ans.unit();
        while(b){
            if(b&1)
                ans=ans*a%p;
            a=a*a%p;
            b>>=1;
        }
        return ans%p;
    }
};
int main(){
	return 0;
}

三、实践

1 斐波那契数列

1.1 题目描述

斐波那契数列是满足如下性质的一个数列：

$$F_n=\{\begin{array}{r}1(n\le 2)\\ F_{n-1}+F_{n-2}(n\ge 3)\end{array}$$

请你求出 $F_n\mathrm{mod}10^9+7$ 的值。

1.2 思路

首先是原始矩阵 $a$ ：
$$\left[\begin{array}{cc}1& 1\\ 0& 0\end{array}\right]$$
然后要找到满足斐波那契数列的单位矩阵 $b$ ：
$$\left[\begin{array}{cc}1& 1\\ 1& 0\end{array}\right]$$
每次使用矩阵快速幂即可。
注意：因为数列前两个都是 $1$ ，所以指数为 $n-2$ ！
最后输出 $a[1][1]$ 。

1.3 代码

#include<bits/stdc++.h>
#define ll long long
#define bug printf("---OK---")
#define pa printf("A: ")
#define pr printf("\n")
#define pi acos(-1.0)
using namespace std;
ll n;
const ll mod=1e9+7;
struct Mat{
	ll a[105][105];
	ll r,c;
	void scan(){//输出
		for(int i=1;i<=r;i++){
			for(int j=1;j<=c;i++){
				scanf("%lld",&a[i][j]);
			}
		}
	}
	Mat(ll _r=0,ll _c=0){
		r=_r,c=_c;
		memset(a,0,sizeof a);
		if(c==0){
			c=r;//构造一个方阵
		}
	}
	void unit(){//赋单位矩阵
        memset(a,0,sizeof a);
        for(int i=1;i<=r;i++){
        	a[i][i]=1;
		}
    }
	Mat operator+(const Mat t)const{//加法
		Mat ans(r,c);
		for(int i=1;i<=r;i++){
			for(int j=1;j<=c;j++){
				ans.a[i][j]=a[i][j]+t.a[i][j];
			}
		}
		return ans;
	}
	Mat operator-(const Mat t)const{//减法
		Mat ans(r,c);
		for(int i=1;i<=r;i++){
			for(int j=1;j<=c;j++){
				ans.a[i][j]=a[i][j]-t.a[i][j];
			}
		}
		return ans;
	}
	Mat operator*(const Mat t)const{//乘法
		Mat ans(r,t.r);
		ll n=r,m=t.c;
		for(int i=1;i<=n;i++){
			for(int j=1;j<=m;j++){
				for(int k=1;k<=c;k++){
					ans.a[i][j]+=a[i][k]*t.a[j][k];
				}
			}
		}
		return ans;
	}
	Mat operator%(const ll t)const{//取余
		Mat ans=*this;
		for(int i=1;i<=r;i++){
			for(int j=1;j<=c;j++){
				ans.a[i][j]%=t;
			}
		}
		return ans;
	}
	Mat qpow(ll b,ll p){
        Mat ans(r,c),a= *this;
        ans.unit();
        while(b){
            if(b&1)
                ans=ans*a%p;
            a=a*a%p;
            b>>=1;
        }
        return ans%p;
    }
	void print(){//输出
		for(int i=1;i<=r;i++){
			for(int j=1;j<=c;i++){
				printf("%lld ",a[i][j]);
			}
			pr;
		}
	}
};
int main(){
	cin>>n;
	if(n<=2){
		cout<<"1";
		return 0;
	}
	Mat a(2,2);
	a.a[1][1]=a.a[1][2]=1;
	Mat b(2,2);
	b.a[1][1]=b.a[1][2]=b.a[2][1]=1;
	b.a[2][2]=0;
	a=a*b.qpow(n-2,mod)%mod;
	cout<<a.a[1][1];
	return 0;
}