coding style - Can one identify the design patterns of Mathematica?

... or are they unnecessary in such a high-level language?

I've been thinking about programming style, coding standards and the like quite a bit lately, the result of my current work on a mixed .Net/Mathematica project.

In object-oriented programming, people talk a lot about "design patterns", essentially, standard ways of approaching particular problems (There's even a conference). The idea has been extended to other domains as well, e.g. for personal productivity.

But some people argue that identified "patterns" are really just a sign that the language in question isn't abstracted enough.

So what, if any, are the patterns of Mathematica programming? Or are they not necessary?

Answer

Leaving the general discussion on design patterns for some later time, I will try to collect a comprehensive information on some specific design patterns / idioms, one pattern per post. This time it will be memoization

Basic memoization

What is memoization

Memoization is a technique which enables a program to reuse parts of computations it has done already, to avoid unnecessary repeated work.

Often (although not always) memoization goes hand in hand with recursion. One reason for that is that many recursive functions have very high (often exponential) complexity when implemented straightforwardly, due to the phenomena of overlapping sub-problems. In many such cases, memoization is able to drastically reduce the complexity.

Types of memoization

This is largely a matter of definitions, but, at least in Mathematica, I will distinguish 3 types of memoization:

• 
  

  "Simple" memoization. This is basically a trick used to remember a number of past values on which some computationally demanding function has been computing. This is typically done as

  
  

  f[x_]:=f[x]= r.h.s.
  
  

  
  

  what is important here is that enabling memoization leads to a speed-up, but not to a drastic change of the complexity of the resulting algorithm. Put in other words, this is a technique of remembering past values for functions without the overlapping sub-problems. A typical example here would be to memoize values of some expensive integral computed on some fixed grid in its parameter space.

  
  


• 
  

  Algorithmic memoization: this may be technically realized in the same way, but the distinctive feature of this type of problems is that they have overlapping sub-problems, so here memoization leads to a drastic change in algorithm's complexity. Perhaps the simplest well-known example here is Fibonacci numbers:

  
  

  fib[0] = fib[1] = 1;
  fib[n_Integer?Positive]:= fib[n] = fib[n-1] + fib[n-2]
  

  
  

  where the complexity of the resulting algorithm is reduced from exponential to linear, as a result of memoization.

  
  



• 
  

  Caching. The main difference between caching and memoization is that memoization is a technique of run-time caching during one and the same function's evaluation (single evaluation process), while general caching is a technique to memorize some computed values for possible future uses, which may be happening for different evaluations. So, we can say that caching is a form of "persistent memoization", whose scope extends to more than one evaluation process.

  
  

Why / how it works in Mathematica

The technical reason why Mathematica makes it rather straightforward to implement memoization is that pattern-based definitions are global rules, and more specific rules are applied before more general ones. The rules we create during memoization, are more specific than the rules for the original general definitions, and therefore they are reordered by the system to automatically be tried first. For example, for a function

f[x_]:= f[x]= x^2

When we call it as

f /@ {1,2,3}

we end up with

?f
Global`f
f[1]=1

f[2]=4

f[3]=9


f[x_]:=f[x]=x^2

What I want to stress here is that the simplicity of memoization in Mathematica is due to the fact that we are reusing powerful core system mechanisms.

Advanced memoization

Here I will mention a few more advanced cases and /or illustrate with some examples

Algorithmic memoization example - longest common subsequence

Here I will show a less trivial example of memoization applied to a longest common subsequence problem. Here is the code:

Clear[lcs];
lcs[{}, _List] := {};
lcs[_List, {}] := {};

lcs[x_List, y_List] /; Last[x] === Last[y] := 
   Append[lcs[Most[x], Most[y]], Last[x]];
lcs[x_List, y_List] :=
   lcs[x, y] = (
     If[Length[#1] > Length[#2], #1, #2] &[
       lcs[Most[x], y], lcs[x, Most[y]]
     ]);

The resulting complexity of this algorithm is quadratic in the length of the lists, as it is also for an imperative version of the algorithm. But here the code is quite simple and self-documenting. As an illustration, here is a longest increasing subsequence for some random list:

lst=RandomInteger[{1,500},100];

lcs[lst,Sort@lst]//AbsoluteTiming

(* {0.129883,{30,89,103,160,221,227,236,250,254,297,307,312,330,354,371,374,374}} *)

The implementation can be improved by using linked lists. It also has an annoying problem that one has to manually clear the memoized values - but this issue will be addressed below.

Memoization for functions of more than one parameter, and using patterns in memoized definitions.

So far, we only looked at memoization of functions depending on a single parameter. Often the situation is more complex. One such example is here. Generally, we may want to memoize functions rather than just values. The general idiom for this is the following (I will show one extra parameter, but this is easy to generalize):

f[x_,y_]:=
   Module[{ylocal},
      f[x,ylocal_] = r.h.s (depends on x and ylocal);

      f[x,y]
   ];

What you see here is that the function is first adding a general definition, which is however computed using a fixed x, and then actually computes the value for f[x,y]. In all calls after the first one, with the same x, the newly added general definition will be used. Sometimes, as in the linked example, this may involve additional complications, but the general scheme will be the same.

There are many more examples of various flavors of this technique. For example, in his very cool solution for selecting minimal subsets, Mr.Wizard was able to use a combination of this technique and Orderless attribute in a very powerful way. The code is so short that I will reproduce it here:

minimal[sets_] :=
  Module[{f},
    f[x__] := (f[x, ___] = Sequence[]; {x});
    SetAttributes[f, Orderless];
    f @@@ Sort @ sets

  ] 

This is a more complex example of memoization, where the newly generated definitions involve patterns and this allows for a much more powerful memoization.

Caching, and selective clearing of memoized definitions

Sometimes, particularly for caching (persistent memoization), one may need to clear all or part of the memoized definitions. The basic idiom for doing so was given in this answer by Simon

ClearCache[f_] :=  DownValues[f] = 
     DeleteCases[DownValues[f], _?(FreeQ[First[#], Pattern] &)]

Sometimes, one may need more complex ways of caching, such as e.g. fixed number of cached results. The two implementation I am aware of, which automate this sort of things, are

In the linked posts there are examples of use.

Additional useful techniques / tricks

Self-blocking and automatic removal of memoized definitions

Self-blocking can be thought of as a separate design pattern. Its application to memoization is for cases when one needs the memoized values to be automatically cleared at the end of a computation. I will show two examples - Fibonacci numbers and longest common subsequence - since both of them were described above. Here is how it may look for the Fibonacci numbers:

ClearAll[fib];
fib[n_Integer]:=
  Block[{fib},         
    fib[0] = fib[1] = 1;
    fib[m_Integer?Positive]:= fib[m] = fib[m-1] + fib[m-2];
    fib[n]

  ]

You can see that the main definition redefines fib in the body of Block, and then calls fib[n] inside Block. This will guarantee that we don't generate new global definitions for fib once we exit Block.

For the Fibonacci numbers, this is less important, but for many memoizing functions this will be crucial. For example, for the lcs function, we can do the same thing:

lcs[fst_List, sec_List]:=
  Block[{lcs},
    definitions-of-lcs-given-before; 
    lcs[fst,sec]
  ]

I have used this technique on a number of occasions and find it generally useful. One related discussion is in my answer here. One particularly nice thing about it is that Block guarantees you the automatic cleanup even in cases when exception or Abort[] was thrown during the execution of its body - without any additional effort from the programmer.

Encapsulation of cached definitions

Often one may need to cache some intermediate result, but not expose it directly on the top-level, because the direct users of that result would be some other functions, not the end user. One can then use another pattern of creation of a mutable state by using Module-generated variables.

The main idiom here is

 Module[{cached},
   cached:=cached = some-computation;
   f[...]:= f-body-involving-cached;
   g[...]:= g-body-involving-cached;
 ]

Some examples involve

Memoization of compiled functions, JIT compilation

I will give here a simple example from my post on meta-programming: create a JIT-compiled version of Select, which would compile Select with a custom predicate:

ClearAll[selectJIT];
selectJIT[pred_, listType_] :=
  selectJIT[pred, Verbatim[listType]] = 
    Block[{lst},
     With[{decl = {Prepend[listType, lst]}},
      Compile @@ 

       Hold[decl, Select[lst, pred], CompilationTarget -> "C", 
          RuntimeOptions -> "Speed"]]];

This is a general method however. More examples can be found in

Symbol's auto-loading / self-uncompression

This technique is logically closely related to memoization, although perhaps is not memoization per se. I will reproduce here one function from this post, which takes a symbol and modifies its definition so that it becomes "self-uncompressing" - meaning that it stores a compressed self-modifying definition of itself:

(* Make a symbol with DownValues / OwnValues self - uncompressing *)
ClearAll[defineCompressed];
SetAttributes[defineCompressed, HoldFirst];

defineCompressed[sym_Symbol, valueType_: DownValues] :=
  With[{newVals = 
     valueType[sym] /.
       Verbatim[RuleDelayed][
         hpt : Verbatim[HoldPattern][HoldPattern[pt_]], rhs_] :>
            With[{eval = Compress@rhs}, hpt :> (pt = Uncompress@ eval)]
        },
        ClearAll[sym];
        sym := (ClearAll[sym]; valueType[sym] = newVals; sym)
];

In the mentioned post there are more explanations on how this works.

Links

Here I will list some memoization-related links which came to my mind (some of them I gave above already, but provide here once again for convenience). This list will surely be incomplete, and I invite the community to add more.

• 
  

  General

  
  
  


• 
  
  

  More complex memoization

  
  
  


• 
  

  Algorithmic memoization

  
  
  


• 
  

  Automation of caching

  
  
  
  


• 
  

  Caching and encapsulation of cached values

  
  
  


• 
  

  Memoization of compiled definitions

  
  
  


• 
  
  

  Automatic clean-up of memoized definitions and the self-blocking trick

  
  
  


• 
  

  Memoization and parallel computations

  
  
  


• 
  

  Memoization in probabilistic inference