Project Euler 2 - Fibonacci sequence

Question

I lately spend a lot of time at Project Euler programming challenges. I have following problem which I should solve:

Each new term in the Fibonacci sequence is generated by adding the previous two terms. By starting with 1 and 2, the first 10 terms will be:

1, 2, 3, 5, 8, 13, 21, 34, 55, 89, ...

By considering the terms in the Fibonacci sequence whose values do not exceed four million, find the sum of the even-valued terms.

I wrote Java code and I submitted it.

Implementation

First, I wrote helper method fibonacci(int n) which calculates a n-th word of fibonacci sequence and returns it. In main method I created a list of numbers and a n variable, that would be incremented in a loop. In the loop I was looking for all numbers that are even-valued that are added to the list. Next, I created result variable and I sum all numbers in the list and then, finally I printed it out to console.

Main.java:

package pl.hubot.projecteuler.problem2;

import java.util.ArrayList;
import java.util.List;

public class Main {
    public static void main(String[] args) {
        List<Integer> numbers = new ArrayList<>();
        int n = 0;
        while (fibonacci(n) < 4_000_000) {
            if (fibonacci(n) % 2 == 0) numbers.add(fibonacci(n));
            n++;
        }
        int result = 0;
        for (Integer number : numbers) {
            result += number;
        }
        System.out.print(result);
    }

    private static int fibonacci(int n) {
        if (n == 0) return 0;
        else if (n == 1) return 1;
        else return fibonacci(n - 1) + fibonacci(n - 2);
    }
}

Answer 1

The recursive definition of Fibonacci is an elegant one, but in practice it's something you would never want to use. Let's say you want to get fibonacci(5). How often does fibonacci get called?

fibonacci(5)
 =    fibonacci(4) + fibonacci(3)
 =   (fibonacci(3) + fibonacci(2)) + (fibonacci(2) + fibonacci(1))
 =   ((fibonacci(2) + fibonacci(1)) + (fibonacci(1) + fibonacci(0)))
   + ((fibonacci(1) + fibonacci(0)) + fibonacci(1))
 =   (((fibonacci(1) + fibonacci(0)) + fibonacci(1)) 
        + (fibonacci(1) + fibonacci(0))
     )
   + ((fibonacci(1) + fibonacci(0)) + fibonacci(1))

That seems a lot for 0 + 1 + 1 + 2 + 3 + 5. Indeed, the recursive variant without memoization has exponential runtime. You have several options to reduce the time:

• use an iterative variant (which leads to \$ \mathcal O(n)\$ instead of \$ \mathcal O(2^n)\$ per call of fibonacci),
• use memoization (which, if done correctly, leads to \$ \mathcal O(n)\$ in worst case and \$ \mathcal O(1)\$ if the number was already calculated).

But that begs the question: do we actually need fibonacci? For a programming challenge, it's enough to calculate the sum in main, while you're calculating the numbers iteratively. If you're looking for performance, that would be the fastest variant (that doesn't use the closed form).

Answer 2

@Zeta already pointed out that constantly calculating fibonacci(n) is a really bad idea.

@kyrill pointed out that every third number will be even. So let's calculate the next 3 fibonacci numbers in a loop and add the third one to the sum as long as it's smaller than 4000000.

private static int fibsum(){
    int f1 = 1;
    int f2 = 1;
    int f3 = 2;
    int sum = 0;
    while(f3 < 4_000_000) {
        sum += f3;
        f1 = f2 + f3;
        f2 = f1 + f3;
        f3 = f1 + f2;
    }
    return sum;
}

Didn't this just solve the question really easy? No recursion, no storing all previously calculated fibonnaci numbers and O(n) computation.

Answer 3

