Software Development

Fibonacci and Lucas Sequences

This posts touches on three of my favorite topics – math, transferring knowledge through experience (tutorial unit tests) and the importance of research.

Most developers are aware of the Fibonacci sequence, mostly through job interviews.

To briefly recap the series is defined a:

F(n) = F(n-1) + F(n-2), n > 2
F(1) = F(2) = 1

 
There’s a variant definition:

F(n) = F(n-1) + F(n-2), n > 1
F(1) = 1
F(0) = 0

There are four well-known solutions to the white-board question “write code to calculate F(n)”.

Recursion – you need to mention this to show that you’re comfortable with recursion but you must also mention that it’s a Really Bad Idea since it requires O(2n) time and space stack since you double the work for each n.

Recursion with memoization – this can be a good approach if you point out it’s a good generalization. Basically it’s recursion but you maintain a cache (the memoization) so you only need to make the recursive call once – subsequent recursive calls just look up the cached value. This is a flexible technique since it can be used for any pure recursive function. (That is, a recursive function that depends solely on its inputs and has no side effects.) The first calls require O(n) time, stack and heap space. I don’t recall if it matters if you do the recursive call on the smaller or larger value first.

If you have a persistent cache subsequent calls require O(1) time and stack space and O(n) heap space.

Iteration – if you can’t cache the values (or just want to efficiently initialize a cache) you can use an iterative approach. It requires O(n) time but only O(1) stack and heap space.

Direct approximation – finally there is a well-known approximation using φ, or a variant using sqrt(5). It is O(1) for time, stack space, and heap space. It’s a good approach if you 1) use a lookup table for the smallest values and 2) make sure n is not too big.

The last point is often overlooked. The approximation only works as long as you don’t exceed the precision of your floating point number. F(100,000) should be good. F(1,000,000,000,000) may not be. The iterative approach isn’t practical with numbers this large.

Research

Did you know there’s two other solutions with performance O(lg(n)) (per Wikipedia) in time and space? (I’m not convinced it’s O(lg(n)) since it’s not a divide-and-conquer algorithm – the two recursive calls do not split the initial work between them – but with memoization it’s definitely less than O(n). I suspect but can’t quickly prove it’s O(lg2(n)).)

Per Wikipedia we know:

F(2n-1) = F2(n) + F2(n-1)
F(2n) = F(n)(F(n) + 2F(n-1))

It is straightforward to rewrite this as a recursive method for F(n).

There is another property that considers three cases – F(3n-2), F(3n-1) and F(3n). See the code for details.

These sites provide many additional properties of the Fibonacci and related Lucas sequences. Few developers will ever need to know these properties but in those rare cases an hour of research can save days of work.

Implementation

We can now use our research to implement suitable methods for the Fibonacci and Lucas sequences.

Fibonacci calculation

(This code does not show an optimization using direct approximation for uncached values for sufficiently small n.)

/**
     * Get specified Fibonacci number.
     * @param n
     * @return
     */
    @Override
    public BigInteger get(int n) {
        if (n < 0) {
            throw new IllegalArgumentException("index must be non-negative");
        }

        BigInteger value = null;

        synchronized (cache) {
            value = cache.get(n);

            if (value == null) {
                int m = n / 3;

                switch (n % 3) {
                case 0:
                    value = TWO.multiply(get(m).pow(3))
                               .add(THREE.multiply(get(m + 1)).multiply(get(m))
                                         .multiply(get(m - 1)));

                    break;

                case 1:
                    value = get(m + 1).pow(3)
                                .add(THREE.multiply(get(m + 1)
                                                        .multiply(get(m).pow(2))))
                                .subtract(get(m).pow(3));

                    break;

                case 2:
                    value = get(m + 1).pow(3)
                                .add(THREE.multiply(get(m + 1).pow(2)
                                                        .multiply(get(m))))
                                .add(get(m).pow(3));

                    break;
                }

                cache.put(n, value);
            }
        }

        return value;
    }

Fibonacci Iterator

/**
     * ListIterator class.
     * @author bgiles
     */
    private static final class FibonacciIterator extends ListIterator {
        private BigInteger x = BigInteger.ZERO;
        private BigInteger y = BigInteger.ONE;

        public FibonacciIterator() {
        }

        public FibonacciIterator(int startIndex, FibonacciNumber fibonacci) {
            this.idx = startIndex;
            this.x = fibonacci.get(idx);
            this.y = fibonacci.get(idx + 1);
        }

        protected BigInteger getNext() {
            BigInteger t = x;
            x = y;
            y = t.add(x);

            return t;
        }

        protected BigInteger getPrevious() {
            BigInteger t = y;
            y = x;
            x = t.subtract(x);

            return x;
        }
    }

