What is in your Mathematica tool bag? [closed]

Question

We all know that Mathematica is great, but it also often lacks critical functionality. What kind of external packages / tools / resources do you use with Mathematica?

I'll edit (and invite anyone else to do so too) this main post to include resources which are focused on general applicability in scientific research and which as many people as possible will find useful. Feel free to contribute anything, even small code snippets (as I did below for a timing routine).

Also, undocumented and useful features in Mathematica 7 and beyond you found yourself, or dug up from some paper/site are most welcome.

Please include a short description or comment on why something is great or what utility it provides. If you link to books on Amazon with affiliate links please mention it, e.g., by putting your name after the link.

---

Packages:

1. LevelScheme is a package that greatly expands Mathematica's capability to produce good looking plots. I use it if not for anything else then for the much, much improved control over frame/axes ticks. Its newest version is called SciDraw, and it will be released sometime this year.
2. David Park's Presentation Package (US$50 - no charge for updates)
3. Jeremy Michelson's grassmannOps package provides resources for doing algebra and calculus with Grassmann variables and operators that have non trivial commutation relations.
4. John Brown's GrassmannAlgebra package and book for working with Grassmann and Clifford algebras.
5. RISC (Research Institute for Symbolic Computation) has a variety of packages for Mathematica (and other languages) available for download. In particular, there is Theorema for automated theorem proving, and the multitude of packages for symbolic summation, difference equations, etc. at the Algorithmic Combinatorics group's software page.

Tools:

1. MASH is Daniel Reeves's excellent Perl script essentially providing scripting support for Mathematica v7. (Now built in as of Mathematica 8 with the -script option.)
2. An alternate Mathematica shell with a GNU readline input (using python, *nix only)
3. ColourMaths package allows you to visually select parts of an expression and manipulate them. http://www.dbaileyconsultancy.co.uk/colour_maths/colour_maths.html

Resources:

1. Wolfram's own repository MathSource has a lot of useful if narrow notebooks for various applications. Also check out the other sections such as

   • Current Documentation,
   • Courseware for lectures,
   • and Demos for, well, demos.

2. The Mathematica Wikibook.

Books:

1. Mathematica programming: an advanced introduction by Leonid Shifrin (web, pdf) is a must read if you want to do anything more than For loops in Mathematica. We have the pleasure of having Leonid himself answering questions here.
2. Quantum Methods with Mathematica by James F. Feagin (amazon)
3. The Mathematica Book by Stephen Wolfram (amazon) (web)
4. Schaum's Outline (amazon)
5. Mathematica in Action by Stan Wagon (amazon) - 600 pages of neat examples and goes up to Mathematica version 7. Visualization techniques are especially good, you can see some of them on the author's Demonstrations Page.
6. Mathematica Programming Fundamentals by Richard Gaylord (pdf) - A good concise introduction to most of what you need to know about Mathematica programming.
7. Mathematica Cookbook by Sal Mangano published by O'Reilly 2010 832 pages. - Written in the well known O'Reilly Cookbook style: Problem - Solution. For intermediates.
8. Differential Equations with Mathematica, 3rd Ed. Elsevier 2004 Amsterdam by Martha L. Abell, James P. Braselton - 893 pages For beginners, learn solving DEs and Mathematica at the same time.

Undocumented (or scarcely documented) features:

7. How to customize Mathematica keyboard shortcuts. See this question.
8. How to inspect patterns and functions used by Mathematica's own functions. See this answer
9. How to achieve consistent size for GraphPlots in Mathematica? See this question.
10. How to produce documents and presentations with Mathematica. See this question.

Answer 1

By popular demand, the code to generate the top-10 SO answerers plot (except annotations) using the SO API.

getRepChanges[userID_Integer] :=
 Module[{totalChanges},
  totalChanges = 
   "total" /. 
    Import["http://api.stackoverflow.com/1.1/users/" <> 
      ToString[userID] <> "/reputation?fromdate=0&pagesize=10&page=1",
      "JSON"];
  Join @@ Table[
    "rep_changes" /. 
     Import["http://api.stackoverflow.com/1.1/users/" <> 
       ToString[userID] <> 
       "/reputation?fromdate=0&pagesize=10&page=" <> ToString[page], 
      "JSON"],
    {page, 1, Ceiling[totalChanges/10]}
    ]
  ]

topAnswerers = ({"display_name", 
      "user_id"} /. #) & /@ ("user" /. ("top_users" /. 
      Import["http://api.stackoverflow.com/1.1/tags/mathematica/top-\
answerers/all-time", "JSON"]))

