*Original Author: Angelo Pesce / C0de517e*

This will appear in a longer, more rambling form, on C0de517e.

It’s main section is inspired by LearnXinYminutes

## Introduction

Games and rendering are becoming increasingly “data driven”, in a sense that we use big amounts of “offline” data, either through acquisition or simulation, that we have to fit into our games via forms of approximation. These can be derived analytically, through simplifying assumptions, numerically, by fitting functions to data (learning), or using both methods.

Mathematica can help with all these steps:

- Being a Computer Algebra System, it aids with experimentation on analytic forms
- It has great support for a wide range of numerical solvers, allowing to write simulations
- It has strong visualization and interactive manipulation abilities, which makes reasoning about the data easier
- It’s a prototyping, experimention language, with an interactive REPL and very coincise, powerful functional primitives

## Mathematica’s language

It stands to reason that a CAS is built on top of a symbolic language that supports programmatic manipulation of its own programs (code as data, a.k.a. homoiconicity), and indeed this is the case here. The most famous homoiconic language is Lisp, and indeed you’re familiar with the Lisp family of languages, Mathematica won’t be too far off, but there are a few notable differences.

While in Lisp everything is a list, in Mathematica, everything is an expression tree. Also, expressions in mathematica can have different forms, that is, input (or display) versions of the same internal expression node. This allows you for example to have equations entered in the standard mathematical notation (TraditionalForm) via equation editors or in a textual form that can be typed without auxiliary graphical symbols (InputForm) and so on. Mathematica’s environment, the notebook, is not a purely textual one, but supports graphics, so even images and graphs can be displayed as output or input, inside equations, while still maintaining the same internal representation.

Mathematica is an interactive environment, but it’s not a standard REPL (read-eval-print loop), instead it relies on the concept of “notebooks” which are a collection of “cells”. Each cell can be evaluated (shift-enter) and it will yield an output cell underneath them, thus allowing for changes and re-evaluation of cells in any order. Cells can also be marked as not containing Mathematica code but just text, thus the notebook is a mix of code and documentation which enables a sort of “literary programming” style.

For completeness it’s worth noticing that Mathematica also has a traditional text-only interface that can be invoked by running the Kernel outside the notebook environment, which has only textual input and output and has only the standard REPL you would expect, but there’s little reason to use it. There is also a more “programming” oriented environment called the Workbench, an optional product that can make your life easier if you write lots of Mathematica code and need to profile, debug and so on.

## Mathematica by example

If you don’t have Mathematica, you can try Mathics which is an OpenSource implementation based on SciPy and Sage. It’s obviously far from being as complete as Wolfram’s solution, but it has a lot of the basic language covered and also has an online interface! I expect most of the code below to work on it.