Lucas calculation

/**
     * Get specified Lucas number.
     * @param n
     * @return
     */
    public BigInteger get(int n) {
        if (n < 0) {
            throw new IllegalArgumentException("index must be non-negative");
        }

        BigInteger value = null;

        synchronized (cache) {
            value = cache.get(n);

            if (value == null) {
                value = Sequences.FIBONACCI.get(n + 1)
                                           .add(Sequences.FIBONACCI.get(n - 1));
                cache.put(n, value);
            }
        }

        return value;
    }

Lucas iterator

/**
     * ListIterator class.
     * @author bgiles
     */
    private static final class LucasIterator extends ListIterator {
        private BigInteger x = TWO;
        private BigInteger y = BigInteger.ONE;

        public LucasIterator() {
        }

        public LucasIterator(int startIndex, LucasNumber lucas) {
            idx = startIndex;
            this.x = lucas.get(idx);
            this.y = lucas.get(idx + 1);
        }

        protected BigInteger getNext() {
            BigInteger t = x;
            x = y;
            y = t.add(x);

            return t;
        }

        protected BigInteger getPrevious() {
            BigInteger t = y;
            y = x;
            x = t.subtract(x);

            return x;
        }
    }

Education

What is the best way to educate other developers about the existence of these unexpected relationships? Code, of course!

What is the best way to educate other developers about the existence of code that demonstrates these relationships? Unit tests, of course!

It is straightforward to write unit tests that simultaneous verify our implementation and inform other developers about tricks they can use to improve their code. The key is to provide a link to additional information.

Fibonacci Sequence

public class FibonacciNumberTest extends AbstractRecurrenceSequenceTest {
    private static final BigInteger MINUS_ONE = BigInteger.valueOf(-1);

    /**
     * Constructor
     */
    public FibonacciNumberTest() throws NoSuchMethodException {
        super(FibonacciNumber.class);
    }

    /**
     * Get number of tests to run.
     */
    @Override
    public int getMaxTests() {
        return 300;
    }

    /**
     * Verify the definition is properly implemented.
     *
     * @return
     */
    @Test
    @Override
    public void verifyDefinition() {
        for (int n = 2; n < getMaxTests(); n++) {
            BigInteger u = seq.get(n);
            BigInteger v = seq.get(n - 1);
            BigInteger w = seq.get(n - 2);
            Assert.assertEquals(u, v.add(w));
        }
    }

    /**
     * Verify initial terms.
     */
    @Test
    @Override
    public void verifyInitialTerms() {
        verifyInitialTerms(Arrays.asList(ZERO, ONE, ONE, TWO, THREE, FIVE, EIGHT));
    }

    /**
     * Verify that every third term is even and the other two terms are odd.
     * This is a subset of the general divisibility property.
     *
     * @return
     */
    @Test
    public void verifyEvenDivisibility() {
        for (int n = 0; n < getMaxTests(); n += 3) {
            Assert.assertEquals(ZERO, seq.get(n).mod(TWO));
            Assert.assertEquals(ONE, seq.get(n + 1).mod(TWO));
            Assert.assertEquals(ONE, seq.get(n + 2).mod(TWO));
        }
    }

    /**
     * Verify general divisibility property.
     *
     * @return
     */
    @Test
    public void verifyDivisibility() {
        for (int d = 3; d < getMaxTests(); d++) {
            BigInteger divisor = seq.get(d);

            for (int n = 0; n < getMaxTests(); n += d) {
                Assert.assertEquals(ZERO, seq.get(n).mod(divisor));

                for (int i = 1; (i < d) && ((n + i) < getMaxTests()); i++) {
                    Assert.assertFalse(ZERO.equals(seq.get(n + i).mod(divisor)));
                }
            }
        }
    }

    /**
     * Verify the property that gcd(F(m), F(n)) = F(gcd(m,n)). This is a
     * stronger statement than the divisibility property.
     */
    @Test
    public void verifyGcd() {
        for (int m = 3; m < getMaxTests(); m++) {
            for (int n = m + 1; n < getMaxTests(); n++) {
                BigInteger gcd1 = seq.get(m).gcd(seq.get(n));
                int gcd2 = BigInteger.valueOf(m).gcd(BigInteger.valueOf(n))
                                     .intValue();
                Assert.assertEquals(gcd1, seq.get(gcd2));
            }
        }
    }

