高效计算Fibonacci数

From: http://blog.csdn.net/jcwkyl/article/details/3677747

以前只知道使用递归或递推的方法，最近在CSDN论坛上学到一种新的解法，在自己所知的几个算法里，它具有最好的运行效率。这种方法使用下面这个关于Fibonacci数的矩阵恒等式： 

这个算法就是根据这个恒等式，通过计算等式右边的那个矩阵的n次方来计算第n个Fibonacci数。n次方的计算使用快速模幂算法，这样计算n次方只用做log(n)次2*2矩阵的乘法，而每次乘法的计算时间是常量。所以，上述方法计算第n个Fibonacci数的时间复杂度是O(log(n))。相比之下，用递归来计算第n个Fibonacci的时间复杂度是O(Fn)，是指数时间复杂度的。用递推来计算，是O(n)，然而当n很大时，n之于log(n)，等于无穷大之于无穷小。
顺便把关于Fibonacci的常用等式再总结一下：

From: http://tieba.baidu.com/p/1890441099

Fibonacci数是组合数学中非常重要的一个数列，它的递推公式是：
　　F(1)=F(2)=1
　　F(n)=F(n-1)+F(n-2)
　　当然，用这个公式来计算F(n)是非常慢的，当计算F(n)时需要从F(1)一直计算到F(n)。Fibonacci数列还满足一些其他的公式，如：
　　F(a+b+1)=F(a+1)*F(b+1)+F(a)*F(b)
　　利用这个公式，可以加速Fibonacci数的计算。我们考虑同时计算F(2n+1)和F(2n)，则按照上面的公式：
　　F(2n+1)=F(n+1)*F(n+1)+F(n)*F(n)
　　F(2n)=F(n+1)*F(n)+F(n)*F(n-1)=F(n+1)*F(n)+F(n)*(F(n+1)-F(n))
　　这样，F(2n+1)和F(2n)的计算变为了F(n+1)和F(n)的计算，即下标变为了原来的一半。重复利用这种方法，可以每次让下标变为原来的一半，总共需要大约log n次计算(以2为底)。
　　当n较大时，后面的方法就比直接的递推要快得多，比如当n=1000000时，后面的方法大概需要20次计算，而直接递推的方法大概需要1000000次计算

From: http://blog.sina.com.cn/s/blog_993d2542010143qw.html

现在的问题转换为求矩阵{1, 1, 1, 0}的乘方。如果简单第从0开始循环，n次方将需要n次运算，并不比前面的方法要快。但我们可以考虑乘方的如下性质：

        /  a^n/2*a^n/2                      n为偶数时
aⁿ=
        /  a(^n-1)/2*a^(n-1)/2            n为奇数时

要求得n次方，我们先求得n/2次方，再把n/2的结果平方一下。如果把求n次方的问题看成一个大问题，把求n/2看成一个较小的问题。这种把大问题分解成一个或多个小问题的思路我们称之为分治法。这样求n次方就只需要logn次运算了。

实现这种方式时，首先需要定义一个2×2的矩阵，并且定义好矩阵的乘法以及乘方运算。当这些运算定义好了之后，剩下的事情就变得非常简单。完整的实现代码如下所示。

#include <cassert>

///
// A 2 by 2 matrix
///
struct Matrix2By2
{
      Matrix2By2
      (
            long long m00 = 0, 
            long long m01 = 0, 
            long long m10 = 0, 
            long long m11 = 0
      )
      :m_00(m00), m_01(m01), m_10(m10), m_11(m11) 
      {
      }

      long long m_00;
      long long m_01;
      long long m_10;
      long long m_11;
};

///
// Multiply two matrices
// Input: matrix1 - the first matrix
//        matrix2 - the second matrix
//Output: the production of two matrices
///
Matrix2By2 MatrixMultiply
(
      const Matrix2By2& matrix1, 
      const Matrix2By2& matrix2
)
{
      return Matrix2By2(
            matrix1.m_00 * matrix2.m_00 + matrix1.m_01 * matrix2.m_10,
            matrix1.m_00 * matrix2.m_01 + matrix1.m_01 * matrix2.m_11,
            matrix1.m_10 * matrix2.m_00 + matrix1.m_11 * matrix2.m_10,
            matrix1.m_10 * matrix2.m_01 + matrix1.m_11 * matrix2.m_11);
}

///
// The nth power of matrix 
// 1  1
// 1  0
///
Matrix2By2 MatrixPower(unsigned int n)
{
      assert(n > 0);

      Matrix2By2 matrix;
      if(n == 1)
      {
            matrix = Matrix2By2(1, 1, 1, 0);
      }
      else if(n % 2 == 0)
      {
            matrix = MatrixPower(n / 2);
            matrix = MatrixMultiply(matrix, matrix);
      }
      else if(n % 2 == 1)
      {
            matrix = MatrixPower((n - 1) / 2);
            matrix = MatrixMultiply(matrix, matrix);
            matrix = MatrixMultiply(matrix, Matrix2By2(1, 1, 1, 0));
      }

      return matrix;
}

///
// Calculate the nth item of Fibonacci Series using devide and conquer
///
long long Fibonacci_Solution3(unsigned int n)
{
      int result[2] = {0, 1};
      if(n < 2)
            return result[n];

      Matrix2By2 PowerNMinus2 = MatrixPower(n - 1);
      return PowerNMinus2.m_00;
}

From: http://bbs.bccn.net/thread-352527-1-1.html

//2ms AC代码
#include <stdio.h>
struct Matrix
{
    int _11,_12;
    int _21,_22;
};
void Ride(const Matrix *src1,const Matrix *src2,Matrix *out)
{
    out->_11 = (src1->_11*src2->_11 + src1->_12*src2->_21)%10000;
    out->_12 = (src1->_11*src2->_12 + src1->_12*src2->_22)%10000;
    out->_21 = (src1->_21*src2->_11 + src1->_22*src2->_21)%10000;
    out->_22 = (src1->_21*src2->_12 + src1->_22*src2->_22)%10000;
}
void Power(Matrix *src,int n,Matrix *out)
{
    Matrix one = {1,1,1,0};
    if(n == 1)
    {
        *out = one;
        return ;
    }
    Power(src,n>>1,out);
    Matrix temp = *out;
    Ride(&temp,&temp,out);
    if(0 == n%2)
        return ;
    else
    {
        temp = *out;
        Ride(&temp,&one,out);
    }
}
int main()
{
    int i,j,n;
    while(scanf("%d",&n) && n >=0)
    {
        if(0 == n || 1 == n)
        {
            printf("%d\n",n);
            continue;
        }
        Matrix m1 = {1,1,1,0},res;
        Power(&m1,n-1,&res);
        printf("%d\n",res._11);
    }
    return 0;
}