題目鏈接POJ307

題意：求第n位斐波那契數(shù)mod 10000的大小。其中n的大小高達(dá)1000000000

由于我對(duì)矩陣快速冪也不是很了解颅停，所以只能淺談谓晌，不敢深談，希望對(duì)新手有點(diǎn)幫助癞揉。

先聊聊快速冪

首先纸肉，如果要你自己寫一個(gè)pow()函數(shù)，需要求2²⁰喊熟，怎么求

#include <iostream>
using namespace std;
int quickPow(int a,int b){   //快速冪
    int ans = 1;
    while(b){
        if(b&1) ans*=a;
        a *= a;
        b /= 2;
    }
    return ans;
}
int pow(int a,int b){   //普通乘
    int ans = 1;
    for(int i=1; i<=b; i++)
        ans *= a;
    return ans;
}
int main(int argc, const char * argv[]) {
    cout<<quickPow(2, 20)<<endl;   //第一種方式
    cout<<pow(2,20)<<endl;    //第二種方式
    return 0;
}

上面的兩種方式柏肪，一種復(fù)雜度是O(logn)，一種復(fù)雜度是O(n)芥牌。當(dāng)你數(shù)量級(jí)很大的時(shí)候预吆，比如10000000000時(shí)，這個(gè)復(fù)雜度的時(shí)間差別就很大了

怎么得到這個(gè)快速冪的代碼的呢胳泉？我們假設(shè)需要求2⁷⁴

74轉(zhuǎn)化為二進(jìn)制為1001010，根據(jù)二進(jìn)制對(duì)應(yīng)位的1岩遗，可以把74拆分為64+8+2 扇商，所以
74 = 2⁶ + 2³ + 2¹
2⁷⁴ = 2²⁶ + 2²³ + 2²¹

 while(b){    //快速冪核心代碼
      if(b&1) ans*=a;
      a *= a;
      b /= 2;
 }

下面看一下上面這段快速冪程序中的ans和a的用法 (74二進(jìn)制1001010為從右到左下面表格從1-7)

這里有兩個(gè)變量，第一個(gè)變量a默認(rèn)用來存2ⁿ宿礁，案铺，變量ans存儲(chǔ)的是最后的結(jié)果，當(dāng)遇到第k位的二進(jìn)制位是1的時(shí)候梆靖，則ans ×= 2^k

舉個(gè)例子控汉，74的二進(jìn)制表示是1001010，因?yàn)槭?code>除2運(yùn)算返吻，所以從右往左遍歷姑子，到右邊第二個(gè)時(shí)，ans ×= 2²测僵，到右邊第四個(gè)時(shí)街佑，ans ×= 2⁸，到右邊第7個(gè)時(shí)捍靠，ans ×= 2⁶⁴沐旨，最后ans = 2⁷⁴，至于ans乘的過程中的那些2²榨婆、2⁸磁携、2⁶⁴怎么來？a變量不是一直在不斷增加嗎良风，ans乘的就是a變量對(duì)應(yīng)位的值谊迄，對(duì)應(yīng)上面的表格可能看的更清楚一點(diǎn)

不知道闷供，你看懂沒……如果還是沒懂，就給個(gè)面子裝懂好嘛……

矩陣鳞上？怎么過渡到斐波那契

• 宏觀上看矩陣

上圖中这吻，乘號(hào)右邊是一系列的矩陣。我們稱為B₁篙议、B₂唾糯、B₃、B₄ …… B_n
該序列中的同行為一個(gè)斐波那契數(shù)列鬼贱。乘號(hào)左邊是一個(gè)轉(zhuǎn)移矩陣A移怯，B_n = A×B_n-1

• 微觀上看矩陣

  
  

原始的序列是 [A,B]，需要得到的序列是[A+B,A]这难，根據(jù)矩陣相乘定律舟误，可以得到轉(zhuǎn)移矩陣A

• 總結(jié)一波
  現(xiàn)在需要求的是B_n%10000，B_n = B₁ × A^n-1姻乓，B₁是[1,1]嵌溢，所以 A^n-1就是需要的答案。那把求B_n改成了求A^n-1有什么好處呢蹋岩？見下文

矩陣快速冪的實(shí)際運(yùn)用

上述矩陣代替斐波那契數(shù)列的方式赖草。已經(jīng)把斐波那契數(shù)列從原先的遞加變成了現(xiàn)在的遞乘。那變成遞乘有什么用?

如果斐波那契第一項(xiàng)是A剪个，第二項(xiàng)是A²秧骑，第三項(xiàng)是A³，第4項(xiàng)是A⁴ ……
那第8項(xiàng)的值是扣囊，A⁸ = A⁴×A⁴

換句話說乎折，第4項(xiàng)的值 × 第4項(xiàng)的值 = 第8項(xiàng)的值

求第16項(xiàng)的值
原來的方式是：1->2->3->4……->15->16需要經(jīng)過15次運(yùn)算
現(xiàn)在的方式是：1->2->4->8->16只需要經(jīng)過4次運(yùn)算