    /**
     * Verify second identity (per Wikipedia): sum(F(i)) = F(n+2)-1
     */
    @Test
    public void verifySecondIdentity() {
        BigInteger sum = ZERO;

        for (int n = 0; n < getMaxTests(); n++) {
            sum = sum.add(seq.get(n));
            Assert.assertEquals(sum, seq.get(n + 2).subtract(ONE));
        }
    }

    /**
     * Verify third identity (per Wikipedia): sum(F(2i)) = F(2n+1)-1 and
     * sum(F(2i+1)) = F(2n)
     */
    @Test
    public void verifyThirdIdentity() {
        BigInteger sum = ZERO;

        for (int n = 0; n < getMaxTests(); n += 2) {
            sum = sum.add(seq.get(n));
            Assert.assertEquals(sum, seq.get(n + 1).subtract(ONE));
        }

        sum = ZERO;

        for (int n = 1; n < getMaxTests(); n += 2) {
            sum = sum.add(seq.get(n));
            Assert.assertEquals(sum, seq.get(n + 1));
        }
    }

    /**
     * Verify fourth identity (per Wikipedia): sum(iF(i)) = nF(n+2) - F(n+3) + 2
     */
    @Test
    public void verifyFourthIdentity() {
        BigInteger sum = ZERO;

        for (int n = 0; n < getMaxTests(); n++) {
            sum = sum.add(BigInteger.valueOf(n).multiply(seq.get(n)));

            BigInteger x = BigInteger.valueOf(n).multiply(seq.get(n + 2))
                                     .subtract(seq.get(n + 3)).add(TWO);
            Assert.assertEquals(sum, x);
        }
    }

    /**
     * Verify fifth identity (per Wikipedia): sum(F(i)^2) = F(n)F(n+1)
     */
    public void verifyFifthIdentity() {
        BigInteger sum = ZERO;

        for (int n = 0; n < getMaxTests(); n += 2) {
            BigInteger u = seq.get(n);
            BigInteger v = seq.get(n + 1);
            sum = sum.add(u.pow(2));
            Assert.assertEquals(sum, u.multiply(v));
        }
    }

    /**
     * Verify Cassini's Identity - F(n-1)F(n+1) - F(n)^2 = -1^n
     */
    @Test
    public void verifyCassiniIdentity() {
        for (int n = 2; n < getMaxTests(); n += 2) {
            BigInteger u = seq.get(n - 1);
            BigInteger v = seq.get(n);
            BigInteger w = seq.get(n + 1);

            BigInteger x = w.multiply(u).subtract(v.pow(2));
            Assert.assertEquals(ONE, x);
        }

        for (int n = 1; n < getMaxTests(); n += 2) {
            BigInteger u = seq.get(n - 1);
            BigInteger v = seq.get(n);
            BigInteger w = seq.get(n + 1);

            BigInteger x = w.multiply(u).subtract(v.pow(2));
            Assert.assertEquals(MINUS_ONE, x);
        }
    }

    /**
     * Verify doubling: F(2n-1) = F(n)^2 + F(n-1)^2 and F(2n) =
     * F(n)(F(n-1)+F(n+1)) = F(n)(2*F(n-1)+F(n).
     */
    @Test
    public void verifyDoubling() {
        for (int n = 1; n < getMaxTests(); n++) {
            BigInteger u = seq.get(n - 1);
            BigInteger v = seq.get(n);
            BigInteger w = seq.get(n + 1);

            BigInteger x = v.multiply(v).add(u.pow(2));
            Assert.assertEquals(seq.get((2 * n) - 1), x);

            x = v.multiply(u.add(w));
            Assert.assertEquals(seq.get(2 * n), x);

            x = v.multiply(v.add(TWO.multiply(u)));
            Assert.assertEquals(seq.get(2 * n), x);
        }
    }

