矩阵乘法略解

矩阵乘法

对于一个 n×mn\times m 的矩阵 AA 和一个 m×km\times k 的矩阵 BB，有且仅有一个 n×kn\times k 的矩阵 CC 使A×B=CA\times B=C且对于 ∀i∈[1,n] & ∀j∈[1,k]\forall i\in [1,n]\ \&\ \forall j\in [1,k]，有 Ci,j=∑k=1mAik×BkjC_{i,j}=\sum^{m}_{k=1}{A_{ik}\times B_{kj}}这就是矩阵乘法。比如：
[123321]×[112233]=[14141010]\left[\begin{array}{ccc}1&2&3\\3&2&1\\\end{array}\right]\times\left[\begin{array}{ccc}1&1\\2&2\\3&3\\\end{array}\right]=\left[\begin{array}{ccc}14&14\\10&10\\\end{array}\right]

题目描述 loj10220\text{loj10220}

大家都知道 Fibonacci 数列吧，fn=fn−1+fn−2f_n=f_{n-1}+f_{n-2}。

现在问题很简单，输入 nn 和 mm，求 fnmod  mf_n\mod m。

输入格式

输入 n,mn,m。

输出格式

输出 fnmod  mf_n\mod m。

样例输入

5 1000

样例输出

5

数据范围与提示

对于 100%100\% 的数据，1≤n≤2×109,1≤m≤109+101\leq n\leq2\times10^9,1\leq m\leq 10^9+10 。

Solution 10220\text{Solution 10220}

考虑使用矩阵快速幂优化时间复杂度。
若矩阵 AA 使[fi−2fi−1]×A=[fi−1fi]\left[ \begin{array}{ccc} f_{i-2}&f_{i-1} \end{array} \right] \times A=\left[ \begin{array}{ccc} f_{i-1} & f_i \end{array} \right]则
[fnfn−1]=[f1f2]×An−1\left[ \begin{array}{ccc} f_n & f_{n-1} \end{array} \right]=\left[ \begin{array}{ccc} f_1&f_2 \end{array} \right]\times A^{n-1}而An−1A^{n-1}可以使用快速幂求出。
下面推导 AA。 由矩阵乘法知 AA 是一个 2×22\times 2 的矩阵，且
fi−1=A11×fi−2+A21×fi−1f_{i-1}=A_{11}\times f_{i-2}+A_{21}\times f_{i-1}∴A11=0,A21=1,\therefore A_{11}=0,A_{21}=1,
同理得 A12=1,A22=1,A_{12}=1,A_{22}=1,
∴A=[0111].\therefore A=\left[ \begin{array}{ccc} 0 & 1\\ 1 & 1 \end{array} \right].

#include<cstdio>
#include<cstdlib>
#include<cstring>

#define reg register

typedef long long ll;

int n,m;

struct node{
	ll a[3][3];
	int x,y;
	node(){
		x=y=0;
		memset(a,0,sizeof(a));
	}
	void pt(){
		printf("%lld",a[1][2]);
	}
}s,t,ans;
node T(node a,node b){
	node c;c.x=a.x;c.y=b.y;
	for(reg int i=1;i<=c.x;++i)
		for(reg int j=1;j<=c.y;++j)
			for(reg int k=1;k<=a.y;++k)
				c.a[i][j]=(c.a[i][j]+a.a[i][k]*b.a[k][j])%m;
	return c;
}
node QT(node a,int b){
	if(b<=1) return a;
	node c=QT(a,b/2);
	if(b%2)
		return T(T(c,c),a);
	else
		return T(c,c);
}

int main(){
	scanf("%d%d",&n,&m);
	s.x=1;s.y=2;s.a[1][1]=s.a[1][2]=1;
	t.x=2;t.y=2;t.a[1][1]=0;t.a[1][2]=t.a[2][1]=t.a[2][2]=1;
	ans=T(s,QT(t,n-2));
	ans.pt();
}

