HDU3306 Another kind of Fibonacci 矩阵乘法

//题目参数 n,x,y
//f(n)=x*f(n-1)+y*f(n-2)
//s(n)=sum(f(i)^2)   0<=i<=n
//分解得到
//f(n)^2=(x^2*f(n-1)^2+y^2*f(n-2)^2+2xy*f(n-1)*f(n-2));
//s(n)=s(n-1)+f(n)^2;
//
//因为要求sn 我们构造矩阵[s(n-1),f(n)^2,f(n-1)^2,f(n)*f(n-1)]
//观察如何能够得到矩阵[s(n),f(n+1)^2,f(n)^2,f(n+1)*f(n)]
//很显然利用矩阵乘法，构造矩阵
//[ [ 1   0    0    0 ]
//  [ 1  x^2   1    x ]
//  [ 0  y^2   0    0 ]
//  [ 0  2xy   0    y ] ];
//因为 f(n+1)*f(n)=x*f(n)^2+y*f(n)*f(n-1);
//     f(n+1)^2=x^2*f(n)^2+y^2*f(n-1)^2+2xy*f(n)*f(n-1);

推公式的方法详见该博客，讲的十分清晰 ：http://blog.csdn.net/abcjennifer/article/details/5302198

Another kind of Fibonacci

Time Limit: 3000/1000 MS (Java/Others)    Memory Limit: 65536/65536 K (Java/Others)
Total Submission(s): 1155    Accepted Submission(s): 437

Problem Description

As we all known , the Fibonacci series : F(0) = 1, F(1) = 1, F(N) = F(N - 1) + F(N - 2) (N >= 2).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) ².

 

Input

There are several test cases.
Each test case will contain three integers , N, X , Y .
N : 2<= N <= 2 ³¹ – 1
X : 2<= X <= 2 ³¹– 1
Y : 2<= Y <= 2 ³¹ – 1

 

Output

For each test case , output the answer of S(n).If the answer is too big , divide it by 10007 and give me the reminder.

 

Sample Input

  
  

   
   2 1 1 
3 2 3 
  
  

 

Sample Output

  
  

   
   6
196
  
  

 

 

#include<iostream>
#include<cstring>
#include<string>
#include<cstdio>

using namespace std;

#define MAXN 5

struct Matrix
{
    int size;
    long modulo;
    long element[MAXN][MAXN];
    void setSize(int);
    void setModulo(long);
};

int n,x,y;
Matrix a,an;

void Matrix::setSize(int a)
{
    for (int i=0; i<a; i++)
        for (int j=0; j<a; j++)
            element[i][j]=0;
    size=a;
}

void Matrix::setModulo(long a)
{
    modulo = a;
}

Matrix multi(Matrix p1,Matrix p2)
{
    Matrix p;
    p.setSize(p1.size);
    p.setModulo(p1.modulo);
    for (int i=0; i<p1.size; i++)
        for (int j=0; j<p1.size; j++)
            for (int k=0; k<p1.size; k++)
            {
                p.element[i][j]+=p1.element[i][k]*p2.element[k][j];
                p.element[i][j]%=p1.modulo;
            }
    return p;
}

Matrix power(Matrix p,long exp)
{
    Matrix res,A;
    A=p;
    res.setModulo(p.modulo);
    res.setSize(p.size);
    for(int i=0;i<p.size;i++)
        res.element[i][i]=1;
    while(exp)
    {
        if(exp&1)
            res=multi(res,A);
        exp>>=1;
        A=multi(A,A);
    }
    return res;
}

void init()
{
    a.setSize(4);
    a.setModulo(10007);
    a.element[0][0]=a.element[1][0]=a.element[1][2]=1;
    a.element[1][1]=x*x%10007;
    a.element[2][1]=y*y%10007;
    a.element[3][1]=2*x*y%10007;
    a.element[1][3]=x;
    a.element[3][3]=y;
}

int main()
{
    int s0[]={1,1,1,1};
    while(~scanf("%ld%ld%ld",&n,&x,&y))
    {
        x%=10007;
        y%=10007;
        init();
        an=power(a,n);
        long tmp=0;
        for(int j=0;j<4;j++)
            tmp+=(an.element[j][0]*s0[j]);
        tmp=tmp%10007;
        printf("%ld\n",tmp);
    }
    return 0;
}