    /**
     * Verify tripling.
     */
    @Test
    public void verifyTripling() {
        for (int n = 1; n < getMaxTests(); n++) {
            BigInteger u = seq.get(n - 1);
            BigInteger v = seq.get(n);
            BigInteger w = seq.get(n + 1);

            BigInteger x = TWO.multiply(v.pow(3))
                              .add(THREE.multiply(v).multiply(u).multiply(w));
            Assert.assertEquals(seq.get(3 * n), x);

            x = w.pow(3).add(THREE.multiply(w).multiply(v.pow(2)))
                 .subtract(v.pow(3));
            Assert.assertEquals(seq.get((3 * n) + 1), x);

            x = w.pow(3).add(THREE.multiply(w.pow(2)).multiply(v)).add(v.pow(3));
            Assert.assertEquals(seq.get((3 * n) + 2), x);
        }
    }
}

Lucas Sequence

public class LucasNumberTest extends AbstractRecurrenceSequenceTest {
    private static final FibonacciNumber fibonacci = new FibonacciNumber();

    /**
     * Constructor
     */
    public LucasNumberTest() throws NoSuchMethodException {
        super(LucasNumber.class);
    }

    /**
     * Get number of tests to run.
     */
    @Override
    public int getMaxTests() {
        return 300;
    }

    /**
     * Verify the definition is properly implemented.
     *
     * @return
     */
    @Test
    @Override
    public void verifyDefinition() {
        for (int n = 2; n < getMaxTests(); n++) {
            BigInteger u = seq.get(n);
            BigInteger v = seq.get(n - 1);
            BigInteger w = seq.get(n - 2);
            Assert.assertEquals(u, v.add(w));
        }
    }

    /**
     * Verify initial terms.
     */
    @Test
    @Override
    public void verifyInitialTerms() {
        verifyInitialTerms(Arrays.asList(TWO, ONE, THREE, FOUR, SEVEN, ELEVEN,
                BigInteger.valueOf(18), BigInteger.valueOf(29)));
    }

    /**
     * Verify Lucas properties.
     */
    @Test
    public void verifyLucas() {
        // L(n) = F(n-1) + F(n+1)
        for (int n = 2; n < getMaxTests(); n++) {
            Assert.assertEquals(seq.get(n),
                fibonacci.get(n - 1).add(fibonacci.get(n + 1)));
        }
    }

    /**
     *  F(2n) = L(n)F(n)
     */
    @Test
    public void verifyLucas2() {
        for (int n = 2; n < getMaxTests(); n++) {
            Assert.assertEquals(fibonacci.get(2 * n),
                seq.get(n).multiply(fibonacci.get(n)));
        }
    }

    /**
     * F(n) = (L(n-1)+ L(n+1))/5
     */
    @Test
    public void verifyLucas3() {
        for (int n = 2; n < getMaxTests(); n++) {
            Assert.assertEquals(FIVE.multiply(fibonacci.get(n)),
                seq.get(n - 1).add(seq.get(n + 1)));
        }
    }

    /**
     * L(n)^2 = 5 F(n)^2 + 4(-1)^n
     */
    @Test
    public void verifyLucas4() {
        for (int n = 2; n < getMaxTests(); n += 2) {
            Assert.assertEquals(seq.get(n).pow(2),
                FIVE.multiply(fibonacci.get(n).pow(2)).add(FOUR));
        }

        for (int n = 1; n < getMaxTests(); n += 2) {
            Assert.assertEquals(seq.get(n).pow(2),
                FIVE.multiply(fibonacci.get(n).pow(2)).subtract(FOUR));
        }
    }
}

Conclusion

Obviously developers rarely need to compute Fibonacci numbers unless they’re working on Project Euler problems or at a job interview. This code isn’t going to have direct utility.

At the same time it’s a powerful demonstration of the value of investing an hour or two in research even if you’re sure you already know everything you need to know. You probably don’t need BigInteger implementation but some people might consider the O(lg(n)) approach preferable to the estimate using powers of φ, or could make good use of the relationships discussed on the MathWorld and Wikipedia pages.

Source Code

The good news is that I have published the source code for this… and the bad news is it’s part of ongoing doodling when I’m doing Project Euler problems. (There are no solutions here – it’s entirely explorations of ideas inspired by the problems. So the code is a little rough and should not be used to decide whether or not to bring me in for an interview (unless you’re impressed): http://github.com/beargiles/projecteuler.

Reference: Fibonacci and Lucas Sequences from our JCG partner Bear Giles at the Invariant Properties blog.
Subscribe
Notify of
guest

This site uses Akismet to reduce spam. Learn how your comment data is processed.

0 Comments
Oldest
Newest Most Voted
Inline Feedbacks
View all comments
Back to top button