题目描述 loj10221\text{loj10221}

大家都知道 Fibonacci 数列吧，fn=fn−1+fn−2f_n=f_{n-1}+f_{n-2}。

求sn=∑k=1nfkmod&ThinSpace;&ThinSpace;ms_n=\sum^{n}_{k=1}{f_k} \mod m的值。

样例输入

5 1000

样例输出

12

数据范围与提示

对于 100%100\% 的数据，1≤n≤2×109,1≤m≤109+101\leq n\leq2\times10^9,1\leq m\leq 10^9+10 。

Solution 10221\text{Solution 10221}

对于Fibonacci数列前 nn 项和 sns_n，显然有 sn=sn−1+fns_n=s_{n-1}+f_n。类似的，我们可以构造一矩阵 BB，使 [fi−2fi−1si−1]×B=[fi−1fisi]\left[ \begin{array}{ccc} f_{i-2}&amp;f_{i-1}&amp;s_{i-1}\\ \end{array} \right] \times B=\left[ \begin{array}{ccc} f_{i-1}&amp;f_i&amp; s_i\\ \end{array} \right]，就能完成此题。

#include<cstdio>
#include<cstdlib>
#include<cstring>

#define reg register

typedef long long ll;

int n,m;

struct node{
	ll a[4][4];
	int x,y;
	node(){
		x=y=0;
		memset(a,0,sizeof(a));
	}
	void pt(){
		printf("%lld",a[1][3]);
	}
}s,t,ans;
node T(node a,node b){
	node c;c.x=a.x;c.y=b.y;
	for(reg int i=1;i<=c.x;++i)
		for(reg int j=1;j<=c.y;++j)
			for(reg int k=1;k<=a.y;++k)
				c.a[i][j]=(c.a[i][j]+a.a[i][k]*b.a[k][j])%m;
	return c;
}
node QT(node a,int b){
	if(b<=1) return a;
	node c=QT(a,b/2);
	if(b%2)
		return T(T(c,c),a);
	else
		return T(c,c);
}

int main(){
	scanf("%d%d",&n,&m);
	s.x=1;s.y=3;s.a[1][1]=s.a[1][2]=1;s.a[1][3]=2;
	t.x=3;t.y=3;t.a[1][1]=t.a[3][1]=t.a[3][2]=0;t.a[1][2]=t.a[1][3]=t.a[2][1]=t.a[2][2]=t.a[2][3]=t.a[3][3]=1;;
	ans=T(s,QT(t,n-2));
	ans.pt();
}

题目描述 loj10222\text{loj10222}

大家都知道 Fibonacci 数列吧，fn=fn−1+fn−2f_n=f_{n-1}+f_{n-2}。

求 Tn=∑k=1nk⋅fkmod&ThinSpace;&ThinSpace;mT_n=\sum^{n}_{k=1}{k·f_k} \mod m 的值。

样例输入

5 5

样例输出

4

数据范围与提示

对于 100%100\% 的数据，1≤n,m≤231−11\leq n,m\leq 2^{31}-1 。

Solution 10222\text{Solution 10222}

考虑构造一 1×41\times 4 矩阵
Ai=[fi−2fi−1Ti−1i−1]A_i=\left[ \begin{array}{ccc} f_{i-2}&amp;f_{i-1}&amp;T_{i-1}&amp;i-1 \end{array} \right]
使
Ai×B=Ai+1A_i\times B=A_{i+1}
其中 BB 是一个 4×44\times 4 的矩阵。同上文，我们有
fi−1=B11⋅fi−2+B21⋅fi−1+B31⋅Ti−1+B41⋅(i−1)f_{i-1}=B_{11}·f_{i-2}+B_{21}·f_{i-1}+B_{31}·T_{i-1}+B_{41}·(i-1)fi=B12⋅fi−2+B22⋅fi−1+B32⋅Ti−1+B42⋅(i−1)f_i=B_{12}·f_{i-2}+B_{22}·f_{i-1}+B_{32}·T_{i-1}+B_{42}·(i-1)Ti= ⋅⋅⋅T_i=\ ···
那么，TiT_i 怎么用形如上式的表达式表示呢？
事实上，我们无法用形如上式的表达式表示 Ti=Ti−1+i⋅fiT_i=T_{i-1}+i·f_i，因为我们无法用方阵 BB 的若干次幂的形式表示 ii。

