**Try it on GitHub:** https://github.com/coffeemug/ts-wolfram

I’ve been using Mathematica to help with learning math (and using math to help learn Mathematica). Too often my Mathematica code has surprising behavior, which means I don’t *really* understand the language or how it’s evaluated. So I thought I’d do two things to get better.

First, write a bunch of Mathematica rules to implement a toy symbolic differentiation operator. Second, write a toy Mathematica interpreter in Typescript that’s good enough to correctly run my custom differentiation code.

Writing a toy differentiator turns out to be shockingly easy. It’s a near verbatim transcription of differentiation rules from any calculus textbook:

```
D[_?NumberQ, x_Symbol] = 0;
D[x_, x_Symbol] = 1;
D[Times[expr1_, expr2_], x_Symbol] =
D[expr1, x] expr2 + D[expr2, x] expr1;
D[Plus[expr1_, expr2_], x_Symbol] = D[expr1, x] + D[expr2, x];
D[Sin[x_], x_Symbol] = Cos[x];
D[Cos[x_], x_Symbol] = -Sin[x];
D[f_Symbol[expr_], x_Symbol] :=
(D[f[x], x] /. x -> expr) * D[expr, x];
D[Power[expr_, p_Integer], x_Symbol] := p expr^(p - 1) * D[expr, x];
```

My first implementation had three more rules (one for each of $dx^k/dx, (cf)’(x)$, and $(1/f)’(x)$), which I later realized you don’t need. These are shortcuts for human differentiators, but they’re automagically covered by the multiplication and exponentiation rules above.

Some problems I encountered while writing these rules: infinite recursion, rules not matching, and rules matching but evaluating to a surprising result. The hobby edition of Mathematica has *some* tools for debugging these problems, but not much. Ultimately I fixed bad rules by staring at the code and thinking really hard. I found `Trace`

impossible to read, and `TracePrint`

(which is supposed to be better) not much better. Also `MatchQ`

is good, but somehow not as useful for debugging as I would have liked.

I first implemented basic parsing of integers, symbols, forms, and arithmetic operators using ts-parsec (which I wrote for this project). In Mathematica `a 2 b 3`

evaluates to `Times[6, a, b]`

because `Times`

has a `Flat`

attribute. To get this behavior I implemented attributes next– `Flat`

, and also `HoldFirst`

, `HoldRest`

, and `HoldAll`

which I’d eventually need. I also exposed `Attributes`

, `SetAttributes`

, and `ClearAttributes`

. These all accept (and `Attributes`

returns) lists, so I added those too. All this was easy enough.

I wanted to implement assignment next so I could say `a = 1`

. In Mathematica even something this simple is implemented by adding `RuleDelayed[HoldPattern[a], 1]`

to `OwnValues`

.^{1} So the next step was to build the skeleton of a term rewriting system.^{2} I first implemented a version of `MatchQ`

that does a deep equal, extended it to handle `Blank[]`

, extended it again to handle `Blank[...]`

, and extended it *again* to handle `Pattern[foo, Blank[...]]`

. The version exposed to the interpreter returns `True`

or `False`

, but under the hood it also returns an environment with matched variables. I built on that next to implement term rewriting.

A really simple rewrite rule is `f[1] /. f[x_] :> x + 1`

. In Mathematica this parses to

`ReplaceAll[f[1], RuleDelayed[f[Pattern[x, Blank[]]], Plus[x, 1]]]`

With my `MatchQ`

implementation there was now enough machinery to get this working. I added operators `/.`

and `:>`

to the grammar, implemented `Replace`

, and built `ReplaceAll`

on top. I tested a bunch of rewrite rules and they all worked! From here it was also easy to add `->`

and `//.`

which I did next.

I had enough machinery to implement assignment. I added `Set`

and `SetDelayed`

, and modified evaluation to apply the rules stored in `OwnValues`

and `DownValues`

. This let me assign values to symbols and define functions! I could now run code like this:

```
fib[1] := 1
fib[2] := 1
fib[x_] := fib[x-2] + fib[x-1]
```

One caveat is that Mathematica inserts rules into value lists in order of specificity, so specific rules are tried before general rules. I initially added rules in order of definition to get assignment working, ~~but then went back and added a simple specificity ordering function~~^{3}.

**EDIT:** I wrote specificity ordering code at late at night, and just realized it was completely broken. I removed it; rules are now processed in the order they’re entered. But everything else still works!

Finally, I added support for `PatternTest`

/`?`

and `NumberQ`

. These were all the features needed to do differentiation!

`ts-wolfram`

I ran `D`

in `ts-wolfram`

on the following examples, and cross-checked with the results in Mathematica:

```
D[1, x]
D[x, x]
D[x^5, x]
D[3 x^2, x]
D[(x + 1) (x + 2), x]
D[x^2 + x^3, x]
D[Cos[x], x]
D[x^3/(x^2 + 1), x]
D[Cos[Cos[x]], x]
D[Cos[Cos[Cos[x]]], x]
D[Cos[x^2 + 1], x]
D[(x + 1)^2, x]
```

Somewhat miraculously, `ts-wolfram`

got correct results on every test! Since I didn’t add any manipulation to simplify expressions, the results weren’t *exactly* the same. For example:

```
D[Cos[Cos[x]], x]
(* ts-wolfram outputs *)
Times[Minus[Sin[Cos[x]]], Minus[Sin[x]]]
(* With `FullForm`, Mathematica outputs *)
Times[Sin[x], Sin[Cos[x]]]
```

This would be easy to fix by adding the rule `Times[Minus[x_], Minus[y_]]=x y`

, but (a) `Times`

is implemented in typescript, and the interpreter doesn’t currently support mixing kernel and userspace rules, and (b) extending the system to simplify algebraic expressions feels like a different undertaking.

Overall, this has been a really fun and instructive project. I built it in four days hacking on it after work, and learned a great deal about Mathematica internals. Of course this is still only scratching the surface, but now I feel a lot less lost when Mathematica doesn’t behave the way I expect. I’m very happy with this outcome!

You can check this by running

`a = 1; OwnValues[a]`

.↩︎If you don’t already understand how this works, reverse engineering it is tricky, even with Mathematica’s great docs. For examples why, see Stack Overflow questions here, here, and here.↩︎

The exact details for how Mathematica handles pattern specificity are not so easy to find. I didn’t try too hard to reverse engineer it; I just did what seemed reasonable and moved on.↩︎

Oct 17, 2024