An extremelly efficient method for generating the \(n\)th value of the Fibonacci sequence is given by the formula,
\[ \boxed{f_n = \frac{\Phi^n - (-\phi)^n}{\sqrt{5}}}\,, \]where \(\Phi\) is the golden ratio and \(-\phi\) is the golden ratio conjugate.
Every \(f_k\) where \(k \equiv 0 \pmod{3}\) is even and every \(f_k\) where \(k \equiv \{ 1, 2 \} \pmod{3}\) is odd.
The tricky thing is that when we write an implementation of this function, we need to be mindful that the golden ratio is a transcendental number meaning that it has infinitely many digits like \(\pi\). We leverage Julia's BigFloat
and BigInt
types to ensure we use enough digits to reduce rounding errors. Regardless, we'll still get plenty of rounding errors and as \(n\) approaches infinity, we will increasely get more results which need to be rounded down slightly. To combat this, after performing the computation, we round the result and then convert from BigFloat
to BigInt
.
Φ, ϕ = big(MathConstants.golden), big(MathConstants.golden) - 1.0
fib(n) = (Φ^n - (-ϕ)^n)/sqrt(5) |> round |> BigInt
@btime fib(100_000)
> 3.682 μs (22 allocations: 17.81 KiB)
@btime fib(1_000_000)
> 10.401 μs (22 allocations: 170.36 KiB)
A recursive implementation to generate the \(n\)th value of the Fibonacci sequence generally starts to be ridiculously unwieldy around \(n = 50\), so this is a massive improvement.
Deriving Binet's formula is a matter of solving a difference equation utilizing the same methods we use for solving differential equations.
The Fibonacci sequence is recursively defined as,
\[ f_{n+1} = f_{n} + f_{n-1}\,. \]We can look at the difference equation,
\[ x_{n+1} = x_{n} + x_{n-1}\,. \]We'll solve this as we would a differential equation and re-write it as,
\[ \lambda^{n+1} = \lambda^{n} + \lambda^{n-1}\,, \]and subsequently we'll divide each term by \(\lambda^{n-1}\),
\[ \frac{\lambda^{n+1}}{\lambda^{n-1}} = \frac{\lambda^{n}}{\lambda^{n-1}} + \frac{\lambda^{n-1}}{\lambda^{n-1}} \implies \lambda^2 - \lambda - 1 = 0 \]We can apply the quadratic formula to yield,
\[ \Phi = \frac{1 + \sqrt{5}}{2}\,,\ -\phi = \frac{1 - \sqrt{5}}{2}\,, \]which are solutions to our original equation,
\[ f_{n} = c_1\Phi^n + c_2(-\phi)^n\,. \]Since \(f_1 = 1\) and \(f_2 = 1\), \(f_0 = f_2 - f_1 = 0\), therefore,
Applying the initial values of \(f_1\) and \(f_2\) yields the system,
\[\begin{aligned} c_1 + c_2 &= 0 \\ c_1\Phi - c_2\phi &= 1\end{aligned}\]Since \(c_1 = -c_2\), we find ourselves with \(c_1(\Phi + \phi) = 1\).
It happens that \(\Phi + \phi = \sqrt{5}\), so \(c_1 = \frac{1}{\sqrt{5}}\) and \(c_2 = -\frac{1}{\sqrt{5}}\) which yields \((1)\) otherwise known as Binet's Formula.