(* This is a comment, if you're entering this in a notebook remember that to evaluate the content of a cell you need to use shift-enter or the numeric pad enter *) (* Basic math is as expected, but it's kept at arbitrary precision unless you use machine numbers *) (1+2)*3/4 (1.+2.)*3./4. (* % refers to the last computed value *) %+2 (* Functions are invoked passing parameters in square braces, all built-in functions start with capitals*) Sin[Pi/3] (* N[] forces evaluation to machine numbers, using machine numbers makes evaluation faster, but will defeat many CAS functions *) N[Sin[Pi/3]] (* Infix and postfix operators all have a functional form, use FullForm to show *) FullForm[Hold[(1+2)*3/4]] (* Each expression in a cell will yield an output in a separate output cell. Expressions can be terminated with ; if we don't want them to emit output, which is useful when doing intermediate assignments that would yield large outputs otherwise *) 1+2; (* Assigning a symbol to an expression *) x = 10 (* If a symbol is not yet defined, it will be kept in its symbolic version as the evaluation can't proceed further *) y = x*w (* This will recursively expand z until it reaches expansion limit and errors out *) z = z+1 (* Clears the previous assignments. It's not wise to assign as globals such common symbols, we use these here for brevity and will clear as needed *) Clear[x,y,z] (* Evaluation is controlled by symbols attributes *) x = 10 (* y will be equal to "x*2", not 20 as := is the infix version of the function SetDelayed, which doesn't evaluate the right hand...*) y := x*2 (* …that's because SetDelayed has attribute HoldAll, which tells the evaluator to not evaluate any of its arguments. HoldAll and HoldFirst attributes are one of the "tricky" parts, and a big difference from Lisp where you should explicitly quote to stop evaluation *) Attributes[SetDelayed] (* As many functions in Mathematica are supposed to deal with symbolic expressions and not their evaluated version, you'll find that many of them have HoldAll or HoldFirst, for example Plot has HoldFirst to not evaluate its first argument, that is the expression that we want to graph *) Plot[Sin[x], {x, 0, 6*Pi}] (* The Hold function can be used to stop evaluation, and the Evaluate function can be used to counter-act HoldFirst or HoldAll *) Hold[x*2] y:=Evaluate[x*2] y (* A neat usage of SetDelayed is for memoization of computations, the following pi2, the first time it will be evaluated, will set itself to the numerical value of Pi*Pi to 50 decimal points *) pi2:=pi2=N[Pi*Pi,50] pi2 (* Defining functions can be done with the Function function, which has attributes HoldAll *) fn=Function[{x,y}, x*y]; fn[10,20] (* As many Mathematica built-ins, Function has multiple input forms, the following is a shorthand with unnamed parameters #1 and #2, ended with the & postfix *) fn2=#1*#2& fn2[10,20] (* Third version, infix notation. Note that [Function] is a textual representation of a graphical symbol that can be more easily entered in Mathematica with the key presses: esc f n esc, many symbols can be similarly entered, try for example esc theta esc *) fn3={x,y}[Function]x*y fn3[10,20] (* A second, very common way of defining functions is to use pattern matching and delayed evaluation, the following defines the fn4 symbol to evaluate the expression x*y when it's encountered with two arguments that will matched to the symbols x and y *) fn4[x_,y_]:=x*y fn4[10,20] (* _ or Blank[] can match any Mathematica expression, _h matches only expressions with the Head[] h *) fn5[x_Integer,y_Integer]:=x+y fn5[10,20] fn5[10,20.] (* A symbol can have multiple matching rules *) fn6[0] = 1; fn6[x_Integer] := x*fn6[x - 1] fn6[3] (* In general pattern matching is more powerful than Function as it's really an evaluation rule, but it's slower to evaluate, thus not the best if a function has to be applied over large datasets *) (*Note that pattern matching can be used also with=,not only:=,but beware that=evaluates RHS,in the following fnWrong will multiply y by 3,not by the value matching test at "call" site,as test*y gets fully evaluated and test doesn't "stay" a symbol,it evaluates to its global value*) test=3; fnWrong[test_,y_]=test*y (* Lists are defined with {} *) a={1,2,3,{4,5},{aa,bb}}

(* Elements are accessed with [[index]], indices are one-based, negative wrap-around *) a[[1]] a[[-1]] (* Ranges are expressed with ;; or Span *) a[[2;;4]] (* From the beginning to the second last *) a[[;;-2]] (* Vectors and matrices are just appropriately sized lists and lists of lists *) b={1,2,3} m={{1,0,0},{0,1,0},{0,0,1}} (* . is the product for vector, matrices, and tensors *) m.b (* Expression manipulation and CAS. ReplaceAll or /. applies rules to an expression *) (x+y)/.{x->2,y->Sin[Pi]} (* Rules can contain patterns, the following will match only the x symbols that appear to a power, match the expression of the power and replace it *) Clear[x]; 1+x+x^2 +x^(t+n)/.{x^p_->f[p]} (* Mathematica has lots of functions that deal with expressions, Integrate, Limit, D, Series, Minimize, Reduce, Refine, Factor, Expand and so on. We'll show only some basic examples. Solve finds solution to systems of equations or inequalities *) Clear[a]; Solve[x^2+a*x+1==0, x] (* It returns results as list of replacement rules that we can replace into the original equation *) eq=x^2+a*x+1 sol=Solve[eq==0, x] neweq=eq/.sol[[1]] (* Simplifying neweq yields true as the equation is satisfied *) Simplify[neweq] (* Assumptions on the variables can be made *) Simplify[Sqrt[x^2], Assumptions -> x < 0] (* fn7 will compute the Integral and Derivative every time it's evaluated, as Function is HoldAll, fn8, using Evaluate, will force the definition to be equal to the simplified version which yields correctly back the original equation *) fn7[x_]:=Function[x,D[Integrate[x^3,x],x]] fn8[x_]:=Function[x,Evaluate[Simplify[D[Integrate[x^3,x],x]]]] (* Many procedural programming primitives are supported *) If[3>2,10,20] For[i = 0,i < 4,i++,Print[i]] n=1; While[n < 4,Print[n];n++] Do[Print[n^2],{n,4}] (* Boolean operators are C-like for the most, only Xor is not ^ which means Power instead *) !((1>2)||(4>3))&&((1==1)&&(5<=6)) (* Equality tests can be chained *) (5>4>3)&&(1!=2!=3) (* == compares the result of the evaluation on both sides, === is true only if the expression are identical *) v1=1;v2=1; v1==v2 v1===v2 (* Boolean values are False and True. No output is Null *) (* With, Block and Module can be used to set symbols to temporary values in an expression *) With[{x = Sin[y]}, x*y] Block[{x = Sin[y]}, x*y] Module[{x = Sin[y]}, x*y] (* The difference is subtle. With acts as a replacement rule. Block temporarily assigns the value to a symbol and the restores the previous definition. Module creates an unique, temporary symbol, which affects only the occurrences in the inner scope. *) m=i^2 Block[{i = a}, i + m] Module[{i = a}, i + m] (* In general prefer Block or With, which are faster than Module. Module implements lexical scoping, Block does dynamic scoping *) (* Data operations. Table generates data from expressions *) Table[i^2,{i,1,10}] (* Table can generate multi-dimensional arrays, i.e. matrices *) Table[10*i+j,{i,1,4},{j,1,3}] MatrixForm[%] (* List elements can be manipulated using functional programming primitives, like Map which applies a function over a list *) squareListElements[list_]:=Map[#^2&,list] (* Short-hand, infix notation of Map[] is /@ *) squareListElements2[list_]:=(#^2&)/@list (* You can use MapIndexed to operate in parallel across two lists, it passes to the mapped function *) addLists[list1_,list2_]:=MapIndexed[Function[{element,indexList},element + list2[[indexList[[1]]]] ], list1] addLists[{1,2,3},{3,4,5}] (* A more complete version of the above that is defined only on lists and asserts if the two lists are not equal size. Note the usage of ; to compound two expressions and the need of parenthesis *) On[Assert] addListsAssert[list1_List,list2_List]:=(Assert[Length[list1]==Length[list2]]; MapIndexed[Function[{element,indexList},element + list2[[indexList[[1]]]] ], list1]) (* Or Thread can be used, which "zips" two or more lists together *) addLists2[list1_,list2_]:=MapThread[#1+#2&,{list1,list2}] (* There are many functional list manipulation primitives, in general, using these is faster than trying to use procedural style programming. Extract from a list of the first 100 integers, the ones divisible by five *) Select[Range[100],Mod[#,5]==0&] (* Group together all integers from 1...100 in the same equivalence class modulo 5 *) Gather[Range[100],Mod[#1,5]==Mod[#2,5]&] (* Fold repeatedly applies a function to each element of a list and the result of the previous fold *) myTotal[list_]:=Fold[#1+#2&,0,list] (* Another way of redefining Total is to use Apply, which calls a function with as arguments, the elements of a list. The infix shorthand of Apply is @@ *) myTotal2[list_]:=Apply[Plus,list] (* Mathematica's CAS abilities also help with numerical algorithms, as Mathematica is able to infer some information from the equations passed in order to select or optimize the numerical methods *) (* NMinimize does constrained and unconstrained minimization, linear and nonlinear, selecting among different algorithms as needed *) Clear[x,y] NMinimize[{x^2-(y-1)^2, x^2+y^2<=4}, {x,y}] (* NIntegrate does numerical definite integrals. Uses Monte Carlo methods for many-dimensional integrands *) NIntegrate[Sin[Sin[x]], {x,0,2}] (* NSum approximates discrete summations, even to infinites *) NSum[(-5)^i/i!,{i,0,Infinity}] (* Many other analytic operators have numerical counterparts, like NLimit, ND and so on... *) NLimit[Sin[x]/x,x->0] ND[Exp[x],x,1] (* Mathematica's plots produce Graphics and Graphics3D outputs, which the notebook shows in a graphical interface *) Plot[Sin[x],{x,0,2*Pi}] (* Graphics are objects that can be further manipulated, Show combines different graphics together into a single one *) g1=Plot[Sin[x],{x,0,2*Pi}]; g2=Plot[Cos[x],{x,0,2*Pi}]; Show[g1,g2] (* GraphicsGrid on the other hand takes a 2d matrix of Graphics objects and displays them on a grid *) GraphicsGrid[{{g1,g2}}] (* Graphics and Graphics3D can also be used directly to create primitives *) Graphics[{Thick,Green,Rectangle[{0,-1},{2,1}],Red,Disk[],Blue,Circle[{2,0}]}] (* Most Mathematica functions accept a list of options as the last argument. For Plots an useful one is to override the automatic range. Show by default uses the range of the first Graphics so it will cut the second plot here: *) Show[Plot[x^2,{x,0,1}],Plot[x^3,{x,1,2}]] (* Forcing to show all the plotted data *) Show[Plot[x^2,{x,0,1}],Plot[x^3,{x,1,2}], PlotRange->All] (* Very handy for explorations is the ability of having parametric graphs that can manipulated. Manipulate allows for a range of widgets to be displayed next to the output of an expression *) Manipulate[Plot[x^p,{x,0,1}],{{p,1},1,10}] Manipulate[Plot3D[x^p[[1]]+y^p[[2]],{x,0,1},{y,0,1}],{{p,{1,1}},{1,1},{5,5}}] (* Manipulate output is a Dynamic cell, which is special as it get automatically re-evaluated if any of the symbols it capture changes. That's why you can see Manipulate output behaving "weirdly" if you change symbols that are used to compute its output. This allows for all kind of "spreadsheet-like" computations and interactive applications. *) (* Debugging functional programs can be daunting. Mathematica offers a number of primitives that to a degree help. Monitor generates a temporary output that shows the computation in progress. Here the temporary output is a ProgressIndicator graphical object. Evaluations can be aborted with Alt+. *) Monitor[Table[FactorInteger[2^(2*n)+1],{n,1,100}], ProgressIndicator[n, {1,100}]] (* Another example, we assign the value of the function to be minimized to a local symbol, so we can display how it changes as the algorithm progresses *) complexFn=Function[{x,y},(Mod[Mod[x,1],Mod[y,1]+0.1])*Abs[x+y]] Plot3D[complexFn[x,y],{x,-2,2},{y,-2,2}] Block[{temp},Monitor[NMinimize[{temp=complexFn[x,y],x+y==1},{x,y}],N[temp]]] (* Print forces an output from intermediate computations *) Do[Print[Prime[n]],{n,5}] (* Mathematica also supports reflection, via Names, Definition, Information and more *) (* Performance tuning. A first common step is to reduce the number of results Mathematica will keep around for % *) $HistoryLength=2 (* Evaluate current memory usage *) MemoryInUse[] (* Share[] can sometimes shrink the memory usage by making Mathematica realize that certain subexpressions can be shared, it prints the amount of bytes saved *) Share[] (* Reflection can be used to know which symbols are taking the most memory *) Reverse@Sort[{ByteCount[Symbol[#]],#}&/@Names["`*"]] (* Timing operations is simple with AbsoluteTiming *) AbsoluteTiming[Pause[3]] (* Mathematica's symbolic evaluation is relatively slow. Machine numbers operations are faster, but slow compared to other languages. In general Mathematica is not made for high-performance, and if that's needed it's best to directly go to one of the ways it supports external compilation: LibraryLink, CudaLink, and OpenCLLink *) (* On the upside, many list-based operations are trivially parallelizable via Parallelize *) Parallelize[Table[Length[FactorInteger[10^50+n]],{n,20}]] (* The downside is that only a few functions seems to be natively parallelized, mostly image-related, and many others require manual parallelization via domain-splitting. E.G. integrals *) sixDimensionalFunction=Function[{a,b,c,d,e,f},Re[(a*b+c)^d/e+f]]; Total[ParallelTable[NIntegrate[sixDimensionalFunction[a,b,c,d,e,f],{a,-1,1},{b,-1,1},{c,-1,1},{d,-1,1},{e,-1,1},{f,-1+i/4,-1+(i+1)/4}],{i,0,7}]] (* Even plotting ca be parallelized, see http://mathematica.stackexchange.com/questions/30391/parallelize-plotting. Intra-thread communication is expensive, beware of the amount of data you move! *) (* There is a Compile functionality that can translate -some- Mathematica expressions into bytecode or C code, even parallelizing, but it'http://mathematica.stackexchange.com/questions/1096/list-of-compilable-functions *)