在C中使用黄金比例计算第n个斐波那契数模m

Question

我正在尝试计算第 n 个斐波那契数模 10^9+7，其中 n 由用户输入。

我用黄金比例计算斐波那契数。

以下代码在 n=43 之前产生正确的结果。但是对于 n>=44，phi 超过 10^9+7，我开始得到意想不到的结果。此外，如果移除模数，n>=44 会给出正确的结果。

#include <stdio.h>
#include <math.h>
long double mod=1000000007;

long double power(long double base, long long int expo)
{
    if(base==1  ||  expo==0)
        return 1;
    if(expo&1)
    {
        long double temp = power(base, expo>>1);
        return fmodl(base * fmodl(temp*temp, mod), mod);
    }
    else
    {
        long double temp=power(base, expo>>1);
        return fmodl(temp*temp,mod);
    }
}


int main(void) {
    // your code goes here
    long double phi = (1+powl(5, 0.5))/2;
    long double phi_cap = (1 - powl(5, 0.5))/2;
    long double root5 = powl(5, 0.5);
    long long int n;
    scanf("%lld",&n);
    long double ans = fmodl(  (power(phi, n) - power(phi_cap, n)) * power(root5,mod-2), mod);
    printf("%.0Lf\n", ans);
    return 0;
}

为什么会这样？使用long double存储无理数是不是错了？

谢谢

Answer 1

closed form expression for the n'th Fibonacci number是

F_n = φⁿ / √5 - ψⁿ / √5

哪里

φ = 1/2 + √5/2 ≅ 1.6180339887
ψ = 1/2 - √5/2 = φ - 1 ≅ 0.6180339887

由于 F_n 中减去的部分总是小于二分之一，我们可以计算 F_n via rounding 趋近于零（或趋向负无穷大，或 floor()，因为 F_n 对于 n ≥ 0 是非负的）：

F_n = ⌊ φⁿ / √5 ⌋ = floor( φⁿ / √5 )

如果我们对F_nmodulo M感兴趣，M≥2，需要观察C中的modulo operation affects the above formula. Note that the expression "expr MOD M" is typically computed using fmod(expr, M)

我们可以应用 property of modular multiplication：对于正实数 a、b 和 m, (a mod m)(b mod m) = (a b mod m).楼层操作不受影响。这里，a = φⁿ，b = 1/√5。这意味着我们可以简化表达式

F_n MOD M = ⌊ φⁿ / √5 ⌋ MOD M

进入

F_n MOD M = ⌊ ( φⁿ MOD (M√5) ) / √5 ⌋

这里，我们可以套用modular exponentiation，注意此时的modulus不是整数M，而是M·√5 .

换句话说，如果我们有一个函数，它使用 mod 具有浮点数 modulus 的平方幂来计算浮点值的整数次幂，比如说

double pow_mod(double base, unsigned int exponent, double modulus);

我们可以使用

计算第n个斐波那契数F_n modulo m

#define PHI    1.61803398874989484820458683436563811772
#define SQRT5  2.236067977499789696409173668731276235441;

double fn_mod_m;
unsigned int n, m;

fn_mod_m = floor(pow_mod(PHI, n, SQRT5 * (double)m) / SQRT5);

right-to-left method 是 mod 平方求幂的极佳候选。

在C中使用黄金比例计算第n个斐波那契数模m

Calculating nth fibonacci number modulo m using golden ratio in C

c

fibonacci

modulo