repChangesTopUsers =
  Monitor[Table[
    repChange = 
     ReleaseHold[(Hold[{DateList[
              "on_date" + AbsoluteTime["January 1, 1970"]], 
             "positive_rep" - "negative_rep"}] /. #) & /@ 
        getRepChanges[userID]] // Sort;
    accRepChange = {repChange[[All, 1]], 
       Accumulate[repChange[[All, 2]]]}\[Transpose],
    {userID, topAnswerers[[All, 2]]}
    ], userID];

pl = DateListLogPlot[
  Tooltip @@@ 
   Take[({repChangesTopUsers, topAnswerers[[All, 1]]}\[Transpose]), 
    10], Joined -> True, Mesh -> None, ImageSize -> 1000, 
  PlotRange -> {All, {10, All}}, 
  BaseStyle -> {FontFamily -> "Arial-Bold", FontSize -> 16}, 
  DateTicksFormat -> {"MonthNameShort", " ", "Year"}, 
  GridLines -> {True, None}, 
  FrameLabel -> (Style[#, FontSize -> 18] & /@ {"Date", "Reputation", 
      "Top-10 answerers", ""})]

Answer 2

Recursive pure functions (#0) seem to be one of the darker corners of the language. Here are a couple of non-trivial examples of their use , where this is really useful (not that they can not be done without it). The following is a pretty concise and reasonably fast function to find connected components in a graph, given a list of edges specified as pairs of vertices:

ClearAll[setNew, componentsBFLS];
setNew[x_, x_] := Null;
setNew[lhs_, rhs_]:=lhs:=Function[Null, (#1 := #0[##]); #2, HoldFirst][lhs, rhs];

componentsBFLS[lst_List] := Module[{f}, setNew @@@ Map[f, lst, {2}];
   GatherBy[Tally[Flatten@lst][[All, 1]], f]];

What happens here is that we first map a dummy symbol on each of the vertex numbers, and then set up a way that, given a pair of vertices {f[5],f[10]}, say, then f[5] would evaluate to f[10]. The recursive pure function is used as a path compressor (to set up memoization in such a way that instead of long chains like f[1]=f[3],f[3]=f[4],f[4]=f[2], ..., memoized values get corrected whenever a new "root" of the component is discovered. This gives a significant speed-up. Because we use assignment, we need it to be HoldAll, which makes this construct even more obscure and more attractive ). This function is a result of on and off-line Mathgroup discussion involving Fred Simons, Szabolcs Horvat, DrMajorBob and yours truly. Example:

In[13]:= largeTest=RandomInteger[{1,80000},{40000,2}];

In[14]:= componentsBFLS[largeTest]//Short//Timing
Out[14]= {0.828,{{33686,62711,64315,11760,35384,45604,10212,52552,63986,  
     <<8>>,40962,7294,63002,38018,46533,26503,43515,73143,5932},<<10522>>}}

It is certainly much slower than a built-in, but for the size of code, quite fast still IMO.

Another example: here is a recursive realization of Select, based on linked lists and recursive pure functions:

selLLNaive[x_List, test_] :=
  Flatten[If[TrueQ[test[#1]],
     {#1, If[#2 === {}, {}, #0 @@ #2]},
     If[#2 === {}, {}, #0 @@ #2]] & @@ Fold[{#2, #1} &, {}, Reverse[x]]];

For example,

In[5]:= Block[
         {$RecursionLimit= Infinity},
         selLLNaive[Range[3000],EvenQ]]//Short//Timing

Out[5]= {0.047,{2,4,6,8,10,12,14,16,18,20,22,24,26,28,30,
 <<1470>>,2972,2974,2976,2978,2980,2982,2984,2986,2988,2990,
  2992,2994,2996,2998,3000}}

It is however not properly tail recursive, and will blow the stack (crash the kernel) for larger lists. Here is the tail-recursive version:

selLLTailRec[x_List, test_] :=
Flatten[
 If[Last[#1] === {},
  If[TrueQ[test[First[#1]]],
   {#2, First[#1]}, #2],
  (* else *)
  #0[Last[#1],
   If[TrueQ[test[First[#1]]], {#2, First[#1]}, #2]
   ]] &[Fold[{#2, #1} &, {}, Reverse[x]], {}]];

For example,

In[6]:= Block[{$IterationLimit= Infinity},
       selLLTailRec[Range[500000],EvenQ]]//Short//Timing
Out[6]= {2.39,{2,4,6,8,10,12,14,16,18,20,22,
       <<249978>>,499980,499982,499984,499986,499988,499990,499992,
        499994,499996,499998,500000}} 

Answer 3

I find it really useful when developing packages to add this keyboard shortcut to my SystemFiles/FrontEnd/TextResources/Windows/KeyEventTranslations.tr file.

(* Evaluate Initialization Cells: Real useful for reloading library changes. *)

Item[KeyEvent["i", Modifiers -> {Control, Command}],
    FrontEndExecute[
        FrontEndToken[
            SelectedNotebook[],
            "EvaluateInitialization"]]],

Next for every Packagename.m I make a PackagenameTest.nb notebook for testing and the first 2 cells of the test notebook are set as initialization cells. In the first cell I put

Needs["PackageManipulations`"]

to load the very useful PackageManipulations library which was written by Leonid. The second cell contains

PackageRemove["Packagename`Private`"]
PackageRemove["Packagename`"]
PackageReload["Packagename`"]

which all do the actual package reloading. Note the first two lines are there only to Remove all symbols as I like to keep the contexts as clean as possible.

Then the workflow for writing and testing a package becomes something like this.

1. Save changes to Packagename.m.
2. Go to PackagenameTest.nb and do CTRL + ALT + i.

This causes the initialization cells to reload the package, which makes testing real simple.

Answer 4

One trick I've used, which allows you to emulate the way most built-in functions work with bad arguments (by sending a message and then returning the whole form unevaluated) exploits a quirk of the way Condition works when used in a defintion. If foo should only work with one argument:

foo[x_] := x + 1;
expr : foo[___] /; (Message[foo::argx, foo, Length@Unevaluated[expr], 1]; 
                    False) := Null; (* never reached *)

If you have more complex needs, it's easy to factor out the argument validation and message generation as an independent function. You can do more elaborate things by using side effects in Condition beyond just generating messages, but in my opinion most of them fall into the "sleazy hack" category and should be avoided if possible.

Also, in the "metaprogramming" category, if you have a Mathematica package (.m) file, you can use the "HeldExpressions" element to get all the expressions in the file wrapped in HoldComplete. This makes tracking things down much easier than using text-based searches. Unfortunately, there's no easy way to do the same thing with a notebook, but you can get all the input expressions using something like the following:

inputExpressionsFromNotebookFile[nb_String] :=
 Cases[Get[nb],
  Cell[BoxData[boxes_], "Input", ___] :>
   MakeExpression[StripBoxes[boxes], StandardForm],
  Infinity]

Lastly, you can use the fact that Module emulates lexical closures to create the equivalent of reference types. Here's a simple stack (which uses a variation the Condition trick for error handling as a bonus):

ClearAll[MakeStack, StackInstance, EmptyQ, Pop, Push, Peek]
 With[{emptyStack = Unique["empty"]},
  Attributes[StackInstance] = HoldFirst;
  MakeStack[] :=
   Module[{backing = emptyStack},
    StackInstance[backing]];

  StackInstance::empty = "stack is empty";

  EmptyQ[StackInstance[backing_]] := (backing === emptyStack);

  HoldPattern[
    Pop[instance : StackInstance[backing_]]] /;
    ! EmptyQ[instance] || (Message[StackInstance::empty]; False) :=
   (backing = Last@backing; instance);

  HoldPattern[Push[instance : StackInstance[backing_], new_]] :=
   (backing = {new, backing}; instance);

  HoldPattern[Peek[instance : StackInstance[backing_]]] /;
    ! EmptyQ[instance] || (Message[StackInstance::empty]; False) :=
   First@backing]

Now you can print the elements of a list in reverse order in a needlessly convoluted way!

With[{stack = MakeStack[], list},
 Do[Push[stack, elt], {elt, list}];

 While[!EmptyQ[stack],
  Print[Peek@stack];
  Pop@stack]]

Answer 5

Internal`InheritedBlock

I have learned recently the existence of such useful function as Internal`InheritedBlock, from this message of Daniel Lichtblau in the official newsgroup.

As I understand, Internal`InheritedBlock allows to pass a copy of an outbound function inside the Block scope:

In[1]:= Internal`InheritedBlock[{Message},
Print[Attributes[Message]];
Unprotect[Message];
Message[x___]:=Print[{{x},Stack[]}];
Sin[1,1]
]
Sin[1,1]
During evaluation of In[1]:= {HoldFirst,Protected}
During evaluation of In[1]:= {{Sin::argx,Sin,2},{Internal`InheritedBlock,CompoundExpression,Sin,Print,List}}
Out[1]= Sin[1,1]
During evaluation of In[1]:= Sin::argx: Sin called with 2 arguments; 1 argument is expected. >>
Out[2]= Sin[1,1]

I think this function can be very useful for everyone who need to modify built-in functions temporarily!

Comparison with Block

Let us define some function:

a := Print[b]

Now we wish to pass a copy of this function into the Block scope. The naive trial does not give what we want:

In[2]:= Block[{a = a}, OwnValues[a]]

During evaluation of In[9]:= b

Out[2]= {HoldPattern[a] :> Null}

Now trying to use delayed definition in the first argument of Block (it is an undocumented feature too):

In[3]:= Block[{a := a}, OwnValues[a]]
Block[{a := a}, a]

Out[3]= {HoldPattern[a] :> a}

During evaluation of In[3]:= b

We see that in this case a works but we have not got a copy of the original a inside of the Block scope.

Now let us try Internal`InheritedBlock:

In[5]:= Internal`InheritedBlock[{a}, OwnValues[a]]

Out[5]= {HoldPattern[a] :> Print[b]}

We have got a copy of the original definition for a inside of the Block scope and we may modify it in the way we want without affecting the global definition for a!

Answer 6

Mathematica is a sharp tool, but it can cut you with its somewhat untyped behaviour and avalanches of cryptic diagnostic messages. One way to deal with this is to define functions following this idiom:

ClearAll@zot
SetAttributes[zot, ...]
zot[a_] := ...
zot[b_ /; ...] := ...
zot[___] := (Message[zot::invalidArguments]; Abort[])

That is a lot of boilerplate, which I'm frequently tempted to skip. Especially when prototyping, which happens a lot in Mathematica. So, I use a macro called define that allows me to stay disciplined, with much less boilerplate.

A basic usage of define is like this:

define[
  fact[0] = 1
; fact[n_ /; n > 0] := n * fact[n-1]
]

fact[5]

120

It doesn't look like much at first, but there are some hidden benefits. The first service that define provides is that it automatically applies ClearAll to the symbol being defined. This ensures that there are no leftover definitions -- a common occurrence during the initial development of a function.

The second service is that the function being defined is automatically "closed". By this I mean that the function will issue a message and abort if it is invoked with an argument list that is not matched by one of the definitions:

fact[-1]

define::badargs: There is no definition for 'fact' applicable to fact[-1].
$Aborted

This is the primary value of define, which catches a very common class of error.

Another convenience is a concise way to specify attributes on the function being defined. Let's make the function Listable:

define[
  fact[0] = 1
; fact[n_ /; n > 0] := n * fact[n-1]
, Listable
]

fact[{3, 5, 8}]

{6, 120, 40320}

In addition to all of the normal attributes, define accepts an additional attribute called Open. This prevents define from adding the catch-all error definition to the function:

define[
  successor[x_ /; x > 0] := x + 1
, Open
]

successor /@ {1, "hi"}

{2, successor["hi"]}

Multiple attributes may be defined for a function:

define[
  flatHold[x___] := Hold[x]
, {Flat, HoldAll}
]

flatHold[flatHold[1+1, flatHold[2+3]], 4+5]

Hold[1 + 1, 2 + 3, 4 + 5]

Without further ado, here is the definition of define:

ClearAll@define
SetAttributes[define, HoldAll]
define[body_, attribute_Symbol] := define[body, {attribute}]
define[body:(_Set|_SetDelayed), attributes_List:{}] := define[CompoundExpression[body], attributes]
define[body:CompoundExpression[((Set|SetDelayed)[name_Symbol[___], _])..], attributes_List:{}] :=
  ( ClearAll@name
  ; SetAttributes[name, DeleteCases[attributes, Open]]
  ; If[!MemberQ[attributes, Open]
    , def:name[___] := (Message[define::badargs, name, Defer@def]; Abort[])
    ]
  ; body
  ;
  )
def:define[___] := (Message[define::malformed, Defer@def]; Abort[])

define::badargs = "There is no definition for '``' applicable to ``.";
define::malformed = "Malformed definition: ``";

The exhibited implementation supports neither up-values nor currying, nor patterns more general than simple function definition. It remains useful, however.