考虑如何将 TiT_i 表示成一次递推式。
(1)Tn=∑k=1nk⋅fkT_n=\sum^{n}_{k=1}{k·f_k}\tag 1(2)nsn=∑k=1nn⋅fkns_n=\sum^{n}_{k=1}{n·f_k}\tag 2
(2)−(1)(2)-(1) 得
(3)nsn−Tn=∑k=1n(n−k)⋅fkns_n-T_n=\sum^{n}_{k=1}(n-k)·f_k\tag3
(3)+sn(3)+s_n 得
nsn−Tn+sn=∑k=1n(n−k)⋅fk+sn=∑k=1n(n−k)⋅fk+∑k=1nfk=∑k=1n(n+1−k)⋅fk=∑k=1n+1(n+1−k)⋅fk=(n+1)sn+1−Tn+1\begin{aligned} ns_n-T_n+s_n&amp;=\sum^{n}_{k=1}{(n-k)·f_k}+s_n\\ &amp;=\sum^{n}_{k=1}{(n-k)·f_k}+\sum^n_{k=1}{f_k}\\ &amp;=\sum^{n}_{k=1}{(n+1-k)·f_k}\\ &amp;=\sum^{n+1}_{k=1}{(n+1-k)·f_k}\\ &amp;=(n+1)s_{n+1}-T_{n+1} \end{aligned}
设 pk=ksk−Tkp_k=ks_k-T_k，于是我们得到 pp 的递推关系式
pk=pk−1+sk−1p_k=p_{k-1}+s_{k-1}所以我们可以构造矩阵
Ai′=[pisififi−1]A_i&#x27;=\left[ \begin{array}{ccc} p_i&amp;s_i&amp;f_i&amp;f_{i-1} \end{array} \right]B′=[1100011000110010]B&#x27;=\left[ \begin{array}{ccc} 1&amp;1&amp;0&amp;0\\ 0&amp;1&amp;1&amp;0\\ 0&amp;0&amp;1&amp;1\\ 0&amp;0&amp;1&amp;0 \end{array} \right]
使
Ai′×B′=Ai+1′A&#x27;_i\times B&#x27;=A&#x27;_{i+1}

#include<cstdio>
#include<cstdlib>
#include<cstring>

#define reg register

typedef long long ll;

int n,m;

struct node{
	ll a[5][5];
	int x,y;
	node(){
		x=y=0;
		memset(a,0,sizeof(a));
	}
	void pt(){
		printf("%lld",(a[1][2]*n%m-a[1][1]+m)%m);
	}
}s,t,ans;
node T(node a,node b){
	node c;c.x=a.x;c.y=b.y;
	for(reg int i=1;i<=c.x;++i)
		for(reg int j=1;j<=c.y;++j)
			for(reg int k=1;k<=a.y;++k)
				c.a[i][j]=(c.a[i][j]+a.a[i][k]*b.a[k][j])%m;
	return c;
}
node QT(node a,int b){
	if(b<=1) return a;
	node c=QT(a,b/2);
	if(b%2)
		return T(T(c,c),a);
	else
		return T(c,c);
}

int main(){
	scanf("%d%d",&n,&m);
	s.x=1;s.y=4;s.a[1][2]=s.a[1][3]=s.a[1][4]=1;
	t.x=4;t.y=4;t.a[1][1]=t.a[2][1]=t.a[2][2]=t.a[3][2]=t.a[3][3]=t.a[4][3]=t.a[3][4]=1;
	ans=T(s,QT(t,n-1));
	ans.pt();
}

