HDOJ 3117 - Fibonacci Numbers

Matrix Multiplication (& Quick Power) 

Description

输入n，输出斐波那契数列f(n)。

关键在于，如果f(n)大于8位，那么就输出“前4位...后4位”这种格式。

Type

Matrix Multiplication

Quick Power

Analysis

这题虽然是矩阵乘法和快速幂。

但是其真正的难点在于f(n)大于8位后，如何精确得算出前4位。

小于8位，我们可以递推得到。

大于8位的后4位，我们可以用矩阵乘法和快速幂得到。

而大于8位的前4位，如果我们用矩阵乘法和快速幂，然后对计算中的数只取前几位的话，

当n接近10^8的时候，精度误差会越来越大，导致WA。

那么如何去求前4位呢……

我们可以用通项公式：

对于后面一项，

因为括号内的式子约等于-0.75，在n >= 40（结果超过8位）的时候，已经近乎于0。

在n趋近正无穷大的时候，更是趋近无穷小，所以可以忽略。

之后我们把这个式子两边对10取对数，为的是把指数n化为乘数n，免去了求高阶幂的过程。

现在可以计算出log10 f(n)了，但是我们仍然会得到一个很大的数字，直接作为指数求10的幂，结果依然不精确。

我们知道log10 f(n)的整数部分，相当于我们的结果最后大概有几个0，或者说是个几位数。

而我们所求的答案并不需要知道所求的斐波那契数有几位。

因此，我们可以把整数部分去掉，只保留小数部分。

换句话说，我们所求的前四位数是什么，跟log10 f(n)的整数部分一点**关系都木有。

这时候我们可以得到我们的结果了，把log10 f(n)的小数部分作为指数求10的幂，然后前四位即是我们要的答案。

可以求幂后乘以1000再取整，即得到前四位的数字。

这样做，在求幂的时候甚至都不用快速幂，直接利用pow( )函数搞定。

而且充分利用了通项公式在计算机计算中，n越大，越准确，数位越高，越准确的特点。

又避免了通项公式在n较小或者数位较低时不准确的缺点(因为这一部分我们是用矩阵乘法和快速幂做的)。

而这个牛逼方法是网上各种解题报告的答案，不知道源自何处，我就不提出处了，反正不是我想的。

Solution

// HDOJ 3117
// Fibonacci Numbers
// by A Code Rabbit

#include <cstdio>
#include <cstring>
#include <cmath>

const int MAXO = 2;
const int MOD = 10000;

enum Need {
    ALL,
    FIRST,
    LAST,
};

template <typename T>
struct Matrix {
    T e[MAXO][MAXO];
    int o;
    Matrix(int ord) { memset(e, 0, sizeof(e)); o = ord; }
    Matrix operator*(const Matrix& one) {
        Matrix res(o);
        for (int i = 0; i < o; i++)
            for (int j = 0; j < o; j++)
                for (int k = 0; k < o; k++)
                    res.e[i][j] += e[i][k] * one.e[k][j];
        return res;
    }
    Matrix operator%(int mod) {
        for (int i = 0; i < o; i++)
            for (int j = 0; j < o; j++)
                e[i][j] %= mod;
        return *this;
    }
    Matrix operator%=(int mod) { return *this = *this % mod; }
};

/* if bo is true, it means that the program is competing the first four numbers of Fibonacci numbers. */
template <typename T>
T QuickPower(T rdx, int exp, int mod, Need need) {
    T res = rdx;
    exp--;
    while (exp) {
        if (exp & 1) res = res * rdx;
        exp >>= 1;
        rdx = rdx * rdx;
        if (need == LAST) {
            res %= mod;
            rdx %= mod;
        }
    }
    return res;
}

int n;

int Fibonacci(int x, Need need);

int main() {
    while (scanf("%d", &n) != EOF) {
        if (n < 40) printf("%d\n", Fibonacci(n, ALL));
        else printf("%d...%04d\n", Fibonacci(n, FIRST), Fibonacci(n, LAST));
    }

    return 0;
}

int Fibonacci(int x, Need need) {
    if (need == FIRST) {
        double ans = x * log10((1.0 + sqrt(5.0)) / 2.0) - 0.5 * log10(5.0);
        ans = ans - (int)ans;
        ans = pow(10, ans);
        return (int)(ans * 1000);
    }
    Matrix<long long> mat_one(2);
    mat_one.e[0][0] = 1;
    mat_one.e[0][1] = 1;
    mat_one.e[1][0] = 1;
    if (x) {
        Matrix<long long> mat_ans = QuickPower(mat_one, x, MOD, need);
        return mat_ans.e[1][0];
    } else {
        return 0;
    }
}