BZOJ 3240: [Noi2013]矩阵游戏 矩阵乘法 费马小定理

3240: [Noi2013]矩阵游戏

Time Limit: 10 Sec  Memory Limit: 256 MB
Submit: 1752  Solved: 800
[Submit][Status][Discuss]

Description

婷婷是个喜欢矩阵的小朋友,有一天她想用电脑生成一个巨大的n行m列的矩阵(你不用担心她如何存储)。她生成的这个矩阵满足一个神奇的性质：若用F[i][j]来表示矩阵中第i行第j列的元素，则F[i][j]满足下面的递推式:
F[1][1]=1
F[i,j]=a*F[i][j-1]+b (j!=1)
F[i,1]=c*F[i-1][m]+d (i!=1)
递推式中a,b,c,d都是给定的常数。
现在婷婷想知道F[n][m]的值是多少,请你帮助她。由于最终结果可能很大，你只需要输出F[n][m]除以1,000,000,007的余数。

Input

一行有六个整数n,m,a,b,c,d。意义如题所述

Output

包含一个整数，表示F[n][m]除以1,000,000,007的余数

Sample Input

3 4 1 3 2 6

Sample Output

85

HINT

样例中的矩阵为：
1 4 7 10
26 29 32 35
76 79 82 85

1<=N,M<=10^1000 000,a<=a,b,c,d<=10^9

这个题。。显然是可以用数学方法求通项的

然鹅

我用了矩乘

之后天真地把n,m模了1e9+7

一顿 狂WA

最后终于想明白

之所以n,m可以模

是因为在矩乘过程中n,m是次数

模数是质数根据费马小定理才可以搞

特别的 a或c==1时就还是模1e9+7

具体自己YY一下就明白了

#include<cmath>
#include<ctime>
#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<complex>
#include<iostream>
#include<algorithm>
#include<iomanip>
#include<vector>
#include<string>
#include<bitset>
#include<queue>
#include<set>
#include<map>
using namespace std;

typedef long long ll;

inline ll read()
{
	ll x=0,f=1;char ch=getchar();
	while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();}
	while(ch<='9'&&ch>='0'){x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
	return x*f;
}
void print(ll x)
{if(x<0)putchar('-'),x=-x;if(x>=10)print(x/10);putchar(x%10+'0');}

const int N=1001000;const ll mod=int(1e9)+7;

struct matrix
{
	ll ma[3][3];
	friend matrix operator *(const matrix &x,const matrix &y)
	{
		register int i,j,k;matrix z;
		for(i=1;i<3;++i)for(j=1;j<3;++j)z.ma[i][j]=0;
		for(i=1;i<3;++i)for(j=1;j<3;++j)for(k=1;k<3;++k)(z.ma[i][j]+=x.ma[i][k]*y.ma[k][j]%mod)%=mod;
		return z;
	}
};

char sn[N],sm[N];

ll n,m,a,b,c,d;

matrix tmp,t;

matrix qpow(matrix tmp,ll x)
{
	matrix res;
	res.ma[1][1]=1;res.ma[1][2]=0;
	res.ma[2][1]=0;res.ma[2][2]=1;
	while(x)
	{
		if(x&1)res=res*tmp;
		tmp=tmp*tmp;
		x>>=1;
	}
	return res;
}

int main()
{
	scanf("%s%s",sn,sm);
	ll modn=mod,modm=mod;register int i;
	a=read();b=read();c=read();d=read();
	if(a^1)modn--;if(c^1)modm--;
	int len=strlen(sn);for(i=0;i<len;++i)n=(n*10+sn[i]-'0')%modn;
	len=strlen(sm);for(i=0;i<len;++i)m=(m*10+sm[i]-'0')%modm;
	tmp.ma[1][1]=a;tmp.ma[1][2]=0;
	tmp.ma[2][1]=b;tmp.ma[2][2]=1;
	t.ma[1][1]=c;t.ma[1][2]=0;
	t.ma[2][1]=d;t.ma[2][2]=1;
	tmp=qpow(tmp,m-1);t=tmp*t;
	t=qpow(t,n-1);t=t*tmp;
	print((t.ma[1][1]+t.ma[2][1])%mod);puts("");
	return 0;
}
/*
3 4 1 3 2 6

85
*/