数学专题（2）、快速幂

嗯，今天是快速幂专题。

相信大家应该都会求幂，即求a^n.那么这个东西有什么可讲的呢？当然有！朴素的幂计算复杂度是O(n)的，而我们的快速幂则能够达到O(logn)级别。可能你会觉得线性复杂度已经够快了，但是在实际应用中，你就会发现快速幂的高效之处。

那什么是快速幂呢？

快速幂就是在计算a的b次方的时候，把b写成二进制和的形式。

比如说计算a^43,我们就可以把43写成二进制数101011，每一位乘以对应的权值，就变成了a^(2^0+2^1+2^3+2^5).

那么这么做有什么好处呢？我们的循环从原来的43次降为log(43)次，每次不再是n-=1而是n>>=1.这样就能大大降低复杂度了。

快速幂一般还可以搭配矩阵使用，成为矩阵快速幂。

我们在做矩阵乘法运算的时候，首先有个前提，就是矩阵A的列与矩阵B的行必须相同才能进行乘法运算。最后的结果C的行数等于A的行数，C的列数等于B的列数。C中某个元素等于A对应行向量与B对应列向量乘积之和。这种算法的时间复杂度为O(n^3).

实际上，在clrs里，介绍了一种O(n^2.81)复杂度的矩阵乘法，即Strassen算法，它的原理是利用分治，将矩阵分块，再计算分块矩阵的乘积，最后合并，在这儿不是重点表过不提，感兴趣的可以去看一下。

我们今天要将的是如何求矩阵的幂。我们可以很清楚地看到，如果用朴素矩阵乘法，复杂度会很令人厌恶，所以我们运用快速幂思想，将这个问题转化成了快速幂问题。矩阵快速幂有很多应用，例如给dp加速之类的，这儿我们讲一个最简单的应用。

想必大家都知道斐波那契数列（Fibonacci sequence）吧！她是数学史上最著名的数列之一，有着许多有趣的性质。那我们怎么计算第n项斐波那契数列是多少呢？我们至少有四种办法可以计算。

最朴素的计算方法就是递归计算。由于斐波那契数列的递推式是F(n)=F(n-1)+F(n-2)(n>2)，所以我们可以递归计算F(n).但是这种方法的时间复杂度高达O(2^n),是无论如何不可能接受的。

我们可以使用非递归形式计算递推式，只要写一个一重循环就可以了。这种方法是线性复杂度，但仍不够好。

我们拥有一个常数复杂度的方法，就是计算它的通项公式。我们用特征方程与特征根求出通项公式直接计算，就可以把复杂度压缩到O(1),但是有个问题，这个整数数列的通项公式里含有大量的无理数，n一旦变大之后会对精度造成影响。

最后一种自然就是我们的矩阵快速幂。我们可以构造如下的矩阵递推式：

|Fn+1 Fn   |      | 1         1  |^n

                    = 

|Fn     Fn-1|      | 1          0 |    

这时候用到我们的矩阵快速幂就可以了。代码如下：

#include<cstdio>
#include<iostream>
#include<cmath>
#include<algorithm>
#include<cstring>
using namespace std;

typedef long long ll;
const int mod=1e4;

struct mat{
    ll a[2][2];
};

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 quickpow(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));
    res.a[0][0]=res.a[1][1]=1;
    while(n){
        if(n&1)
            res=mul(res,c);
        c=mul(c,c);
        n>>=1;
    }//快速幂
    cout<<res.a[0][1]<<endl;
}

int main(){
    int n;
    while(cin>>n){
        if(n==-1) break;
        quickpow(n);
    }
    return 0;
}