Leonardo Numbers

This post is part of a series on Mohammad Anwar’s excellent Weekly Challenge, where hackers submit solutions in Perl, Raku, or any other language, to two different challenges every week. (It’s a lot of fun, if you’re into that sort of thing.)

Happy new year! We are on Week 41, and this is Challenge #2.

The Leonardo Numbers (A001595) are a simple recursively defined sequence:

\(
L(n) = \begin{cases}
1 & \text{if } n \lt 2 \\
1 + L(n – 1) + L(n – 2) & \text{if } n \geq 2
\end{cases}\)

The sequence starts: 1, 1, 3, 5, 9, 15, 25, 41, 67, 109, 177, 287, 465, 753, 1219, 1973, 3193, 5167, 8361, 13529, 21891, …

You’ll note this sequence is very similar to the well-known Fibonacci sequence, which differs only in that the Fibonacci sequence does not have the + 1 term, and starts at F(0) = 0.

These are classic problems for beginning students in computer science, as they introduce recursion, and how incredibly simple algorithms can have very high complexity. A naïve implementation of either sequence will quickly bring a computer to its knees. L(35) takes my VM around 10 seconds to compute. L(36) takes 16 seconds. By L(40) we’re into minutes per number.

In general, a recursive implementation of either sequence has a big-O complexity of O(2n), and has a good chance of running your computer out of memory with the deep recursion. So how do we improve that?

Analytic Solution

For the truly cheeky among us, we could solve this problem with the closed-form expression for Leonardo numbers, which gives us an O(1) formula:

\(L(n) = 2 \frac{\varphi^{n+1} \,-\, \psi^{n+1}}{\varphi \,-\, \psi} – 1\)

Where the golden ratio \(\varphi = \frac{1 + \sqrt{5}}{2}\) and \(\psi = \frac{1 – \sqrt{5}}{2}\). See Wolfram Mathworld for a very detailed discussion.

However, I’m not going to use the closed-form solution, because this is a programming challenge, not a number theory challenge.

Memoization

Both my Perl 5 and Raku solutions are going to use memoization. Memoization, informally, means caching results of previous function calls. Memoization only works on pure functions (that is, functions with no side effects), and the Leonardo sequence function is a great example of that. All of the complexity from the recursive solution comes from repeatedly calculating the same results over and over again. By storing results we’ve already calculated, the complexity goes from O (2n) to O(n). Huge win, tiny effort.

Perl 5

Perl 5 has a core module aptly named Memoize. To use it, we start with the slow recursive function, and add the call to memoize:

use Memoize;
sub leo {
    my $n = shift;
    return 1 if $n < 2;
    1 + leo($n - 1) + leo($n - 2);
}
memoize 'leo';

What memoize does is replaces leo in the symbol table with a wrapper that saves each result to a hash, keyed on the argument(s) to the function. The speedup is dramatic indeed.

If you fancy doing the memoization yourself, it’s really easy in this case:

{
    my @leo = (1, 1);
    sub leo_my_memo {
        my $n = shift;
        $leo[$n] //= 1 + leo_my_memo($n - 1) + leo_my_memo($n - 2);
    }
}

The Memoize module has a lot of interesting features, and that introduces some overhead, which is why leo_my_memo runs ~750% faster than the version using Memoize, even though the big-O complexity is the same. Since the code is equally simple either way, in this case, I’d go with leo_my_memo.

Iterative (Pre-Compute) Approach

There is still another fast way to solve this challenge. Since we are asked to print out the Leonardo sequence in order, we can simply iterate through the list, remembering the last two values. This makes the O(n) complexity very apparent, and is even faster than the memoized versions, thanks to the lower overhead.

The iterative approach works for this challenge, but may not be a convenient option if you need arbitrary Leonardo numbers in random order. However, even then, we can pre-compute all of the values we need and store them in an array:

sub leo_to_n {
    my @r = (1, 1);
    push @r, $r[-1] + $r[-2] + 1 for 1..$_[0]-1;
    @r;
}

On my system, for L(1..20), this performs about 50% faster than leo_my_memo, but performance is often not the most important consideration.

Performance Summary

Here is how all of the methods so far stack up in a benchmark:

            Rate no_memo Memoize my_memo    to_n
no_memo   52.0/s      --   -100%   -100%   -100%
Memoize  21151/s  40581%      --    -88%    -92%
my_memo 180696/s 347442%    754%      --    -31%
to_n    263463/s 506630%   1146%     46%      --

Keep in mind this only looks at L(0) through L(20). no_memo‘s performance will be astonishingly bad for larger L(n).

Raku Solution

Raku’s syntax once again provides some very expressive ways to solve this problem. For starters, we can emulate Perl5’s Memoize functionality with Raku’s experimental :cached feature:

use experimental :cached;
sub leo( Int $n ) is cached { $n < 2 ?? 1 !! 1 + leo($n - 1) + leo($n - 2) }

Note that, again, the subroutine code simply a recursive implementation, but we’ve added the is cached trait to the routine to cache the results. That does help a bit, but the results were somewhat disappointing: it’s only about twice as fast as the pure recursive definition. Please keep in mind :cached is still an experimental feature; I trust the Raku team still have some tricks up their collective sleeves to improve performance. For now, let’s roll our own again:

sub leo_my_memo( Int $n ) {
    state @leo = (1, 1);
    @leo[$n] //= 1 + leo_my_memo($n - 1) + leo_my_memo($n - 2);
}

Much better. This runs about 3000% faster than the naïve method.

However, both of these solutions are really just Raku translations of the Perl 5 code. There’s nothing wrong with them, but we’d be remiss if we didn’t explore the full range of Raku features. Raku gives us powerful lazy evaluation with the sequence operator (...), so we can define an array of all Leonardo numbers as follows:

my @leo = 1, 1, 1+*+* ... ∞;

It’s hard not to like lazy sequences such as this. The 1+*+* is the bit that adds up consecutive terms (plus 1). Now if we evaluate, say, @leo[20], Raku will evaluate only what it needs to in order to give us the result. We could of course wrap a subroutine around this if we want, but I opted not to, since it was already so simple.

Wrapping Up

Leonardo numbers (and their close friend, the Fibonacci numbers) have a beautifully simple definition that gives rise to an often surprising amount of interesting math and algorithmics. The hard part of this challenge was deciding what to include in this article, of the reams of interesting information on the topic. I’ve tried to focus on the interesting Perl 5 and Raku bits that emerge from this challenge. I hope you’ve enjoyed the ride.

Leave a Reply

Your email address will not be published. Required fields are marked *