hdu 3306 Another kind of Fibonacci 矩阵快速幂

最新推荐文章于 2025-09-21 17:57:55 发布

原创最新推荐文章于 2025-09-21 17:57:55 发布 · 1k 阅读

CC 4.0 BY-SA版权

本文详细介绍了使用快速幂算法解决特定Fibonacci数列平方和问题的方法，包括矩阵快速幂的原理及实现过程。通过定义新的Fibonacci数列规则，并利用矩阵运算简化计算步骤，实现高效求解。

Now we define another kind of Fibonacci : A(0) = 1 , A(1) = 1 , A(N) = X * A(N - 1) + Y * A(N - 2) (N >= 2).And we want to Calculate S(N) , S(N) = A(0)² +A(1)²+……+A(n)².

s( n ) = s(n-1) + An^2

An^2 = X^2*An-1 ^2 + Y^2*An-2^2 + 2*X*Y*An-1*An-2

An*An-1 = X*An-1^2 + Y*An-1 *An-2

矩阵【Sn-1 ， An-1*An-2 An-1^2 An-2 ^2】

矩阵2

{ 1 0 0 0 }

{2XY Y 2XY 0}

{X^2 X X^2 0 }

{Y^2 0 Y^2 0}

通过矩阵就可以递推。

这样就可以矩阵快速幂。

#include <vector>
#include <list>
#include <map>
#include <set>
#include <deque>
#include <stack>
#include <cstring>
#include <bitset>
#include <algorithm>
#include <functional>
#include <numeric>
#include <utility>
#include <sstream>
#include <iostream>
#include <iomanip>
#include <cstdio>
#include <cmath>
#include <cstdlib>
#include <ctime>
#include <assert.h>
#include <queue>
#define REP(i,n) for(int i=0;i<n;i++)
#define TR(i,x) for(typeof(x.begin()) i=x.begin();i!=x.end();i++)
#define ALLL(x) x.begin(),x.end()
#define SORT(x) sort(ALLL(x))
#define CLEAR(x) memset(x,0,sizeof(x))
#define FILLL(x,c) memset(x,c,sizeof(x))
using namespace std;
const double eps = 1e-9;
#define LL long long 
#define pb push_back
const int  maxn = 4;
const int size = 4;
const LL MOD =  10007 ;
LL X, Y;
LL n;
struct Mat{
	long long mat[maxn][maxn];
	void init(){
		 memset(mat,0,sizeof(mat));
		 mat[0][0] = 1;
		 mat[1][0] = 2*X%MOD*Y%MOD;  mat[1][1] = Y ; mat[1][2] = 2*X%MOD*Y%MOD;
		 mat[2][0] = X*X%MOD;        mat[2][1] = X%MOD;  mat[2][2] = X*X%MOD; mat[2][3] =1;
		 mat[3][0] = Y*Y%MOD;        mat[3][2] =  Y*Y%MOD;       
	}
	void show(){
		for(int i=0;i<=4;i++){
			for(int j=0;j<=4;j++){
				 cout << mat[i][j]<< " ";
			}
			cout <<endl;
		}
	}
	LL get_ans(){
		return (mat[1][0] + mat[2][0]+ mat[3][0])%MOD;
	}
}E,g;
void init_E(){
   memset(E.mat,0,sizeof(E.mat));
   for(int i=0;i<size ;i++){
      E.mat[i][i] =1 ;
   }
}
Mat operator*(const Mat &a,const Mat &b){
	Mat c;
	for(int i=0;i<size;i++)
		for(int j=0;j<size;j++){
			c.mat[i][j] = 0;
			for(int k=0;k<size;k++){
				if(a.mat[i][k] && b.mat[k][j])
					c.mat[i][j]=(c.mat[i][j]+a.mat[i][k]*b.mat[k][j]%MOD)%MOD;
			}
			
		}
	return c;
}

Mat operator^(Mat a,LL x){
	Mat c=E;
	for(;x;x>>=1){
		if(x&1)
			c=c*a;
		a=a*a;
	}
	return c;
}
Mat pow2(Mat a,LL x){
	Mat c = E;
	while(x){
		if(x&1){
			c = c*a;
		}
		a= a*a;
		x= x>>1;
	}
	return c;
}
void solve(){
	g.init();
	g = g^(n-1);
	LL ans = g.get_ans();
	ans +=2;
	ans %= MOD;
	cout << ans<<endl;
}
int main(){
	init_E();
    while(~scanf("%lld",&n)){
    	scanf("%lld%lld",&X,&Y);
    	
    	solve();
    }
    return 0;
}