Mathematical patterns with Hy

I recently read a blog posting that among other things touched patterns in mathematics, how they seem to be obvious to some people while others strugle and if there would be something that schools should do differently in teaching mathematics. Those of you who know me, can guess that this is close to my heart (after all, I keep bringing up Lockhart’s Lament often in day to day discussions). So I dug out my trusty Hy REPL and started playing with sequences just to see what one could do with them.

For a very simple sequence, I would imagine something along the lines:

```=> (seq [n] 1)
[1, 1, 1, 1, 1, 1, 1, 1, 1, 1, ...]
```

This creates an infinite sequence of ones. No matter what point of the sequence I choose, there’s always number 1 there. Such a sequence is sort of boring, so lets spicy up things a bit and try following:

```=> (seq [n] n)
[0, 1, 2, 3, 4, 5, 6, 7, 8, 9, ...]
```

What would it look like, if we were to define something in term of itself (ie. recursively)? Not much different really, we just need a way to bind our sequence to a symbol to give it a name and we’re good to go:

```=> (defseq numbers [n]
...  (cond [(= n 0) 1]
...        [true (+ 1 (get numbers (- n 1)))]))
=> numbers
[1, 2, 3, 4, 5, 6, 7, 8, 9, 10, ...]
```

Previous example defines a sequence called numbers. First element is always 1 (the first line in cond form). The subsequent elements are defined as the previous element increased by one. And thus, we now have all (well, to some rather large value of all) positive integers at our disposal.

How would we create a sequence of even numbers? We could repeat the code that was used to create numbers and increment by 2 instead of by 1. Or we could get creative and just add two numbers streams together:

```=> (sum-of numbers numbers)
[2, 4, 6, 8, 10, 12, 14, 16, 18, 20, ...]
```

Even numbers can be built of top of the previous definition:

```=> (defseq ones [n] 1)
=> (difference-of (sum-of numbers numbers) ones)
[1, 3, 5, 7, 9, 11, 13, 15, 17, 19, ...]
```

Lets build something more complicated based on numbers. Imagine that we are building triangles from little balls. Small triangle is just a single ball, slightly bigger is three balls, next one is 6 balls and so on:

First triangle has only a single ball and it will be our special case. Subsequent triangles have as many balls as the previous one, plus as many balls as there are rows in triangle. So we can just use our previously defined numbers stream and add corresponding element to previous value of triangles:

```=> (defseq triangles [n]
...  (cond [(= n 0) 1]
...        [true (+ (get triangles (- n 1))
...                 (get numbers n))]))
=> triangles
[1, 3, 6, 10, 15, 21, 28, 36, 45, 55, ...]
```

Lets play a bit more with our triangle-numbers and define a new sequence that is formed by subtracting two subsequent triangle-numbers from each other:

```=> (seq [n] (- (get triangles (+ n 1))
...            (get triangles n)))
[2, 3, 4, 5, 6, 7, 8, 9, 10, 11, ...]
```

That seems to be a sequence of each integer again, but starting from 2 this time. What would happen if we were to sum instead of subtract:

```=> (seq [n] (+ (get triangles (+ n 1))
...            (get triangles n)))
[4, 9, 16, 25, 36, 49, 64, 81, 100, 121, ...]
```

Looks awfully same as if we had squared our numbers stream (and skipped first element):

```=> (squared numbers)
[1, 4, 9, 16, 25, 36, 49, 64, 81, 100, ...]
```

Lets do one more sequence, this time rather famous Fibonacci numbers one. The Wikipedia article does much better job explaining various places where it can be found and how it can be used, so I’m not even going to try this time. The main point is that the first two elements are 1 and 1 and every subsequent element is sum of two previous elements:

```=> (defseq fibonacci-numbers [n]
...  (cond [(< n 2) 1]
...        [true (+ (get fibonacci-numbers (- n 1))
...                 (get fibonacci-numbers (- n 2)))]))
=> fibonacci-numbers
[1, 1, 2, 3, 5, 8, 13, 21, 34, 55, ...]
```

How does this work then? There’s a particularly gnarly piece of code that defines a new class called Sequence. Its initializer takes a function as parameter that is used to calculate any given element in sequence. The actual beef is in –getitem– method. This is the method that Hy (and Python too) calls when some code tries to access an element in list. Since we define this method, Sequence can behave like a list:

```(defclass Sequence []
[[--init-- (fn [self func]
(setv (. self func) func)
(setv (. self cache) {0 (.func self 0)})
(setv (. self high-water) 0)
nil)]
[--getitem-- (fn [self n]
(if (hasattr n "start")
(if n.step
(genexpr (get self x) [x (range n.start n.stop n.step)])
(genexpr (get self x) [x (range n.start n.stop 1)]))
(do (assert (>= n 0))
(if (<= n (. self high-water))
(get (. self cache) n)
(do (for [x (range (. self high-water) (inc n))]
(assoc (. self cache) x (.func self x)))
(setv (. self high-water) n)
(get self n))))))]
[--str-- (fn [self]
(.format "[{0}, ...]"
(.join ", " (list-comp (str (get self x)) [x (range 0 10)]))))]
[--repr-- (fn [self]
(.--str-- self))]])
```

Also, in order to speed up things a bit and make recursive definitions easier to use, we have caching logic. When ever –getitem– method is called, it will first check the cache and return value from there if found. If value hasn’t been asked before for that particular index, it’s calculated and placed in cache. Before that, we calculate all the previous values that haven’t been calculated yet. This helps mitigating running over recursion limit in some cases.

–str– and –repr– methods are used to turn 10 first elements of our sequence into textual string. This is used by REPL and str function. Without these methods, we would only get type information, which isn’t particularly informative.

Since I want to shave a little bit of typing, I defined two macros to help defining sequences. One defines a sequence, while other defines a sequence and binds it to a symbol:

```(defmacro seq [param seq-code]
`(Sequence (fn ~param ~seq-code)))

(defmacro defseq [seq-name param seq-code]
`(def ~seq-name (seq ~param ~seq-code)))
```

Mathematical operations for sequences aren’t too complicated to code (albeit it took me several tries to get the interfaces to look pretty enough):

```(defn combine [op seqs]
(seq [n] (reduce (fn [a b] (op a b))
(map (fn [x] (get x n)) seqs))))

(defn sum-of [&rest seqs]
"sum of sequences"
(combine + seqs))

(defn difference-of [&rest seqs]
"difference of sequences"
(combine - seqs))

(defn product-of [&rest seqs]
"product of sequences"
(combine * seqs))

(defn squared [sequence]
"square sequence"
(seq [n] (pow (get sequence n) 2)))
```

The important thing to remember is that trying to process whole sequence in one go is going to take very, very long time (and you’re bound to run out of memory), so all operations need to be lazy (ie. they process only as much as the calling code requests). combine function captures this general pattern. It takes a function that can combine two elements and arbitrary amount of sequences. Then it creates a new sequences which element is result of reducing corresponding elements of all sequences with the given function.

Initially I created the operations using map and genexpr, but that had a nasty side effect: results weren’t sequences anymore. If I wanted to combine result with other sequence again, I had to resort writing little conversions here and there. Now that sequence operations return sequences, combining them is much more easier and there are less surprises involved (sequences are closed under sequence operations is a fancy way of saying this).

For those interested, SICP lectures 6A and 6B talk about streams and how to implement them in Scheme.