Some theory (closed-form solution & Binet's formula)

Fibonacci sequence is defined by the recurrence relation $$F_n = F_{n-1} + F_{n-2}$$

This is a second-order homogeneous linear recurrence with constant coefficients.

• second-order = you need the last two terms
• homogeneous = each term in the sequence depends only on previous terms, not some other function. Eg. \$F_n = 2F_{n-1} + sin(n)\$ is not homogeneous because of the \$sin\$.
• linear = terms are in first power (eg. \$F_n = F_{n-1}^2\$ is not linear)
• recurrence relation — I guess that's clear
• with constant coefficients = you don't multiply \$F_{n-1}\$ and \$F_{n-2}\$ by something that depends on \$n\$. Eg. \$F_n = sin(n)*F_{n-1}\$ doesn't have constant coefficients.

A sequence defined by a linear recurrence with constant coefficients is called a constant-recursive sequence or a C-finite sequence [1]. The name is not really important, what's important is that every C-finite sequence has a closed form solution.

---

A closed-form solution basically allows you to compute the \$n\$th term of a sequence in constant time. This is quite similar to an analytical solution, although according to wikipedia, there are a few differences. The concept of closed-form or analytical solution is not specifically related to computing terms of a sequence; it generally denotes a solution which can be calculated in a finite number of operations (more info in the wikipedia article).

That means the \$n\$th term of the Fibonacci sequence can be calculated in constant time, using its closed-form solution. This is a big improvement compared to linear time required by the recurrent solution (not the same as recursive!).

The closed-form solution for the Fibonacci sequence is known as Binet's formula, and is of form $$F_n = \frac{\phi^n - \psi^n}{\phi - \psi}$$ where \$\phi = \frac{1+\sqrt5}{2}\$ and \$\psi = \frac{1-\sqrt5}{2}\$.

The number \$\phi\$ (Greek letter phi) is also known as golden ratio and has some interesting properties (see the linked wikipedia article).

Given the values of \$\phi\$ and \$\psi\$, the formula can be simplified to $$F_n = \frac{\phi^n - \psi^n}{\sqrt5}$$

We can also make use of the fact that \$\left| \frac{\psi^n}{\sqrt5} \right| \lt 0.5\$ if \$n \ge 0\$. That means we can calculate just \$\frac{\phi^n}{\sqrt5}\$ and then round it. The rounding can be accomplished by adding \$0.5\$ and taking just the integral part. $$F_n = \left\lfloor \frac{\phi^n}{\sqrt5}+\frac{1}{2} \right\rfloor$$

Here's how it would look like in Java:

static int fib(int n) {
    double sqrt_5 = Math.sqrt(5);
    double phi = (1 + sqrt_5) / 2;

    double fib_n = Math.pow(phi, n) / sqrt_5;
    return (int) (fib_n + 0.5);
}

---

EDIT

@Zeta has pointed out one fact, which I did not consider relevant and thus I didn't mention it. If you implemented the original Binet's formula \$\frac{\phi^n-\psi^n}{\sqrt5}\$, and your function returned double or float instead of an integer type such as int or long, you would see that the result is incorrect due to floating-point errors. Here's an illustration of how incorrect the result would be:

actual value [long]   Binet's formula [double]
3                     3.0000000000000004
5                     5.000000000000001
8                     8.000000000000002
⋮                    ⋮
498454011879264       498454011879265.2

\$498\ 454\ 011\ 879\ 264\ \$ is the first number for which the integral part of both results is not equal. If the rounding version is used, the result is incorrect starting with the previous term \$308\ 061\ 521\ 170\ 129\$.

This is not a problem in our scenario, since we only need Fibonacci numbers which are less or equal to \$4\ 000\ 000\$. Moreover the function given above returns an int, for which the maximum possible value is \$2^{31}-1\$, which is \$2\ 147\ 483\ 647\$ — well out of danger.

I would post a link to a now-classical paper by David Goldberg, What Every Computer Scientist Should Know About Floating-Point Arithmetic. But since it's not just some easy wikipedia-style reading, I doubt anybody would actually read it. It's no problem to google it if you're interested.

---

Some more theory (even Fibonacci numbers)

• The Fibonacci sequence begins with \$0\$ (even) and \$1\$ (odd).
• Each next number is the sum of the two previous numbers.
• Adding even and odd number yields an odd number, ergo the 3rd term should be odd. Indeed, \$1\$ is odd.
• Adding two odd numbers yields an even number, ergo the 4th term should be even. Voila, \$2\$ is even.
• Then we have \$1+2=3\$, an odd number. And we're right where we started – even number (\$2\$) followed by odd (\$3\$).
• Via induction it's clear that every 3rd Fibonacci number is even.

Even more theory (exponential sums)

We know how to calculate \$n\$th Fibonacci number efficiently, and we know which Fibonacci numbers are even. Now we just need to sum them.

In a naive approach, we could calculate the terms using the recurrence formula and add them to the sum as we calculate them. For that, we would need something like this pseudocode:

sum := 0
a, b := 0, 1

while a < 4_000_000:
    if a is even:
        sum := sum + a
    a, b := b, a+b

We could use the fact that every 3rd Fibonacci number is even, but there would be little gain. We would save a modulo operation, but on most modern architectures integer division is quite fast anyway. Moreover, modulo by a power of 2 will be probably optimized into an AND by the compiler.

Now if we were to sum every 3rd Fibonacci number, and the numbers were to be calculated using the closed-form solution, it would be most likely less efficient than the simple code above. This is because the simple code performs just addition. That's is. The closed-form solution, however, requires raising \$\phi\$ to the \$n\$th power and dividing the result. And that's in the optimal case, when we save all the constants statically and don't need to recalculate them on each call to the function.

Instead, we can calculate the entire sum in constant time.

To calculate the sum, we must first write it down: $$\sum_{n=0}^{\frac{N}{3}}{ \frac{(\phi^3)^n}{\sqrt5} }$$

\$N\$ is the index of the last Fibonacci number less or equal to 4,000,000. The sum goes up to \$\frac{N}{3}\$ because we only want every 3rd number, eg. \$\frac{N}{3}\$ numbers total. How to obtain \$N\$ will be discussed later.

The sum can be also written as $$\frac{1}{\sqrt5} \; \sum_{n=0}^{\frac{N}{3}}{(\phi^3)^n}$$

Now we have a sum of an exponential sequence. This sum can be calculated using the formula $$\sum_{n=0}^{N-1}{r^n} = \frac{1 - r^N}{1 - r}$$

In our case that is $$\sum_{n=0}^{\frac{N}{3}} {(\phi^3)^n} = \frac{ 1 - (\phi^3)^{\frac{N}{3} + 1} } { 1 - \phi^3 } $$

The entire sum then is $$\frac{1}{\sqrt5} \frac{ 1 - (\phi^3)^{\frac{N}{3} + 1} } { 1 - \phi^3 } $$

Written In Java:

public static int sumEvenFibs(int numFibs) {
    int numTerms = numFibs / 3;

    double sqrt_5 = Math.sqrt(5);
    double phi_to3 = Math.pow((1 + sqrt_5) / 2, 3);
    double phi_to3_toN = Math.pow(phi_to3, numTerms+1);
    double sum_phi = (1 - phi_to3_toN) / (1 - phi_to3);

    return (int) (sum_phi / sqrt_5);
}

The last bit of theory (given fib(n) determine n)

We have to sum all even Fibonacci numbers \$F_n\$ such that \$F_n \le 4\ 000\ 000\$. In order to use the exponential sum formula, we have to know how many terms we want to sum, not just "all that are less or equal 4,000,000".

Given a number \$x\$ such that \$F_n \le x \lt F_{n+1}\$ for two successive Fibonacci numbers \$F_n\$ and \$F_{n+1}\$, we have to determine \$n\$.

In other words, given some number \$x\$, we have to find the index \$n\$ of the closest Fibonacci number \$F_n\$, where \$F_n\$ is less or equal \$x\$.

For that, we have to take a look at the closed-form formula and do some math. $$F_n = \frac{\phi^n - \psi^n}{\sqrt5}$$ $$\sqrt5 \, F_n = \phi^n - \psi^n$$ $$\sqrt5 \, F_n + \psi^n = \phi^n$$ $$\log(\sqrt5 \, F_n + \psi^n) = \log \phi^n$$ $$\log(\sqrt5 \, F_n + \psi^n) = n \log \phi$$ $$\frac{\log(\sqrt5 \, F_n + \psi^n)}{\log \phi} = n $$ $$\log_\phi(\sqrt5 \, F_n + \psi^n) = n$$ $$\left\lfloor \log_\phi(\sqrt5 \, F_n + \frac{1}{2}) \right\rfloor = n$$

So now we know how to calculate the index of the last Fibonacci number in our sum.

Note that \$1\$ occurs in the sequence twice, but obviously this method can produce only one index.

Java code:

public static int getNumOfFib(int fib) {
    double sqrt_5 = Math.sqrt(5);
    double phi = (1 + sqrt_5) / 2;
    double naturalLog = Math.log(sqrt_5 * fib + 1/2);

    return (int) (naturalLog / Math.log(phi));
}

Code

So that's how you sum all even Fibonacci numbers up to 4,000,000 in constant time.

Here's the full code.

public class Main {

    static final double sqrt_5 = Math.sqrt(5);
    static final double phi = (1 + sqrt_5) / 2;

    public static void main(String[] args) {
        System.out.print(sumEvenFibs(getNumOfFib(4_000_000)));
    }

    public static int getNumOfFib(int fib) {
        final double naturalLog = Math.log(sqrt_5 * fib + 1/2);
        return (int) (naturalLog / Math.log(phi));
    }

    public static int sumEvenFibs(int numFibs) {
        final int numTerms = numFibs / 3;

        final double
          phi_to3 = Math.pow(phi, 3),
          phi_to3_toN = Math.pow(phi_to3, numTerms+1),
          sum_phi = (1 - phi_to3_toN) / (1 - phi_to3);

        return (int) (sum_phi / sqrt_5);
    }
}