轉載於：演算法學習筆記(4)：快速冪 - 知乎 (zhihu.com)

快速冪（Exponentiation by squaring，平方求冪）是一種簡單而有效的小演算法，它可以以的時間複雜度計算乘方。快速冪不僅本身非常常見，而且後續很多演算法也都會用到快速冪。

引入-例題

讓我們先來思考一個問題：7的10次方，怎樣算比較快？

方法1：最樸素的想法，7*7=49，49*7=343，... 一步一步算，共進行了9次乘法。

這樣算無疑太慢了，尤其對計算機的CPU而言，每次運算只乘上一個個位數，無疑太屈才了。這時我們想到，也許可以拆分問題。

方法2：先算7的5次方，即7*7*7*7*7，再算它的平方，共進行了5次乘法。

但這並不是最優解，因為對於“7的5次方”，我們仍然可以拆分問題。

方法3：先算7*7得49，則7的5次方為49*49*7，再算它的平方，共進行了4次乘法。

模仿這樣的過程，我們得到一個在時間內計算出冪的演算法，也就是快速冪。

遞迴快速冪

剛剛我們用到的，無非是一個二分的思路。我們很自然地可以得到一個遞迴方程：

計算a的n次方，如果n是偶數（不為0），那麼就先計算a的n/2次方，然後平方；如果n是奇數，那麼就先計算a的n-1次方，再乘上a；遞迴出口是a的0次方為1。

遞迴快速冪的思路非常自然，程式碼也很簡單（直接把遞迴方程翻譯成程式碼即可）：

//遞迴快速冪
int qpow(int a, int
 
 n)
{
    if (n == 0)
        return 1;
    else if (n % 2 == 1)
        return qpow(a, n - 1) * a;
    else
    {
        int temp = qpow(a, n / 2);
        return temp * temp;
    }
}

注意，這個temp變數是必要的，因為如果不把記錄下來，直接寫成qpow(a, n /2)*qpow(a, n /2)，那會計算兩次，整個演算法就退化為了。

在實際問題中，題目常常會要求對一個大素數取模，這是因為計算結果可能會非常巨大，但是在這裡考察高精度又沒有必要。這時我們的快速冪也應當進行取模，此時應當注意，原則是步步取模，如果MOD較大，還應當開long long。

//遞迴快速冪（對大素數取模）
#define MOD 1000000007
typedef long long ll;
ll qpow(ll a, ll n)
{
    if (n == 0)
        return 1;
    else if (n % 2 == 1)
        return qpow(a, n - 1) * a % MOD;
    else
    {
        ll temp = qpow(a, n / 2) % MOD;
        return temp * temp % MOD;
    }
}

大家知道，遞迴雖然簡潔，但會產生額外的空間開銷。我們可以把遞迴改寫為迴圈，來避免對棧空間的大量佔用，也就是非遞迴快速冪。

非遞迴快速冪

我們換一個角度來引入非遞迴的快速冪。還是7的10次方，但這次，我們把10寫成二進位制的形式，也就是。

現在我們要計算，可以怎麼做？我們很自然地想到可以把它拆分為. 實際上，對於任意的整數，我們都可以把它拆成若干個的形式相乘。而這些，恰好就是、、……我們只需不斷把底數平方就可以算出它們。

我們先看程式碼，再來仔細推敲這個過程：

//非遞迴快速冪
int qpow(int a, int n){
    int ans = 1;
    while(n){
        if(n&1)        //如果n的當前末位為1
            ans *= a;  //ans乘上當前的a
        a *= a;        //a自乘
        n >>= 1;       //n往右移一位
    }
    return ans;
}

最初ans為1，然後我們一位一位算：

1010的最後一位是0，所以a^1這一位不要。然後1010變為101，a變為a^2。

101的最後一位是1，所以a^2這一位是需要的，乘入ans。101變為10，a再自乘。

10的最後一位是0，跳過，右移，自乘。

然後1的最後一位是1，ans再乘上a^8。迴圈結束，返回結果。

這裡的位運算子，>>是右移，表示把二進位制數往右移一位，相當於/2；&是按位與，&1可以理解為取出二進位制數的最後一位，相當於%2==1。這麼一等價，是不是看出了遞迴和非遞迴的快速冪的關係了？雖然非遞迴快速冪因為牽扯到二進位制理解起來稍微複雜一點，但基本思路其實和遞迴快速冪沒有太大的出入。

快速冪的拓展

上面所述的都是整數的快速冪，但其實，在算時，只要a的資料型別支援乘法且滿足結合律，快速冪的演算法都是有效的。矩陣、高精度整數，都可以照搬這個思路。下面給出一個模板：

//泛型的非遞迴快速冪
template <typename T>
T qpow(T a, ll n)
{
    T ans = 1; // 賦值為乘法單位元，可能要根據建構函式修改
    while (n)
    {
        if (n & 1)
            ans = ans * a; // 這裡就最好別用自乘了，不然過載完*還要過載*=，有點麻煩。
        n >>= 1;
        a = a * a;
    }
    return ans;
}

注意，較複雜型別的快速冪的時間複雜度不再是簡單的，它與底數的乘法的時間複雜度有關。

例如，矩陣快速冪的一個經典應用是求斐波那契數列：

（洛谷P1962） 斐波那契數列

題目背景
大家都知道，斐波那契數列是滿足如下性質的一個數列：

題目描述
請你求出的值。

（以下內容涉及到基本的線性代數知識）

設矩陣，我們有，於是 :

這樣，我們把原來較為複雜的問題轉化成了求某個矩陣的冪的問題，這就可以應用快速冪求解了。

#include <cstdio>
#define MOD 1000000007
typedef long long ll;

struct matrix
{
    ll a1, a2, b1, b2;
    matrix(ll a1, ll a2, ll b1, ll b2) : a1(a1), a2(a2), b1(b1), b2(b2) {}
    matrix operator*(const matrix &y)
    {
        matrix ans((a1 * y.a1 + a2 * y.b1) % MOD,
                   (a1 * y.a2 + a2 * y.b2) % MOD,
                   (b1 * y.a1 + b2 * y.b1) % MOD,
                   (b1 * y.a2 + b2 * y.b2) % MOD);
        return ans;
    }
};

matrix qpow(matrix a, ll n)
{
    matrix ans(1, 0, 0, 1); //單位矩陣
    while (n)
    {
        if (n & 1)
            ans = ans * a;
        a = a * a;
        n >>= 1;
    }
    return ans;
}

int main()
{
    ll x;
    matrix M(0, 1, 1, 1);
    scanf("%lld", &x);
    matrix ans = qpow(M, x - 1);
    printf("%lld\n", (ans.a1 + ans.a2) % MOD);
    return 0;
}