题目描述 loj10225\text{loj10225}

原题来自：SCOI 2009

Windy 在有向图中迷路了。 该有向图有 nn 个节点，Windy 从节点 00 出发，他必须恰好在 TT 时刻到达节点 。

现在给出该有向图，你能告诉 Windy 总共有多少种不同的路径吗？

注意：Windy 不能在某个节点逗留，且通过某有向边的时间严格为给定的时间。

输入格式

第一行包含两个整数 nn，TT；
接下来有 nn 行，每行一个长度为 nn 的字符串。第 ii 行第 jj 列为 00 表示从节点 ii 到节点 jj 没有边，为 11 到 99 表示从节点 ii 到节点 jj 需要耗费的时间。

输出格式

包含一个整数，可能的路径数，这个数可能很大，只需输出这个数除以 20092009 的余数。

样例输入 1

2 2
11
00

样例输出 1

1

样例输入 2

5 30
12045
07105
47805
12024
12345

样例输出 2

852

数据范围与提示

对于 100%100\% 的数据，满足 2≤n≤10,1≤T≤1092\leq n\leq 10,1\leq T\leq 10^9。

Solution 10225\text{Solution 10225}

nn 的规模和边权都很小，考虑拆点，则时间复杂度为 O(nT)O(n^T)。
如何把 TT 优化成 log⁡T\log T 呢？观察输入格式，我们很自然地想到矩阵，尝试使用矩阵快速幂优化。
设从 ii 到 jj 的边权为 uiju_{ij}。考虑当边权为 00 或 11 时，从 ii 到 jj 路径长为 TT 的路径数
sijT=∑k=1n(∑u=T)s_{ijT}=\sum^{n}_{k=1}{(\sum u=T)}
如果构造矩阵U={uij∣i,j∈[1,n]}U=\{u_{ij}|i,j\in [1,n]\}S=UTS=U^T
则Sij=sijT (∀i,j)S_{ij}=s_{ijT}\ (\forall i,j)

#include<cstdio>
#include<cstdlib>
#include<cstring>

#define reg register

typedef long long ll;

int n,t;
const int m=2009;
char str[20][20];

struct node{
	ll a[110][110];
	int x,y;
	node(){
		x=y=0;
		memset(a,0,sizeof(a));
	}
}s,ans;
node T(node a,node b){
	node c;c.x=a.x;c.y=b.y;
	for(reg int i=1;i<=c.x;++i)
		for(reg int j=1;j<=c.y;++j)
			for(reg int k=1;k<=a.y;++k)
				c.a[i][j]=(c.a[i][j]+a.a[i][k]*b.a[k][j])%m;
	return c;
}
node QT(node a,int b){
	node c;c.x=n;c.y=n;
	for(reg int i=1;i<=n;++i)
		c.a[i][i]=1;
	while(b){
		if(b&1) c=T(c,a);
		a=T(a,a);b>>=1;
	}
	return c;
}

int main(){
	scanf("%d%d",&n,&t);
	for(reg int i=1;i<=n;++i){
		scanf("%s",str[i]+1);
		for(reg int j=1;j<=n;++j)
			str[i][j]-='0';
	}
	for(reg int i=1;i<=n;++i){
		for(reg int j=1;j<=8;++j)
			s.a[(i-1)*9+j][(i-1)*9+j+1]=1;
		for(reg int j=1;j<=n;++j)
			if(str[i][j])
				s.a[(i-1)*9+str[i][j]][(j-1)*9+1]=1;
	}
	n*=9;s.x=s.y=n;
	ans=QT(s,t);
	printf("%d",ans.a[1][n-8]);
}

题目描述 情书