為什么只需要經(jīng)過4次？從A->A²得到第2項(xiàng)侵歇，我要求第4項(xiàng)不需要經(jīng)過第3項(xiàng)骂澄，直接A² × A² -> A⁴就可以得到第4項(xiàng)。這又回到了剛開始的快速冪

具體怎么得到惕虑，下面請(qǐng)看代碼解析

代碼拆分

代碼不長(zhǎng)酗洒，先看下面這個(gè)函數(shù)。這里使用了引用傳遞枷遂，即a數(shù)組發(fā)生變化后樱衷，傳入的實(shí)參也會(huì)改變。例如傳入數(shù)組a和數(shù)組b酒唉，得到的數(shù)組a = 數(shù)組a × 數(shù)組b矩桂。即Matrix(cot,temp)執(zhí)行的結(jié)果是cot *= temp，且cot和temp都表示的是矩陣數(shù)組

void Matrix(int (&a)[2][2],int b[2][2]){   //矩陣乘法
    int tmp[2][2] = {0};
    for(int i = 0; i < 2; ++i)
        for(int j = 0; j < 2; ++j)
            for(int k = 0; k < 2; ++k)
                tmp[i][j] = (tmp[i][j] + a[i][k] * b[k][j]) % N;  //矩陣相乘，注意mod
    for(int i = 0; i < 2; ++i)
        for(int j = 0; j < 2; ++j)  
            a[i][j] = tmp[i][j]; //賦值返回?cái)?shù)據(jù)
}

A¹侄榴，也就是轉(zhuǎn)移矩陣雹锣，題中已經(jīng)給出了，用數(shù)組表示是[1,1,1,0]癞蚕。這里還需要一個(gè)A⁰蕊爵，A⁰其實(shí)就是單位矩陣E（E × 任意矩陣A = A）。單位矩陣用數(shù)組表示就是[1,0,0,1]（對(duì)角線都是1的矩陣）

while(n){
       if(n & 1) Matrix(cot,temp);   //如果是奇數(shù)
       Matrix(temp,temp);   
       n /= 2;  //不斷除2
}

上述5行代碼是矩陣快速冪關(guān)鍵（至于怎么得到這5行代碼的桦山，在上述快速冪的解釋中應(yīng)該已經(jīng)可以得到了）攒射。cot數(shù)組一開始是A⁰，temp數(shù)組一開始是A¹恒水。我們來模擬下求第13項(xiàng)的斐波那契值是怎么求的会放。

第一步:n = 13，執(zhí)行Matrix(cot,temp)和Matrix(temp,temp)钉凌，得到cot數(shù)組為A¹咧最，t數(shù)組為A²，n = 6

第二步:n = 6御雕，執(zhí)行Matrix(temp,temp)矢沿，得到cot數(shù)組為A¹，t數(shù)組為A⁴酸纲，n = 3

第三步:n = 3捣鲸，執(zhí)行Matrix(cot,temp)和Matrix(temp,temp)，得到cot數(shù)組為A⁵福青，t數(shù)組為A⁸，n = 1

第四步:n = 1脓诡，執(zhí)行Matrix(cot,temp)和Matrix(temp,temp)无午，得到cot數(shù)組為A¹³，t數(shù)組為A¹⁶祝谚，n = 0

退出while循環(huán)并輸出cot數(shù)組對(duì)應(yīng)的值宪迟。一共進(jìn)行了4步。

完整代碼

#include <iostream>
using namespace std;
int N = 10000,n;
void Matrix(int (&a)[2][2],int b[2][2]){
    int tmp[2][2] = {0};
    for(int i = 0; i < 2; ++i)
        for(int j = 0; j < 2; ++j)
            for(int k = 0; k < 2; ++k)
                tmp[i][j] = (tmp[i][j] + a[i][k] * b[k][j]) % N;
    for(int i = 0; i < 2; ++i)
        for(int j = 0; j < 2; ++j)
            a[i][j] = tmp[i][j];
}
int main(){
    while(cin>>n && n!=-1){
        int temp[2][2] = {1, 1, 1, 0},cot[2][2] = {1, 0, 0, 1};
        while(n){
            if(n & 1) Matrix(cot,temp);
            Matrix(temp,temp);
            n /= 2;
        }
        cout<<cot[0][1]<<endl;
    }
    return 0;
}

這篇文章主要利用斐波那契數(shù)列來淺談矩陣快速冪的用法交惯，矩陣快速冪還有更大的使用空間次泽，后續(xù)有空會(huì)補(bǔ)上矩陣快速冪的其他用法，另外送上一篇文章 十個(gè)利用矩陣乘法解決的經(jīng)典題目

后續(xù)學(xué)習(xí)的傳送門

矩陣快速冪——進(jìn)階
矩陣快速冪——實(shí)戰(zhàn)