4

In some cases, recursion in Mathematica (without storing values) is not very efficient.

For example, the recursive implementation of Fibonacci numbers (version which doesn't store the values) :

fib1[1] = fib1[2] = 1;
fib1[n_Integer /; n > 2] := fib1[n - 1] + fib1[n - 2]

fib1[30] // AbsoluteTiming
(* {1.500157, 832040} *)

Of course, the version which stores the values is faster (but consumes memory) :

fib2[1] = fib2[2] = 1;
fib2[n_Integer /; n > 2] := fib2[n] = fib2[n - 1] + fib2[n - 2]

fib2[30] // AbsoluteTiming
(* {0., 832040} *)


I try to optimize (without storing values and using LinearRecurrence) recursive implementations of number sequences $u_{n}$ defined by $u_{n} = f(u_{n-1},\;\ldots,\;u_{n-p})$ where $f(x_{1},\;\ldots,\;x_{p})$ is a function with $p$ arguments.

My idea is to use the //. (ReplaceRepeated) operator instead of standard recursion.

I defined the generic function transformRecurrence by :

transformRecurrence[function_, init_][n_Integer /; n > 0] := Block[
    {index, length = Length@init, terms = Reverse@init},
    ReplaceRepeated[
        {n, init}, 
        {index_ /; index > length, terms_List} :> {index - 1, Flatten@{Rest@terms, function@@terms}},
        MaxIterations -> Infinity
    ] // Last // Last]

I use this function to redefine the Fibonacci numbers :

fib3 = transformRecurrence[Plus, {1, 1}]

Some execution times :

fib3[30] // AbsoluteTiming
(* {0., 832040} *)

fib3[10000]; // AbsoluteTiming
(* {0.031265, Null} *)

fib3[100000]; // AbsoluteTiming
(* {0.328130, Null} *)

fib3[1000000]; // AbsoluteTiming
(* {6.813015, Null} *)

An improved implementation (specialized without list manipulations) is :

fib4[n_Integer /; n > 0] := Block[
    {index, x, y},
    ReplaceRepeated[
        {n, 1, 1}, 
        {index_ /; index > 2, x_, y_} -> {index - 1, y, x + y},
        MaxIterations -> Infinity
    ] // Last]

fib4[30] // AbsoluteTiming
(* {0., 832040} *)

fib4[10000]; // AbsoluteTiming
(* {0.015622, Null} *)

fib4[100000]; // AbsoluteTiming
(* {0.171880, Null} *)

fib4[1000000]; // AbsoluteTiming
(* {4.954125, Null} *)

In comparaison, the built-in function always wins :

Fibonacci[30] // AbsoluteTiming
(* {0., 832040} *)

Fibonacci[10000]; // AbsoluteTiming
(* {0., Null} *)

Fibonacci[100000]; // AbsoluteTiming
(* {0., Null} *)

fib3[1000000]; // AbsoluteTiming
(* {0.031244, Null} *)


And you, have you ever used the //. operator to do such things (even for fun) ?

Olivier
  • 131
  • 2
  • For linear recurrences, matrix power-based methods are very efficient. Re: ReplaceRepeated[], I used it here instead of recursion. – J. M.'s missing motivation Aug 30 '15 at 18:11
  • 1
    I did not, but I have a question. As far as I know, internally any function definition is the application of a definition rule. Consequently, expr //. rule amounts to find the FixedPoint[rule, expr] (where rule is considered a function) and this is exactly what the evaluation loop does (always as far as I understand). So, what would be the conceptual benefit? Perhaps you can clarify. – user8074 Aug 30 '15 at 19:36
  • @user8074 You're right about the evaluation loop, and the correspondence between //. and FixedPoint. The latter is just written in terms of functions (which themselves are of course also rules). So ReplaceRepeated is (in the Mathematica sense) more "low-level", I think. I did a toy example for //. in the context of evaluation here: Showing steps for TrigExpand. – Jens Aug 30 '15 at 20:18

0 Answers0