矩陣快速冪求斐波那契數列（總結）

 

第一部分：矩陣的基礎知識

1.結合性 （AB)C=A(BC）．

2.對加法的分配性 （A+B)C=AC+BC，C(A+B）=CA+CB ．

3.對數乘的結合性 k(AB）=（kA)B =A(kB）．

4.關於轉置 (AB)'=B'A'．

一個矩陣就是一個二維數組，為了方便聲明多個矩陣，我們一般會將矩陣封裝一個類或定義一個矩陣的結構體，我采用的是后者。（弱雞的我也直只會用結構體實現）

 

第二部分：矩陣相乘

若A為n×k矩陣，B為k×m矩陣，則它們的乘積AB(有時記做A·B)將是一個n×m矩陣。前一個矩陣的列數應該等於后一個矩陣的行數，得出的矩陣行數等於前一個矩陣的行數，列數等於后一個矩陣的行數。

其乘積矩陣AB的第i行第j列的元素為：

 

舉例：A、B均為3*3的矩陣：C=A*B，下面的代碼會涉及到兩種運算順序，第一種就是直接一步到位求，第二種就是每次求一列，比如第一次，C00+=a00*b00,C01+=a00*b01……第二次C00+=a00*b10，C01+=a01*b11……以此類推。。。

C00 = a00*b00 + a01*b10 + a02*b20
C01 = a00*b01 + a01*b11 + a02*b21 
C02 = a00*b02 + a01*b12 + a02*b22
C10 = a10*b00 + a11*b10 + a12*b20
C11 = a10*b00 + a11*b11 + a12*b21
C12 = a10*b02 + a11*b12 + a12*b22
C20 = a20*b00 + a21*b10 + a22*b20
C21 = a20*b01 + a21*b11 + a22*b21
C22 = a20*b02 + a21*b12 + a22*b22
C00 = a00*b00 + a01*b10 + a02*b20
C01 = a00*b01 + a01*b11 + a02*b21 
C02 = a00*b02 + a01*b12 + a02*b22

 

下面先來實現一個矩陣相乘的函數吧。

const int MOD=10000;
struct mat
{
    int a[2][2];//這里數據范圍就用小的示范
};
mat mat_mul(mat x,mat y)//實現兩個矩陣相乘，返回的還是一個矩陣。
{
    mat res;//用來表示得到的新的矩陣；
    memset(res.a,0,sizeof(res.a));
    for(int i=0;i<2;i++)
        for(int j=0;j<2;j++)
        for(int k=0;k<2;k++)
    {
        res.a[i][j]+=x.a[i][k]*y.a[k][j];
        res.a[i][j]%=MOD;//這一步看題目具體需要了
    }
    return res;
}

 

 

學了現代的話這個是很好理解的（個人認為）。

 

 

第三部分：矩陣快速冪   //其實和普通快速冪類似，只不過這里需要得到的是一個矩陣

神馬是冪？【很多時候會被高大上的名字嚇到。。。導致學習效率降低。。。其實沒辣么可怕，很簡單！！！】

冪又稱乘方。表示一個數字乘若干次的形式,如n個a相乘的冪為a^n ，或稱a^n為a的n次冪。a稱為冪的底數，n稱為冪的指數。——引自.度娘百科

這類題，指數都是很大很大很大很大很大很大很大的。。。霸王硬上弓的話，很容易超時的 T_T 。。。所以得快速冪→_→

學過之后發現，其實矩陣快速冪 的核心思想跟 以前學過的快速冪取模非常非常相似，只是矩陣乘法需要另外寫個函數，就是上面那個代碼。。。

 

快速冪的思路就是：

設A為矩陣，求A的N次方，N很大，1000000左右吧。。。

先看小一點的，A的9次方

A^9

= A*A*A*A*A*A*A*A*A  【一個一個乘，要乘9次】

= A*(A*A)*(A*A)*(A*A)*(A*A)【保持格式的上下統一，所以加上這句】

 = A*(A^2)^4 【A平方后，再四次方，還要乘上剩下的一個A，要乘6次】

= A*((A^2)^2)^2【A平方后，再平方，再平方，還要乘上剩下的一個A，要乘4次】

 

也算是一種二分思想的應用吧，1000000次冪，暴力要乘1000000次，快速冪就只要（log2底1000000的對數） 次，大約20次。。。這。。。我沒錯吧。。。

 

 單位矩陣： n*n的矩陣 mat ( i , i )=1; 任何一個矩陣乘以單位矩陣就是它本身 n*單位矩陣=n， 可以把單位矩陣等價為整數1。（單位矩陣用在矩陣快速冪中）

 

例如下圖就是一個7*7的單位矩陣：

 

下面來實現一個矩陣快速冪：

int pow(int n)//還是小范圍數據來說吧,要不然返回值的類型自己定義
{
    mat c,res;
    memset(res.a,0,sizeof(res.a));
    c.a[0][0]=1;//給矩陣賦初值
    c.a[0][1]=1;
    c.a[1][0]=1;
    c.a[1][1]=0;
    for(int i=0;i<n;i++) res.a[i][i]=1;//單位矩陣；
    while(n)
    {
        if(n&1) res=mat_mul(res,c);//這里看就要用到上面的矩陣相乘了；
        c=mat_mul(c,c);
        n=n>>1;
    }
    return res.a[0][1];
}//時間復雜度log(n)

 

         

      但是矩陣如何與斐波那契聯系在一起呢？？？

找了很多博客，看了第二位大神的博客才理解。

 

對於矩陣乘法與遞推式之間的關系：

 

如：在斐波那契數列之中

f[i] = 1*f[i-1]+1*f[i-2]  f[i-1] = 1*f[i-1] + 0*f[i-2];

即

所以

就這兩幅圖完美詮釋了斐波那契數列如何用矩陣來實現。

 

  下面一POJ3070/NYOJ148為例

  給出了矩陣相乘的定義，要你求出斐波那契的第n項對1e4取余。

  代碼一：在網上看到的很簡潔的代碼

#include <iostream>
#include <cstddef>
#include <cstring>
#include <vector>
using namespace std;
typedef long long ll;
const int mod=10000;
typedef vector<ll> vec;
typedef vector<vec> mat;
mat mul(mat &a,mat &b)//表示不會這樣用，，，，
{
    mat c(a.size(),vec(b[0].size()));
    for(int i=0; i<2; i++)
    {
        for(int j=0; j<2; j++)
        {
            for(int k=0; k<2; k++)
            {
                c[i][j]+=a[i][k]*b[k][j];
                c[i][j]%=mod;
            }
        }
    }
    return c;
}
mat pow(mat a,ll n)
{
    mat res(a.size(),vec(a.size()));
    for(int i=0; i<a.size(); i++)
        res[i][i]=1;//單位矩陣；
    while(n)
    {
        if(n&1) res=mul(res,a);
        a=mul(a,a);
        n/=2;
    }
    return res;
}
ll solve(ll n)
{
    mat a(2,vec(2));
    a[0][0]=1;
    a[0][1]=1;
    a[1][0]=1;
    a[1][1]=0;
    a=pow(a,n);
    return a[0][1];//也可以是a[1][0];
}
int main()
{
    ll n;
    while(cin>>n&&n!=-1)
    {
        cout<<solve(n)<<endl;
    }
    return 0;
}

代碼二：自己寫的代碼

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int MOD=10000;
struct mat
{
    ll a[2][2];
};
mat mat_mul(mat x,mat y)
{
    mat res;
    memset(res.a,0,sizeof(res.a));
    for(int i=0;i<2;i++)
        for(int j=0;j<2;j++)
        for(int k=0;k<2;k++)
        res.a[i][j]=(res.a[i][j]+x.a[i][k]*y.a[k][j])%MOD;
    return res;
}
void mat_pow(int n)
{
    mat c,res;
    c.a[0][0]=c.a[0][1]=c.a[1][0]=1;
    c.a[1][1]=0;
    memset(res.a,0,sizeof(res.a));
    for(int i=0;i<2;i++) res.a[i][i]=1;
    while(n)
    {
        if(n&1) res=mat_mul(res,c);
        c=mat_mul(c,c);
        n=n>>1;
    }
    printf("%I64d\n",res.a[0][1]);
}
int main()
{
    int n;
    while(~scanf("%d",&n)&&n!=-1)
    {
        mat_pow(n);
    }
    return 0;
}

 

 

  感覺不是很難，裸模板的話沒什么問題，但如果靈活度高一點我就不會構造矩陣了。還是繼續加油吧！

 

參考文章：

                 http://blog.csdn.net/u013795055/article/details/38599321

                 http://blog.csdn.net/g_congratulation/article/details/52734306