Can one identify the design patterns of Mathematica?
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
- Szabolcs's implementation
- My own implementation
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
- Caching file import
- Compiling and memoizing
Outerfor an arbitrary number of input arguments - Caching data chunk import in
definePartAPIfunction in my large data framework
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
- shaving the last 50 ms off NMinimize
- Compiling and memoizing
Outerfor an arbitrary number of input arguments - Implementing a function which generalizes the merging step in merge sort
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
- General discussion of memoization in Mathematica
- How to memoize a function with options
More complex memoization
- Dynamic programming for functions with more than one argument
- Dynamic Programming with delayed evaluation
- selecting minimal subsets
Algorithmic memoization
- selecting minimal subsets
- Happy 2K Prime question (Rojo's solution)
- Count number of sub-lists not greater than a given maximum
- Performance tuning for game solving - peg solitaire
Automation of caching
- Clearing the cache
- Caching automation using deques (my version)
- Caching automation (Szabolcs's version)
Caching and encapsulation of cached values
- Accessing list elements by name
- Caching file import
- Compiling and memoizing
Outerfor an arbitrary number of input arguments - Caching data chunk import in
definePartAPIfunction in my large data framework - make-mathematica-treat-e-i2-as-numeric
- Built-in Mathematica data - are they cached?How to speed up the loading
Memoization of compiled definitions
- JIT Compilation with memoization
- shaving the last 50 ms off NMinimize
- Compiling and memoizing
Outerfor an arbitrary number of input arguments - Implementing a function which generalizes the merging step in merge sort
Automatic clean-up of memoized definitions and the self-blocking trick
- How to automatically localize and clear memoized defintions
Memoization and parallel computations
- Parallelize evaluation of function with memoization
Memoization in probabilistic inference
- Probabilistic Programming with Stochastic Memoization
Linked lists: general discussion
I would like to describe here a pattern which is IMO quite important yet seems to be generally under-appreciated and under-used in our community (and probably in Mathematica programming in general) - namely, linked lists.
Linked lists are important for several reasons. Here is a (partial) list of benefits they can bring to the table:
Idiomatic recursive functional programming. Linked lists allow one to use recursion in Mathematica in ways idiomatic for functional programming. Many people actually consider linked lists and recursion to be one of the corner stones of FP. Since main Mathematica data structures - lists, are actually arrays rather than linked lists, recursion using lists is not a very practical tool in Mathematica, because of performance reasons (we will discuss that in more detail below).
Transparent, easy to guess complexity of algorithms implemented using linked lists. This is quite important particularly in Mathematica, where a lot of internal details may be hidden behind concise programming constructs and built-in functions. In general, Mathematica is a two (or multi) scale language in terms of performance, where one performance scale is roughly given by built-in functions, while another one by the best-written top-level code (the real situation is yet more complex, but this schematic view will be good enough for us). Many algorithms which may be concisely implemented may be fast for smaller samples but then change their behavior, since their speed may be due to small time constants associated with built-in functions, rather than due to the asymptotic behavior of the resulting algorithm itself. And in any case, my point is that an reasonably efficient Mathematica solution may require an expert to quickly determine the complexity of the resulting algorithm.
The algorithms using linked lists typically don't use the "standard" machinery such as vectorization and working with large amounts of data at once, so they are much easier to understand since they use much less of the internal machinery.
Predictable and usually quite reasonable speed. While solutions based on linked lists may be slower than vectorized versions, they are still among the fastest solutions one can achieve purely with the top-level code. And sometimes, they can actually lead to some of the fastest solutions for a problem (one example is a Boggle problem).
Memory efficiency. The solutions using linked lists are typically quite memory efficient, often much more so that vectorized solutions.
Generality. Because solutions based on linked lists don't usually require to have or know all the data at once, they can be used for a more general class of problems than the techniques such as vectorization or structural operations. In particular, there are several classes of problems where linked lists provide by far the most natural and clear solutions.
One final comment here is that, to an extent, declaring linked lists as a design pattern means that here we face the language limitations. Were linked lists so natural to use within the language, and we would not feel the need to describe it as a pattern. To some extent also, it is because the dominant paradigm of practical Mathematica programming is to work with as much data at once as possible and use vectorized and structural operations, an approach which is not quite compatible with linked lists and the way one works with them.
Basic linked lists
Construction
A linked list is a construct of this sort:
llist = {1,{2,{3,{4,{5,{}}}}}}
It can be conveniently constructed from a flat list, using Fold:
toLinkedList[l_List]:=Fold[{#2,#1}&,{},Reverse@l]
so that llist could have been constructed as toLinkedList[Range[5]]. The reverse operation of converting to a flat list can be also done very conveniently using Flatten:
Flatten[llist, Infinity]
(* {1, 2, 3, 4, 5} *)
When elements are lists themselves, one needs to use generalized linked lists, which I will cover later.
The main structural elements of a linked lists are head (the first element), and tail (a list containing the rest of the elements). Note that tail is also a linked list. Extracting one or the other is very easy - just use First for head and Last or tail:
First@llist
(* 1 *)
Last@llist
(* {2,{3,{4,{5,{}}}}} *)
The defining algorithmic property and an illustration - accumulation of results.
In the algorithmic aspect, the main, defining property of a linked list, which makes it very different from usual Mathematica list, is that it is cheap to copy. The reason is that any linked list is a list of two elements only - head and tail. And because expressions in Mathematica are immutable, the copying of tail does not involve a deep copy, but rather (internally) a pointer to the list is copied. This is in contrast to normal Mathematica lists, whose copying will in most cases have a linear complexity in the size of the list.
Everything else follows from this property. For example, collecting of intermediate results is easy with linked lists. The main idiom here is
list = {new-elem, list}
One practical example of using this idiom can be found here. To illustrate this here, consider a toy implementation of Select based in linked lists:
selectLL[l_List, pred_] :=
Module[{result = {}},
Do[If[pred[el], result = {el, result}], {el, l}];
Reverse@Flatten[result]];
Here is how it is used (note that it will only work correctly for list elements which are not themselves lists):
selectLL[Range[10], OddQ]
(* {1, 3, 5, 7, 9} *)
Now, the analog we may write using normal lists will be
selectNaive[l_List, pred_] :=
Module[{result = {}},
Do[If[pred[el], AppendTo[result, el]], {el, l}];
result];
But in this case, entire list is copied every time when a new element is added, which results in much worse performance even for relatively small lists:
selectLL[Range[20000],OddQ];//AbsoluteTiming
selectNaive[Range[20000],OddQ];//AbsoluteTiming
(*
{0.037109,Null}
{0.561523,Null}
*)
Note by the way that the built-in Select is not dramatically faster than our implementation based on linked lists:
Select[Range[20000], OddQ]; // AbsoluteTiming
(* {0.014648, Null} *)
One can argue that Sow and Reap can be used for accumulation of results, and in fact may be more efficient in that than linked lists. This is true, but the Sow - Reap model is more limited because partial lists of results accumulated "so far" are not available to the user.
Generalized linked lists
The general construction for linked lists is similar to the one described above, but the List head is replaced by some (arbitrary) head (say,ll). For certain reasons, it is better in many cases if ll is HoldAllComplete (basically, to prevent search for UpValues, which may be expensive and slow things down). So, we end up with
ClearAll[ll];
SetAttributes[ll,HoldAllComplete];
toLinkedListGen[l_List]:= Fold[ll[#2,#1]&,ll[],Reverse@l]
for example
llgen = toLinkedListGen[Range[5]]
(* ll[1, ll[2, ll[3, ll[4, ll[5, ll[]]]]]] *)
and the code to convert this back to the flat list may look like
List @@ Flatten[llgen, Infinity, ll]
(* {1, 2, 3, 4, 5} *)
This construction is valid for lists with arbitrary elements (except those having themselves the head ll). If you need to use linked lists inside other linked lists, you can define several such heads.
One thing to keep in mind here is that you have to make sure that the element you are adding has been fully evaluated (if that's essential) before writing
llist = ll[elem, llist]
if ll is HoldAllComplete.
Typical use cases and examples
Here I will try to show some examples and give some links for the typical use cases of linked lists, to achieve various goals.
Implementing algorithms requiring complex list iteration, possibly with look-ahead analysis.
A typical such problem is the merge sort algorithm. I will show here the implementation I posted here:
Clear[toLinkedList];
toLinkedList[x_List] := Fold[{#2, #1} &, {}, Reverse[x]];
Module[{h,lrev},
mergeLinked[x_h,y_h]:=
Last[
{x,y,h[]}//.
{
{fst:h[hA_,tA_h],sec:h[hB_,tB_h],e_h}:>
If[hA>hB,{tA,sec,h[hA,e]},{fst,tB,h[hB,e]}],
{fst:h[hA_,tA_h],h[],e_h}:>{tA,h[],h[hA,e]},
{h[],sec:h[hB_,tB_h],e_h}:>{h[],tB,h[hB,e]}
}
];
lrev[set_]:=
Last[h[set,h[]]//.
h[h[hd_,tl_h],acc_h]:>h[tl,h[hd,acc]]];
sort[lst_List]:=
Flatten[
(h[#1,h[]]&)/@lst//.
x_List:>
Flatten[
{toLinkedList[x],{}}//.
{{hd1_,{hd2_,tail_List}},accum_List}:>
{tail,{accum,lrev[mergeLinked[hd1,hd2]]}}
],
\[Infinity],
h
]
]
Here, two different linked lists are used. One is for a list of already merged lists, which at every iteration gets half its previous lengths, since already ordered sub-lists are merged pairwise. For this one, a normal List head is used, and this "higher-level" merging procedure is implemented in sort. Another one is used for individual sub-lists which are merged - they are first converted to linked lists with head h. The lower-level merging procedure is implemented in mergeLinked, while lrev reverses such linked list.
This implementation is a bit unusual in that it uses exclusively linked lists, local rules and ReplaceRepeated, instead of recursion. It was originally posted on MathGroup, where it is explained in detail, and the corrected version I use in this post was posted here. The main idea was to illustrate that ReplaceRepeated can be just as efficient as other means, when used with appropriately constructed patterns. My implementations were certainly motivated by a similar discussion in the book of David Wagner.
In any case, the naive implementation using normal lists, which can look e.g. like
Clear[merge];
merge[x_List,y_List]:=
Block[{merge},
Flatten[
merge[x,y]//.
{
merge[{a_,b___},{c_,d___}]:>
If[a<c,{a,merge[{b},{c,d}]},{c,merge[{a,b},{d}]}],
merge[{},{a__}]:>{a},
merge[{a__},{}]:>{a}
}
]
]
Will have an unacceptable complexity due to the list copying at every stage.
Notes on recursion
Note that while in the previous example recursion was not used explicitly, this was rather artificial, and normally the solutions involving linked lists are recursive. To finish with the merge sort example and illustrate this point, here is the recursive implementation of merge, taken from here
ClearAll[merge];
merge[a_ll, ll[], s_, _] := List @@ Flatten[ll[s, a], Infinity, ll];
merge[ll[], b_ll, s_, _] := List @@ Flatten[ll[s, b], Infinity, ll];
merge[ll[a1_, atail_], b : ll[b1_, _], s_, f_: LessEqual] /;f[a1, b1] :=
merge[atail, b, ll[s, a1], f];
merge[a : ll[a1_, _], ll[b1_, brest_], s_, f_: LessEqual] :=
merge[a, brest, ll[s, b1], f];
merge[a_List, b_List, f_: LessEqual] :=
merge[toLinkedListGen@a, toLinkedListGen@b, ll[], f];
(where I replaced toLinkedList with toLinkedListGen, the latter defined previously in this post).
Four additional important things to note here: first, note how we can effectively combine pattern-matching on linked lists and function overloading (more than one definition / global rule), to end up with code which is quite expressive and self-documenting.
Second, note that the last definition is for normal lists, and automatically translates them to linked lists first, while the first two definitions automatically translate the result into the flat list again. This is a common technique, it allows the function to be self-contained, in that linked lists are only used internally, while both input and output is in terms of normal lists.
Third, you can see that the list of accumulated results (which is also a linked list), is being passed to the function as its third argument. This is also a common /standard technique - make the function pass all additional information that it changes in the process, as extra parameter(s). Not only does this allow one to avoid introducing mutable state into the function (to hold this changing information), but also, and perhaps more importantly, it allows the function to be tail-recursive in the Mathematica sense. This is an absolute requirement for recursive solutions to be practical. The main reason why the function won't be tail-recursive otherwise is explained in the link I just gave, but basically you will have something like
f[x_]:= (
Do-something-with-mutable-state;
f[x-1]
)
And the presence of CompoundExpression around is what makes it non-tail-recursive.
Last thing to note is that, for large enough lists, one may need to increase the setting of $IterationLimit, or set it to infinity altogether. This can always be done on the level of function's definitions. For example, for merge, we could change the last definition to
merge[a_List, b_List, f_: LessEqual] :=
Block[{$IterationLimit = Infinity},
merge[toLinkedListGen@a, toLinkedListGen@b, ll[], f]
];
Note also that if the system complains about $RecursionLimit being exceeded, this means that your function is not properly tail-recursive (as far as Mathematica is concerned). All properly tail-recursive Mathematica functions can only exceed $IterationLimit, but not $RecursionLimit.
General list processing based on some look-ahead analysis
A typical class of use cases for linked lists consists of problems of splitting lists into various groups, according to their structure. I will show a simple example, takesn from here, where the built-in Split would do, but in general Split will not suffice, since the splitting condition may involve more than just two adjacent elements.
This function will count the lengths of sequence of adjacent elements, using linked lists:
Clear[cValsRec];
cValsRec[l_List] := cValsRec[toLinkedList@l, {{}, 0}];
cValsRec[{}, counts_] := Most@Flatten[counts];
cValsRec[{val_, tail : {val_, _List}}, {prev_List, count_}] :=
cValsRec[tail, {prev, count + 1}];
cValsRec[{val_, tail_List}, {prev_List, count_}] :=
cValsRec[tail, {{prev, count + 1}, 0}];
and it can be used e.g. as
test = Sort@RandomInteger[5, 15]
(* {0, 0, 2, 2, 2, 3, 3, 3, 4, 4, 4, 4, 5, 5, 5} *)
cValsRec[test]
(* {2, 3, 3, 4, 3} *)
As before, you can see the same pattern: the first definition is a "starting" one and converts initial list to a linked list. The last definition describes the final stage and conversion the result (which is also accumulated in a linked list) back to normal list. And the results are accumulated in a separate linked list which is passed to the next function call as a separate argument.
Sometimes, instead of splitting the list, one has to do something else, like e.g. deleting or adding some elements, but also based on near-by elements. All these cases are problems for which functional iterators such as Fold are generally ill-suited. I also add such cases to the set of cases I cover in this section, since they are conceptually similar.
Some further links for similar problems solved using linked lists:
- Splitting a list by specifying section headers
- Splitting words into specific fragments
- Performance collapsing repeated contiguous rows cols from a matrix (note that in this example the condition is more complex. The code there illustrates well the powerful combination of linked lists and pattern-matching, in a less trivial setting).
- Split a string at specific positions
- Simpler way to fill date gaps with zero values
- Leveling peaks in list
- Searching from the end of a list in Mathematica
- Splitting up delimited data in lists
- Implementing a function which generalizes the merging step in merge sort
- rule-based implementation of an algorithm
- Partitioning a list by recursion (Programming Paradigms via Mathematica (A First Course))
I particularly recommend to look into the code of the last two links, since those use a number of interesting techniques in combination.
Some general remarks for this section: first, I'd like to stress one clear advantage of linked lists for this type of problems: they are completely straight-forward to use. You don't have to figure out some clever vectorizing scheme or sequence of list-manipulation steps. You don't have to know all subtleties of various built-in functions involved, and their combinations. This means that linked lists provide a straight-forward implementation path for users who are not Mathematica experts.
Second, the code using linked lists and recursion tends to be readable and self-documenting. This is quite important, since cryptic code is hard to understand and maintain, and it is easier to create bugs when the code is cryptic.
Third, code based on linked lists can be completely generic, since it does not rely on vectorization, packed arrays etc. At the same time, it remains completely idiomatic, since linked lists are immutable, and no mutable state is used. Admittedly, for numerical problems where the structure of the data may allow (much) more efficient vectorized solutions, linked lists may not be the fastest. But in those cases, one can dispatch to such special cases automatically, as it has been done e.g. for the merge function here.
Using linked lists for tree traversals and tree-building
This is a more advanced use of linked lists, and here they can fully reveal their potential. The thing is, problems involving trees are usually naturally recursive, and while for flat list creation / processing it is very often possible to get away by using standard Mathematica arsenal for list manipulation, here itis often not enough.
While I plan to extend this part of the post soon, right now I will illustrate this with one relatively simple example - a solution for a variant of the cutting stock problem in Mathematica, and then provide links to a few more involved cases.
So, the problem there was to list the ways in which one could sum elements from a list (with repetitions), so the sum is less-than-or-equal to some value. So, here are the values of the list and the maximum value, from the original question:
i = {7.25, 7.75, 15, 19, 22};
m = 22;
The solution using linked lists is fairly simple:
ClearAll[getCombs];
getCombs[l_List,lim_Integer]:=
Cases[
getCombs[{},toLinkedList[Reverse[Sort[l]]],lim],
getCombs[acc:Except[{}],__]:>Flatten[acc],
\[Infinity]
];
getCombs[accum_,l:{head_,tail_List},lim_?Positive]:=
{
If[lim>=head,getCombs[{accum,head},l,lim-head],Sequence@@{}],
getCombs[accum,tail,lim]
};
The first definition is both the "starting" and "ending" one, and it realizes the "wrap-around definition" pattern - that is, preprocess data (in this case, convert to linked list the original list), call the main function, and then post-process the data (Cases here, I will explain what it does in a moment).
The second definition is the one doing the real work. Borrowing the explanation from the original answer, the way it works is to create a tree by going backwards, subtracting a given element from the current maximal for the sum, to get the maximal for the sums of remaining elements. The case of repetitions is accounted for by the line If[lim>=head,...]. The resulting numbers are also accumulated in linked lists (accum parameter), and finally they are found with Cases and converted back to normal lists using Flatten. In this way, we avoid unnecessary array copying.
As you can see, the most important difference here is that a tree is created as a result, rather than a list. That tree is then post-processed, and the results collected from it, using Cases. This construction allows us to take the best from both worlds: linked lists + recursion, and pattern-based expression destructuring. And since the patterns are "syntactic" See my comment to a post here, it is really fast. Because the algorithm is essentially recursive, and the result is a tree rather than a more regular list, the standard list-manipulation approaches are at a loss here.
Here are some more links to applications of linked lists involving trees:
- Can a trie be implemented efficiently
- How can I use Mathematicas graph functions to cheat at boggle/
- Performance tuning for game solving: peg solitaire
In some of the examples linked, the use of linked lists or associated data structures (and trie can be viewed as a generalization of a linked list) led to some quite dramatic performance boosts. You can read explanations in the linked posts for more details.
Remarks
I tried to outline the usefulness of linked lists and illustrate the main ideas behind their use, and techniques associated with them. One of the main advantages they offer is the possibility of direct and straight-forward implementation of recursive algorithms in Mathematica, with quite reasonable performance.
This post is a work in progress. I plan to add more links, and make a special separate link section where all links will be collected. I also may add a couple of more sections, such as a section of using linked lists to trade speed for memory.
All remarks, suggestions and additions are welcome and encouraged - this is a community post as are all answers here.
"Design patterns" are a standard approach to a class of problems. They may become idioms if the community adopts them. They may be simple, but non-trivial to figure out and identify (for oneself) as a pattern worthy of re-use.
Here's one particular pattern I personally like to use:
Example: the "TransformBy" approach
There are functions which transform lists based on the properties of list elements, e.g. Sort, Split, Gather. Often it's advantageous to do the transformation based on a transform of the list elements, so we have SortBy, SplitBy, GatherBy. It's particularly useful when we deal with equivalence classes which are easier to find based on the comparison f[x] == f[y] than a completely custom comparison function, or for ordering-related problems (put < in place of ==).
This approach becomes a "design pattern" once I internalise it enough that I start approaching problems in terms of it.
There's no built-in MaxBy or MinBy (probably because the result of these wouldn't necessarily be unique), but they're often very useful, so I define them
MaxBy[list_, fun_] := list[[First@Ordering[fun /@ list, -1]]]
MinBy[list_, fun_] := list[[First@Ordering[fun /@ list, 1]]]
Then we can use them in so many ways. I commonly do things like MaxBy[dataTable, Last], MinBy[dataList, Variance], MaxBy[lists, Length], MinBy[points, Norm], etc.
What makes this a "design pattern" is that I internalised the approach enough that it comes as naturally to me as using built-ins such as Outer. Or, to put it less subjectively, that the pattern turned out to be general enough that it was possible and useful for me to internalise it.
Other candidates for a ...By pattern (with caveats): Union, Intersection, Tally, DeleteDuplicates etc. All these deal with equivalence classes and some of them speed up operations using sorting. The ...By version would typically be faster than giving a custom comparison function (which precludes using sorting or binary search for speedup).
Some other general patterns:
Memoization in the form of
f[k_] := f[k] = ...The Villegas-Gayley trick
Vectorization: instead of processing a list (or lists) element wise, use operations that work with whole lists and keep them packed. I'm not sure if I'd call this a design pattern or a paradigm, but in Mathematica, the pattern of vectorizing things is rather standard: Use arithmetic operations on lists; use
UnitStepfor comparisons (>,>=,<,<=); use1and0forTrue/Falseand add/multiply them forandandoroperations; usePickfor conditional extraction. Here's a (very simple) example. Doing this usually leads to a significant speedup, but it also creates code that is not easy and quick to modify. Some other languages, such as MATLAB, provide much better support for this style of programming, but it's possible in Mathematica too (with the same performance benefits)....Blockfor resource management, analogous toJavaBlockorNETBlock. This is a bit specialized an advanced, and it's useful for package authors. It's a very useful pattern for simplifying resource management when one needs to deal with objects that have to be created and destroyed explicitly. (I extendedTetGenLinkto do this as well by wrappingTetGenCreateusing Villegas-Gayley to keep track of created objects.) See also here.One could argue that the standard package structure is a design pattern as well:
BeginPackage[...], mention public symbols (coupling this conveniently with setting usage messages),Begin["`Private`"], definitions,End[]; EndPackage[]. Theoretically it could be done differently. There's nothing forcing me to write packages using this pattern. But this is a time tested pattern for a package structure that makes things easy to manage and robust.Creating closures (i.e. functions with a persistent state) using
Module, as shown here.(will update as I come across them)