矩阵-斐波那契数列

利用矩阵来求解斐波那契数列的有关问题是ACM题中一个比较常见的题型。例:NYOJ 148(斐波那契数列2)

有关斐波那契树列的规律详见这里

(1)、对于n>1,都有f(n)与f(n-1)互质。

(2)、f(n)=f(i)*f(n-i-1)+f(i+1)*f(n-i)。

现在说说怎么利用矩阵来求解斐波那契数列。

我们可以先保存b=f(1),a=f(0),然后每次设: b'=a+b a'=b。后利用a'和b'一直循环即可。同时我们可以将b a看做一个向量[b a],前面的操作就可以乘以矩阵:

|1 1|*[b a]=[a+b b]。

|1 0|

也就是说,如果我们要求第100个fibonacci数,只需要将矩阵[1 0]乘上 
1 1 
1 0 
的一百次方,再取出第二项即可。就像题目中的描述一样:设F0 = 0,F1 = 1,那么有:

.

这个矩阵的n次方的形式是:

F(n+1)  F(n)

F(n)   F(n-1),

现在的问题是如何求出这个矩阵的n次方是多少?

可以使用的方法有两种:

(1)、可以利用二分的思想:假设要求矩阵的N次方为A(N),设i=N/2若N%2==1,则 A(N)=A(i)*A(i)*A(1)若N%2==0,则A(N)=A(i)*A(i)。

(2)、利用二进制的思想:将N拆为二进制数,譬如13=(1101)2那么 A^13= A^8 * A^4 * A^1。这个在求快速幂模的时候经常用。具体见代码,可做模板,写的比较的乱:

#include<iostream>
#include<cstring>
using namespace std;
const int MAX=10;
#define __int64 long long
#define Bit(n) 1<<n
#define CLR(arr,val) memset(arr,val,sizeof(arr))
class Matrix{
public:
    Matrix(int r,int c):row(r),col(c){}
    void Init()  
    {   CLR(map,0);
        map[0][0]=map[0][1]=map[1][0]=1;
    }
    void Unit() //初始化为单位矩阵 
    {   CLR(map,0);
        for(int i=0;i<row;i++)
            map[i][i]=1;
    }
    int Result() const{return map[0][1]%10000;}
    friend Matrix operator*(const Matrix& ,const Matrix&);
    int Pow(int);
private:
    __int64 map[MAX][MAX];    
    int row,col;       
};
Matrix operator*(const Matrix& M1,const Matrix& M2) //矩阵相乘模板
{   Matrix M(M1.row,M2.col); //相乘之后矩阵的行和列会变化   
    for(int i=0;i<M1.row;i++)
        for(int j=0;j<M2.col;j++)
        {   M.map[i][j]=0;
            for(int k=0;k<M1.col;k++)
                M.map[i][j]+=M1.map[i][k]*M2.map[k][j]; 
            M.map[i][j]%=10000;
        }  
    return M;    
}
Matrix M(2,2);
int Matrix::Pow(int n) //矩阵快速幂 
{   Matrix temp(2,2); 
    temp.Init(); 
    for(int i=0;Bit(i)<=n;i++) //利用二进制的思想求解 
    {   if(Bit(i)&n) M=M*temp; 
        temp=temp*temp;     
    }
    return M.Result();    
}
int main()
{   __int64 num;
    while(cin>>num,num!=-1)
    {   M.Unit();
        cout<<M.Pow(num)<<endl;
    }
    return 0;
}

在这里顺便贴下斐波那契数列的通项公式:

 我们可以利用这个公式来简化很多运算:例HDU 2855,意思是输入n,m,求%m的结果。

具体推导过程如下:

 

 

 

 

 

 

 

 

 

 

 

 上个题目中求F(n)%m的值,这个是求F(2n)%m的值,你懂的~代码就不多写了~

调试了几次10^9这个数就是过不了,矩阵快速幂的模板改下就过了~具体改成下面的:

int Matrix::Pow(int n) //矩阵快速幂 
{   Matrix temp(2,2); 
    temp.Init(); 
    while(n) //利用二进制的思想求解 
    {   if(n&1) M=M*temp; 
        temp=temp*temp;     
        n>>=1;
    }
    return M.Result();    
}

例题2:NYOJ 507(斐波那契数列),只要写出那个递推关系式即可~其余的直接套模板:

S(n)=f(0)2+f(1)2+f(2)2+...+f(n)2

S(n-1)=f(0)2+f(1)2+f(2)2+...+f(n-1)2。 

所以:S(n)-S(n-1)=f(n)2=x2f(n-1)2+y2f(n-2)2+2xyf(n-1)f(n-2)。

所以可以构造下面的矩阵:

这个时候就简单了,具体见代码:

#include<iostream>
#include<cstring>
#include<cstdio>
using namespace std;
const int MAX=10;
#define __int64 long long 
#define CLR(arr,val) memset(arr,val,sizeof(arr))
int x,y,num;
class Matrix{
public:
    Matrix(int r,int c):row(r),col(c){}
    void Init()  
    {   CLR(map,0);
        x%=10007,y%=10007;//一定要把这个条件加上~~
        map[0][0]=map[0][1]=map[2][1]=1; 
        map[1][1]=x*x%10007;
        map[1][2]=y*y%10007;
        map[1][3]=(2*x*y)%10007;
        map[3][1]=x%10007,map[3][3]=y%10007; 
    }
    void Unit() //初始化为单位矩阵 
    {   CLR(map,0);
        for(int i=0;i<row;i++)
            map[i][i]=1;
    }
    void Clear()
    {   for(int i=0;i<4;i++)
            map[i][0]=1;
    }
    int Result() const {return map[0][0]%10007;}
    friend Matrix operator*(const Matrix& ,const Matrix&);
    friend ostream& operator<<(ostream& ,const Matrix&); 
    int Pow(int);
private:
    __int64 map[MAX][MAX];    
    int row,col;       
};
Matrix operator*(const Matrix& M1,const Matrix& M2) 
{   Matrix M(M1.row,M2.col);
    for(int i=0;i<M1.row;i++)
        for(int j=0;j<M2.col;j++)
        {   M.map[i][j]=0;
            for(int k=0;k<M1.col;k++)
                M.map[i][j]+=M1.map[i][k]*M2.map[k][j]; 
            M.map[i][j]%=10007;
        }  
    return M;    
}
Matrix M(4,4);
int Matrix::Pow(int n) 
{   Matrix temp(4,4),Start(4,1);
    temp.Init();
    Start.Clear();
    //cout<<temp; 
    while(n)
    {   if(n&1) M=M*temp;
        temp=temp*temp;
        n>>=1;
    }
    Start=M*Start;
    return Start.Result();    
}
ostream& operator<<(ostream &cout,const Matrix& M)
{   for(int i=0;i<M.row;i++)
    {   for(int j=0;j<M.col;j++)
            cout<<M.map[i][j]<<" ";
        cout<<endl;    
    }
    return cout;
}
int main()
{   while(scanf("%d%d%d",&num,&x,&y)!=EOF)
    {   M.Unit();
        cout<<M.Pow(num)<<endl;
    }
    return 0;
}

 

 


 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

已标记关键词 清除标记
©️2020 CSDN 皮肤主题: 大白 设计师:CSDN官方博客 返回首页