Re: idiom for recurrence relations

• To: mathgroup at smc.vnet.net
• Subject: [mg26746] Re: [mg26714] idiom for recurrence relations
• From: Daniel Lichtblau <danl at wolfram.com>
• Date: Fri, 19 Jan 2001 02:14:20 -0500 (EST)
• References: <200101180557.AAA16636@smc.vnet.net>
• Sender: owner-wri-mathgroup at wolfram.com

```maharrix at my-deja.com wrote:
>
> I have a problem which is mostly esthetic, but could have efficiency
> repercussions. First, the general question: What is the best way to
> implement recurrence relations (e.g. Fibonacci[n])?
>
> The obvious and simplest way is to define base cases and an inductive
> case (I'll use Fibonnaci as the generic example:
>
>   F[0]:=0; F[1]:=1; F[n_]:=F[n-1]+F[n-2]
>
> ). This can be sped up using memoization F[n_]:=F[n]=... Great, but this
> fills up internal memory (nothing can get garbage collected since all
> the computed values must hang around forever; this can be bad when you
> have lots of recurrence relations, graphics and a long session).
>
> So the next solution would be to implement it as an imperative program
> keeping only the necessary values to compute the next in the sequence
> (e.g.
>
>   While[i<=n,temp=f2; f2=f2+f1; f1=temp;i++]
>
> ). Again, great, but  the Mathematica docs repeatedly say that list operations
> are much faster than imperative ones.
>
> So my solution is to use Nest:
>
> F[n_] := First@Nest[{#[[2]], #[[2]]+#[[1]]} &, {0, 1}, n]];
>
> Is this the best (fastest, simplest, most appropriate, or most memory
> efficient, etc.) way to do it? It seems somehow ugly (too many [[, #,
> etc.) but it is shorter in number of characters to encode it (not what I
> would usually consider highest on my list of things to optimize (but
> that's not to say terribly low either))
>
> The whole point is to be able to create recurrence relations easily in a
> single function efficiently. I'd also like to create them as lambdda
> functions (as above, as Function[pair, {pair[[2]],pair[[2]]+pair[[1]]}]
> )
>
> So, how do the above rate and what are the reasonable alternatives?
>
> Thanks for any help,
> Mitch
>
> Sent via Deja.com
> http://www.deja.com/

One good way might be via matrix powers. Below I show several variations
of the Fibonacci example, including one that uses such powers.

fib1[0]:=0; fib1[1]:=1; fib1[n_]:=fib1[n-1]+fib1[n-2];

fib2[0]=0; fib2[1]=1; fib2[n_]:=fib2[n]=fib2[n-1]+fib2[n-2];

fib3[0]=0; fib3[1]=1;
fib3[n_] := Module[{f1=1,f2=1,temp},
Do[temp=f2; f2=f2+f1; f1=temp, {n-2}];
f2]

fib4[n_] := First@Nest[{#[[2]], #[[2]]+#[[1]]} &, {0, 1}, n];

fib5[n_] := MatrixPower[{{1,1},{1,0}},n-1][[1,1]]

fib6[n_] := Fibonacci[n]

\$RecursionLimit=Infinity;

The first two are basically useless for serious work. First is
exponential in speed, second is O(n^2) (see below for reason) but seems
to have recursion stack problems at modest values.

In[9]:= Timing[fib1[25];]
Out[9]= {7.33 Second, Null}

In[10]:= Timing[fib2[2000];]
Out[10]= {0.41 Second, Null}

Third and fourth are O(n^2) algorithms. One factor is for the #steps,
the other is for adding numbers with that many bits (which is what we do
because fib[n]~GoldenRatio^n has O(n) bits).

In[12]:= Timing[fib3[100000];]
Out[12]= {20.04 Second, Null}

In[13]:= Timing[fib4[100000];]
Out[13]= {20.42 Second, Null}

Fifth method uses a divide-and-conquer approach in MatrixPower. So it
will be O(log(n)*M(n)) where M(n) is cost of multiplying a pair of n-bit
integers. The built-in Fibonacci code does something similar, I believe.

In[14]:= Timing[fib5[100000];]
Out[14]= {0.32 Second, Null}

In[15]:= Timing[fib6[100000];]
Out[15]= {0.23 Second, Null}

Daniel Lichtblau
Wolfram Research

```

• Prev by Date: Re: Q: Factor with Polynominals?
• Next by Date: gear drawing
• Previous by thread: idiom for recurrence relations
• Next by thread: Re: idiom for recurrence relations