求
∑k=1nfk+∑k=1nck(mod&ThinSpace;&ThinSpace;m)\sum_{k=1}^{n}f_k+\sum_{k=1}^{n}{c^k}\quad(\mod m)
的值，其中 ∀i&gt;2\forall i&gt;2 有
fi=a⋅fi−1+b⋅fi−2f_i=a·f_{i-1}+b·f_{i-2}

输入描述

一行 77 个整数，分别为 a,b,c,n,m,f1,f2a,b,c,n,m,f_1,f_2，相邻两数以一个空格隔开。

输出描述

一行，一个整数，表示答案。

输入样例

    1 1 1 3 5 1 1\ \ \ \ \text{1 1 1 3 5 1 1}

输出样例

    2\ \ \ \ 2

数据规模与约定

对于 100%100\% 的数据，有 0≤a,b,c,n,m,fi≤2310\leq a,b,c,n,m,f_i\leq2^{31}，其中 i∈[1,n]i\in[1,n]。

Solution 5\text{Solution 5}

Pt.1 求 sks_k 的值

仿照例2知，a,ba,b 是矩阵 BB 的两个元素，请读者自己尝试推导。

Pt.2 求 ∑k=1nck\sum_{k=1}^{n}{c^k}

简单等比数列求和。
设
S=∑k=1nck...(4)S=\sum_{k=1}^{n}{c^k}\quad\text{...(4)}
则
cS=∑k=1nck+1...(5)cS=\sum_{k=1}^{n}{c^{k+1}}\quad\text{...(5)}
(5)−(4)(5)-(4) 得
(c−1)S=ck+1−c(c-1)S=c^{k+1}-c
S=ck+1−cc−1S=\frac{c^{k+1}-c}{c-1}
使用乘法逆元和快速幂即可完成求解。

#include<cstdio>
#include<cstdlib>
#include<cstring>

#define reg register
typedef long long ll;

struct node{
	int x,y;
	ll a[5][5];
	node(){
		memset(a,0,sizeof(a));
		x=y=0;
	}
}s,t;

ll a,b,c;
ll n,m,x,y,son;
ll sum=0;

node T(node a,node b){
	node c;c.x=a.x;c.y=b.y;
	for(reg int i=1;i<=c.x;++i)
		for(reg int j=1;j<=c.y;++j)
			for(reg int k=1;k<=a.y;++k)
				c.a[i][j]=(c.a[i][j]+a.a[i][k]*b.a[k][j])%m;
	return c;
}
node QT(node a,int b){
	if(b<=1) return a;
	node c=QT(a,b/2);
	if(b%2)
		return T(T(c,c),a);
	else
		return T(c,c);
}
void exgcd(ll A,ll B,ll&x,ll&y){
	if(!B){
		x=1;y=0;
		return;
	}
	exgcd(B,A%B,x,y);
	ll t=x;x=y;y=t-A/B*y;
}
ll qt(ll a,int b){
	if(b==1) return a;
	ll e=qt(a,b/2)%m;
	if(b%2)
		return (e*e)%m*a%m;
	else
		return (e*e)%m;
}
int main(){
	scanf("%lld%lld%lld%lld%lld%lld%lld",&a,&b,&c,&n,&m,&s.a[1][1],&s.a[1][2]);c%=m;
	s.a[1][3]=s.a[1][1]+s.a[1][2];s.a[1][4]=s.a[1][1];
	s.x=1;s.y=4;t.x=t.y=4;
	t.a[2][1]=t.a[3][3]=t.a[3][4]=t.a[4][4]=1;
	t.a[1][2]=t.a[1][3]=b;t.a[2][2]=t.a[2][3]=a;
	s=T(s,QT(t,n-2));
	exgcd(c-1,m,x,y);son=(qt(c,n+1)-c+m)%m;x=(x%m+m)%m;
	printf("%lld",(((s.a[1][3]*n)%m-s.a[1][4]+m)%m+(x*son)%m+m)%m);
}
posted @ 2019-04-02 09:41  TeacherDai  阅读(289)  评论(0)    收藏  举报