This post is about a simple algebraic structure that I have found useful for algorithms that involve searching or sorting based on some ordered weight. I used it a bit in a pair of papers on graph search (2021; 2025), and more recently I used it to implement a version of the Phases type (Easterly 2019) that supported arbitrary keys, inspired by some work by Blöndal (2025a; 2025b) and Visscher (2025).

The algebraic structure in question is a monus, which is a kind of monoid that supports a partial subtraction operation (that subtraction operation, denoted by the symbol ∸, is itself often called a “monus”). However, before giving the full definition of the structure, let me first try to motivate its use. The context here is heap-based algorithms. For the purposes of this post, a heap is a tree that obeys the “heap property”; i.e. every node in the tree has some “weight” attached to it, and every parent node has a weight less than or equal to the weight of each of its children. So, for a tree like the following:

   ┌d
 ┌b┤
a┤ └e
 └c

The heap property is satisfied when $a \leq b$, $a \leq c$, $b \leq d$, and $b \leq e$.

Usually, we also want our heap structure to have an operation like popMin :: Heap k v -> Maybe (v, Heap k v) that returns the least-weight value in the heap paired with the rest of the heap. If this operation is efficient, we can use the heap to efficiently implement sorting algorithms, graph search, etc. In fact, let me give the whole basic interface for a heap here:

popMin :: Heap k v -> Maybe (v, Heap k v)
insert :: k -> v -> Heap k v -> Heap k v
empty  :: Heap k v

Using these functions it’s not hard to see how we can implement a sorting algorithm:

sortOn :: (a -> k) -> [a] -> [a]
sortOn k = unfoldr popMin . foldr (\x -> insert (k x) x) empty

The monus becomes relevant when the weight involved is some kind of monoid. This is quite a common situation: if we were using the heap for graph search (least-cost paths or something), we would expect the weight to correspond to path costs, and we would expect that we can add the costs of paths in a kind of monoidal way. Furthermore, we would probably expect the monoidal operations to relate to the order in some coherent way. A monus (Amer 1984) is an ordered monoid where the order itself can be defined in terms of the monoidal operations¹:

$$x \leq y \iff \exists z. \; y = x \bullet z$$

I read this definition as saying “$x$ is less than $y$ iff there is some $z$ that fits between $x$ and $y$”. In other words, the gap between $x$ and $y$ has to exist, and it is equal to $z$.

Notice that this order definition won’t work for groups like $(\mathbb{Z},+,0)$. For a group, we can always find some $z$ that will fit the existential (specifically, $z = (- x) \bullet y$). Monuses, then, tend to be positive monoids: in fact, many monuses are the positive cones of some group ($(\mathbb{N},+,0)$ is the positive cone of $(\mathbb{Z},+,0)$).

We can derive a lot of useful properties from this basic structure. For example, if the order above is total, then we can derive the binary subtraction operator mentioned above:

$$x ∸ y = \begin{cases} z, & \text{if } y \leq x \text{ and } x = y \bullet z \\ 0, & \text{otherwise.} \end{cases}$$

If we require the underlying monoid to be commutative, and we further require the derived order to be total and antisymmetric, we get the particular flavour of monus I worked with in a pair of papers on graph search (2021; 2025). In this post I will actually be working with a weakened form of the algebra that I will define shortly.

Getting back to our heap from above, with this new order defined, we can see that the heap property actually tells us something about the makeup of the weights in the tree. Instead of every child just having a weight equal to some arbitrary quantity, the heap property tells us that each child weight has to be made up of the combination of its parent’s weight and some difference.

   ┌d              ┌b•(d∸b)
 ┌b┤       ┌a•(b∸a)┤
a┤ └e  =  a┤       └b•(e∸b)
 └c        └a•(c∸a)

This observation gives us an opportunity for a different representation: instead of storing the full weight at each node, we could instead just store the difference.

     ┌d∸b
 ┌b∸a┤
a┤   └e∸b
 └c∸a

Just in terms of data structure design, I prefer this version: if we wanted to write down a type of heaps using the previous design, we would first define the type of trees, and then separately write a predicate corresponding to the heap property. With this design, it is impossible to write down a tree that doesn’t satisfy the heap property.

More practically, though, using this algebraic structure when working with heaps enables some optimisations that might be difficult to implement otherwise. The strength of this representation is that it allows for efficient relative and global computation: now, if we wanted to add some quantity to every weight in the tree, we can do it just by adding the weight to the root node.

Monuses in Haskell

To see some examples of how to use this pattern, let’s first write a class for Haskell monuses:

class (Semigroup a, Ord a) => Monus a where
  (∸) :: a -> a -> a

You’ll notice that we’re requiring semigroup here, not monoid. That’s because one of the nice uses of this pattern actually works with a weakening of the usual monus algebra; this weakening only requires semigroup, and the following two laws.

$$x \leq y \implies x \bullet (y ∸ x) = y \quad \quad \quad \quad \quad x \leq y \implies z \bullet x \leq z \bullet y$$

A straightforward monus instance is the following:

instance (Num a, Ord a) => Monus (Sum a) where
  (∸) = (-)

Pairing Heaps in Haskell

Next, let’s look at a simple heap implementation. I will always go for pairing heaps (Fredman et al. 1986) in Haskell; they are extremely simple to implement, and (as long as you don’t have significant persistence requirements) their performance seems to be the best of the available pointer-based heaps (Larkin, Sen, and Tarjan 2013). Here is the type definition:

data Root k v = Root !k v [Root k v]
type Heap k v = Maybe (Root k v)

A Root is a non-empty pairing heap; the Heap type represents possibly-empty heaps. The key function to implement is the merging of two heaps; we can accomplish this as an implementation of the semigroup <>.

instance Monus k => Semigroup (Root k v) where
  Root xk xv xs <> Root yk yv ys
    | xk <= yk  = Root xk xv (Root (yk ∸ xk) yv ys : xs)
    | otherwise = Root yk yv (Root (xk ∸ yk) xv xs : ys)

The only difference between this and a normal pairing heap merge is the use of ∸ in the key of the child node (yk ∸ xk and xk ∸ yk). This difference ensures that each child only holds the difference of the weight between itself and its parent.

It’s worth working out why the weakened monus laws above are all we need in order to maintain the heap property on this structure.

The rest of the methods are implemented the same as their implementations on a normal pairing heap. First, we have the pairing merge of a list of heaps, here given as an implementation of the semigroup method sconcat:

  sconcat (x1 :| []) = x1
  sconcat (x1 :| [x2]) = x1 <> x2
  sconcat (x1 :| x2 : x3 : xs) = (x1 <> x2) <> sconcat (x3 :| xs)

merges :: Monus k => [Root k v] -> Heap k v
merges = fmap sconcat . nonEmpty

The pattern of this two-level merge is what gives the pairing heap its excellent performance.

The Heap type derives its monoid instance from the monoid instance on Maybe (instance Semigroup a => Monoid (Maybe a)), so we can implement insert like so:

insert :: Monus k => k -> v -> Heap k v -> Heap k v
insert k v hp = Just (Root k v []) <> hp

And popMin is also relatively simple:

delay :: Semigroup k => k -> Heap k v -> Heap k v
delay by = fmap (\(Root k v xs) -> Root (by <> k) v xs)

popMin :: Monus k => Heap k v -> Maybe (v,Heap k v)
popMin = fmap (\(Root k v xs) -> (v, k `delay` merges xs))

Notice that we delay the rest of the heap, because all of its entries need to be offset by the weight of the previous root node. Thankfully, because we’re only storing the differences, we can “modify” every weight by just increasing the weight of the root, making this an $\mathcal{O}(1)$ operation.

Finally, we can implement heap sort like so:

sortOn :: Monus k => (a -> k) -> [a] -> [a]
sortOn k = unfoldr popMin . foldl' (\xs x -> insert (k x) x xs) Nothing

And it does indeed work:

>>> sortOn Sum [3,4,2,5,1]
[1,2,3,4,5]

Here is a trace of the output:

Input           Heap:
list:

[3,4,2,5,1]


[4,2,5,1]       3


[2,5,1]         3─1


[5,1]           2─1─1


[1]              ┌3
                2┤
                 └1─1


[]                 ┌3
                1─1┤
                   └1─1


Output          Heap:
list:

[]                 ┌3
                1─1┤
                   └1─1


[1]              ┌3
                2┤
                 └1─1


[1,2]            ┌2
                3┤
                 └1


[1,2,3]         4─1


[1,2,3,4]       5


[1,2,3,4,5]

While the heap implementation presented here is pretty efficient, note that we could significantly improve its performance with a few optimisations: first, we could unpack all of the constructors, using a custom list definition in Root instead of Haskell’s built-in lists; second, in foldl' we could avoid the Maybe wrapper by building a non-empty heap. There are probably more small optimisations available as well.

Retrieving a Normal Heap

A problem with the definition of sortOn above is that it requires a Monus instance on the keys, but it only really needs Ord. It seems that by switching to the Monus-powered heap we have lost some generality.

Luckily, there are two monuses we can use to solve this problem:

instance Ord a => Monus (Max a) where
  x ∸ y = x

instance Ord a => Monus (Last a) where
  x ∸ y = x

The Max semigroup uses the max operation, and the Last semigroup returns its second operand.

instance Ord a => Semigroup (Max a) where
  Max x <> Max y = Max (max x y)

instance Semigroup (Last a) where
  x <> y = y

While the Monus instances here might seem degenerate, they do actually satisfy the Monus laws as given above.

x <= y ==> x <> (y ∸ x) = y
Max x <= Max y ==> Max x <> (Max y ∸ Max x) = Max y
Max x <= Max y ==> Max x <> Max y = Max y
Max x <= Max y ==> Max (max x y) = Max y
x <= y ==> max x y = y

x <= y ==> z <> x <= z <> y
Max x <= Max y ==> Max z <> Max x <= Max z <> Max y
Max x <= Max y ==> Max (max z x) <= Max (max z y)
x <= y ==> max z x <= max z y

x <= y ==> x <> (y ∸ x) = y
Last x <= Last y ==> Last x <> (Last y ∸ Last x) = Last y
Last x <= Last y ==> Last x <> Last y = Last y
Last x <= Last y ==> Last y = Last y

x <= y ==> z <> x <= z <> y
Last x <= Last y ==> Last z <> Last x <= Last z <> Last y
Last x <= Last y ==> Last x <= Last y

Either Max or Last will work; semantically, there’s no real difference. Last avoids some comparisons, so we can use that:

sortOn' :: Ord k => (a -> k) -> [a] -> [a]
sortOn' k = sortOn (Last . k)

Phases as a Pairing Heap

The Phases applicative (Easterly 2019) is an Applicative transformer that allows reordering of Applicative effects in an easy-to-use, high-level way. The interface looks like this:

phase     :: Natural -> f a -> Phases f a
runPhases :: Applicative f => Phases f a -> f a

instance Applicative f => Applicative (Phases f)

And we can use it like this:

phased :: IO String
phased = runPhases $ sequenceA
  [ phase 3 $ emit 'a'
  , phase 2 $ emit 'b'
  , phase 1 $ emit 'c'
  , phase 2 $ emit 'd'
  , phase 3 $ emit 'e' ]
  where emit c = putChar c >> return c

>>> phased
cbdae
"abcde"

The above computation performs the effects in the order dictated by their phases (this is why the characters are printed out in the order cbdae), but the pure value (the returned string) has its order unaffected.

I have written about this type before, and in a handful of papers (2021; Gibbons et al. 2022; Gibbons et al. 2023), but more recently Blöndal (2025a) started looking into trying to use the Phases pattern with arbitrary ordered keys (Visscher 2025; Blöndal 2025b). There are a lot of different directions you can go from the Phases type; what interested me most immediately was the idea of implementing the type efficiently using standard data-structure representations. If our core goal here is to order some values according to a key, then that is clearly a problem that a heap should solve: enter the free applicative pairing heap.

Here is the type’s definition:

data Heap k f a where
  Pure :: a -> Heap k f a
  Root :: !k -> (x -> y -> a) -> f x -> Heaps k f y -> Heap k f a

data Heaps k f a where
  Nil :: Heaps k f ()
  App :: !k -> f x -> Heaps k f y -> Heaps k f z -> Heaps k f (x,y,z)

We have had to change a few aspects of the original pairing heap, but the overall structure remains. The entries in this heap are now effectful computations: the fs. The data structure also contains some scaffolding to reconstruct the pure values “inside” each effect when we actually run the heap.

The root-level structure is the Heap: this can either be Pure (corresponding to an empty heap: notice that, though this constructor has some contents (the a), it is still regarded as “empty” because it contains no effects (f)); or a Root, which is a singleton value, paired with the list of sub-heaps represented by the Heaps type. We’re using the usual Yoneda-ish trick here to allow the top-level data type to be parametric and a Functor, by storing the function x -> y -> a.

The Heaps type then plays the role of [Root k v] in the previous pairing heap implementation; here, we have inlined all of the constructors so that we can get all of the types to line up. Remember, this is a heap of effects, not of pure values: the pure values need to be able to be reconstructed to one single top-level a when we run the heap at the end.

Merging two heaps happens in the Applicative instance itself:

instance Functor (Heap k f) where
  fmap f (Pure x) = Pure (f x)
  fmap f (Root k c x xs) = Root k (\a b -> f (c a b)) x xs

instance Monus k => Applicative (Heap k f) where
  pure = Pure
  Pure f <*> xs = fmap f xs
  xs <*> Pure f = fmap ($ f) xs

  Root xk xc xs xss <*> Root yk yc ys yss
    | xk <= yk  = Root xk (\a (b,c,d) -> xc a d (yc b c)) xs (App (yk ∸ xk) ys yss xss)
    | otherwise = Root yk (\a (b,c,d) -> xc b c (yc a d)) ys (App (xk ∸ yk) xs xss yss)

To actually run the heap we will use the following two functions:

merges :: (Monus k, Applicative f) => Heaps k f a -> Heap k f a
merges Nil = Pure ()
merges (App k1 e1 t1 Nil) = Root k1 (,,()) e1 t1
merges (App k1 e1 t1 (App k2 e2 t2 xs)) =
   (Root k1 (\a b cd es -> (a,b, cd es)) e1 t1 <*> Root k2 (,,) e2 t2) <*> merges xs

runHeap :: (Monus k, Applicative f) => Heap k f a -> f a
runHeap (Pure x) = pure x
runHeap (Root _ c x xs) = liftA2 c x (runHeap (merges xs))

And we can lift a computation into Phases like so:

phase :: k -> f a -> Heap k f a
phase k xs = Root k const xs Nil

Stabilising Phases

There’s a problem. A heap sort based on a pairing heap isn’t stable. That means that the order of effects here can vary for two effects in the same phase. If we look back to the example with the strings we saw above, that means that outputs like cdbea would be possible (in actual fact, we don’t get any reordering in this particular example, but that’s just an accident of the way the applicative operators are associated under the hood).

This is problematic because we would expect effects in the same phase to behave as if they were normal applicative effects, sequenced according to their syntactic order. It also means that the applicative transformer breaks the applicative laws, because effects might be reordered according to the association of the applicative operators, which should lawfully be associative.

To make the sort stable, we could layer the heap effect with some state effect that would tag each effect with its order. However, that would hurt efficiency and composability: it would force us to linearise the whole heap sort procedure, where currently different branches of the tree can compute completely independently of each other. The solution comes in the form of another monus: the key monus.

data Key k = !k :* {-# UNPACK #-} !Int deriving (Eq, Ord)

A Key k is some ordered key k coupled with an Int that represents the offset between the original position and the current position of the key. In this way, when two keys compare as equal, we can cascade on to compare their original positions, thereby maintaining their original order when there is ambiguity caused by a key collision. However, in contrast to the approach of walking over the data once and tagging it all with positions, this approach keeps the location information completely local: we never need to know that some key is in the $n$th position in the original sequence, only that it has moved $n$ steps from its original position.

The instances are as follows:

instance Semigroup (Key k) where
  (xk :* xi) <> (yk :* yi) = yk :* (xi + yi)

instance Ord k => Monus (Key k) where
  (xk :* xi) ∸ (yk :* yi) = xk :* (xi - yi)

This instance is basically a combination of the Last semigroup and the $(\mathbb{Z}, +, 0)$ group. We could make a slightly more generalised version of Key that is the combination of any monus and $\mathbb{Z}$, but since I’m only going to be using this type for simple sorting-like algorithms I will leave that generalisation for another time.

The stable heap type is as follows:

data Stable k f a
  = Stable { size :: {-# UNPACK #-} !Int
           , heap :: !(Heap (Key k) f a) }

We need to track the size of the heap so that we can supply the right-hand operand with their offsets. Because we’re storing differences, we can add an offset to every entry in a heap in $\mathcal{O}(1)$ time by simply adding to the root:

delayKey :: Int -> Heap (Key k) f a -> Heap (Key k) f a
delayKey _ hp@(Pure _) = hp
delayKey n (Root (k :* m) c x xs) = Root (k :* (n + m)) c x xs

Finally, using this we can implement the Applicative instance and the rest of the interface:

instance Ord k => Applicative (Stable k f) where
  pure = Stable 0 . pure
  Stable n xs <*> Stable m ys = Stable (n+m) (xs <*> delayKey n ys)

runStable :: (Applicative f, Ord k) => Stable k f a -> f a
runStable = runHeap . heap

stable :: Ord k => k -> f a -> Stable k f a
stable k fa = Stable 1 (phase (k :* 0) fa)

This is a pure, optimally efficient implementation of Phases ordered by an arbitrary total-ordered key.

Local Computation in a Monadic Heap

In (2021), I developed a monadic heap based on the free monad transformer.

newtype Search k a = Search { runSearch :: [Either a (k, Search k a)] }

This type is equivalent to the free monad transformer over the list monad and (,) k functor (i.e. the writer monad).

Search k a ≅ FreeT ((,) k) [] a

In the paper (2021) we extended the type to become a full monad transformer, replacing lists with ListT. This let us order the effects according to the weight k; however, for this example we only need the simplified type, which lets us order the values according to k.

This Search type follows the structure of a pairing heap (although not as closely as the version above). However, this type is interesting because semantically it needs the weights to be stored as differences, rather than absolute weights. As a free monad transformer, the Search type layers effects on top of each other; we can later interpret those layers by collapsing them together using the monadic join. In the case of Search, those layers are drawn from the list monad and the (,) k functor (writer monad). That means that if we have some heap representing the tree from above:

Search [ Right (a, Search [ Right (b, Search [ Right (d, Search [Left x])
                                             , Right (e, Search [Left y])])
                          , Right (c, Search [Left z])])]

When we collapse this computation down to the leaves, the weights we will get are the following:

[(a <> b <> d, x), (a <> b <> e, y), (a <> c, z)]

So, if we want the weights to line up properly, we need to store the differences.

mergeS :: Monus k => [(k, Search k a)] -> Maybe (k, Search k a)
mergeS [] = Nothing
mergeS (x:xs) = Just (mergeS' x xs)
  where
    mergeS' x1 [] = x1
    mergeS' x1 [x2] = x1 <+> x2
    mergeS' x1 (x2:x3:xs) = (x1 <+> x2) <+> mergeS' x3 xs

    (xw, Search xs) <+> (yw, Search ys)
      | xw <= yw  = (xw, Search (Right (yw ∸ xw, Search ys) : xs))
      | otherwise = (yw, Search (Right (xw ∸ yw, Search xs) : ys))

popMins :: Monus k => Search k a -> ([a], Maybe (k, Search k a))
popMins = fmap mergeS . partitionEithers . runSearch

Conclusion

The technique of “don’t store the absolute value, store the difference” seems to be generally quite useful; I think that monuses are a handy algebra to keep in mind whenever that technique looks like it might be needed. The Key monus above is closely related to the factorial numbers, and the trick I used in this post.

References

New paper: “Hyperfunctions: Communicating Continuations”, by myself and Nicolas Wu, will be published at POPL 2026.

The preprint is available here.

The work contained in the paper started with a post on this blog in 2021. I had read a paper by Launchbury, Krstić, and Sauerwein (2000) and I recognised that their hyperfunction construction was quite similar to some types I had used to implement breadth-first traversal (in particular, the Queue in this post). After that, I started seeing hyperfunctions in lots of different settings, rediscovered by different authors, and almost always accompanied by some remark about how difficult it was to understand the type.

My hope with this paper is to clarify and explain what hyperfunctions can do, and where they might be useful. Ideally, the paper will save some future programmer from having to reinvent the type.

References

New paper: “Formalising Graph Algorithms with Coinduction”, by myself and Nicolas Wu, will be published at POPL 2025.

The preprint is available here.

The talk is here.

The paper is about representing graphs (especially in functional languages). We argue in the paper that graphs are naturally coinductive, rather than inductive, and that many of the problems with graphs in functional languages go away once you give up on induction and pattern-matching, and embrace the coinductive way of doing things.

Of course, coinduction comes with its own set of problems, especially when working in a total language or proof assistant. Another big focus of the paper was figuring out a representation that was amenable to formalisation (we formalised the paper in Cubical Agda). Picking a good representation for formalisation is a tricky thing: often a design decision you make early on only looks like a mistake after a few thousand lines of proofs, and modern formal proofs tend to be brittle, meaning that it’s difficult to change an early definition without also having to change everything that depends on it. On top of this, we decided to use quotients for an important part of the representation, and (as anyone who’s worked with quotients and coinduction will tell you) productivity proofs in the presence of quotients can be a real pain.

All that said, I think the representation we ended up with in the paper is quite nice. We start with a similar representation to the one we had in our ICFP paper in 2021: a graph over vertices of type a is simply a function a -> [a] that returns the neighbours of a supplied vertex (this is the same representation as in this post). Despite the simplicity, it turns out that this type is enough to implement a decent number of search algorithms. The really interesting thing is that the arrow methods (from Control.Arrow) work on this type, and they define an algebra on graphs similar to the one from Mokhov (2017). For example, the <+> operator is the same as the overlay operation in Mokhov (2017).

That simple type gets expanded upon and complicated: eventually, we represent a possibly-infinite collection as a function that takes a depth and then returns everything in the search space up to that depth. It’s a little like representing an infinite list as the partial application of the take function. The paper spends a lot of time picking an algebra that properly represents the depth, and figuring out coherency conditions etc.

One thing I’m especially proud of is that all the Agda code snippets in the paper are hyperlinked to a rendered html version of the code. Usually, when I want more info on some code snippet in a paper, I don’t really want to spend an hour or so downloading some artefact, installing a VM, etc. What I actually want is just to see all of the definitions the snippet relies on, and the 30 or so lines of code preceding it. With this paper, that’s exactly what you get: if you click on any Agda code in the paper, you’re brought to the source of that code block, and every definition is clickable so you can browse without having to install anything.

I think the audience for this paper is anyone who is interested in graphs in functional languages. It should be especially interesting to people who have dabbled in formalising some graphs, but who might have been stung by an uncooperative proof assistant. The techniques in the second half of the paper might help you to convince Agda (or Idris, or Rocq) to accept your coinductive and quotient-heavy arguments.

References

New paper: “Algebraic Effects Meet Hoare Logic in Cubical Agda”, by myself, Zhixuan Yang, and Nicolas Wu, will be published at POPL 2024.

The preprint is available here.

Here’s a cool trick:

minimum :: Ord a => [a] -> a
minimum = head . sort

This is $\mathcal{O}(n)$ in Haskell, not $\mathcal{O}(n \log n)$ as you might expect. And this isn’t because Haskell is using some weird linear-time sorting algorithm; indeed, the following is $\mathcal{O}(n \log n)$:

maximum :: Ord a => [a] -> a
maximum = last . sort

No: since the implementation of minimum above only demands the first element of the list, and since sort has been carefully implemented, only a linear amount of work will be done to retrieve it.

It’s not easy to structure programs to have the same property as sort does above: to be maximally lazy, such that unnecessary work is not performed. Today I was working on a maximally lazy implementation of the following program:

groupOn :: Eq k => (a -> k) -> [a] -> [(k,[a])]
groupOn = ...

>>> groupOn (`rem` 2) [1..5]
[(1,[1,3,5]),(0,[2,4])]

>>> groupOn (`rem` 3) [5,8,3,6,2]
[(2,[5,8,2]),(0,[3,6])]

This function groups the elements of a list according to some key function. The desired behaviour here is a little subtle: we don’t want to just group adjacent elements, for instance.

groupOn (`rem` 3) [5,8,3,6,2] ≢ [(2,[5,8]),(0,[3,6]),(2,[2])]

And we don’t want to reorder the elements of the list by the keys:

groupOn (`rem` 3) [5,8,3,6,2] ≢ [(0,[3,6]),(2,[5,8,2])]

These constraints make it especially tricky to make this function lazy. In fact, at first glance, it seems impossible. What should, for instance, groupOn id [1..] return? It can’t even fill out the first group, since it will never find another 1. However, it can fill out the first key. And, in fact, the second. And it can fill out the first element of the first group. Precisely:

groupOn id [1..] ≡ [(1,1:⊥), (2,2:⊥), (3,3:⊥), ...

Another example is groupOn id (repeat 1), or groupOn id (cycle [1,2,3]). These each have partially-defined answers:

groupOn id (repeat 1)      ≡ (1,repeat 1):⊥

groupOn id (cycle [1,2,3]) ≡ (1,repeat 1):(2,repeat 2):(3,repeat 3):⊥

So there is some kind of well-defined lazy semantics for this function. The puzzle I was interested in was defining an efficient implementation for these semantics.

The Slow Case

The first approximation to a solution I could think of is the following:

groupOn :: Ord k => (a -> k) -> [a] -> [(k, [a])]
groupOn k = Map.toList . Map.fromListWith (++) . map (\x -> (k x, [x]))

In fact, if you don’t care about laziness, this is probably the best solution: it’s $\mathcal{O}(n \log n)$, it performs well (practically as well as asymptotically), and it has the expected results.

However, there are problems. Primarily this solution cares about ordering, which we don’t want. We want to emit the results in the same order that they were in the original list, and we don’t necessarily want to require an ordering on the elements (for the efficient solution we will relax this last constraint).

Instead, let’s implement our own “map” type that is inefficient, but more general.

type Map a b = [(a,b)]

insertWith :: Eq a => (b -> b -> b) -> a -> b -> Map a b -> Map a b
insertWith f k v [] = [(k,v)]
insertWith f k v ((k',v'):xs)
  | k == k'   = (k',f v v') : xs
  | otherwise = (k',v') : insertWith f k v xs

groupOn :: Eq k => (a -> k) -> [a] -> [(k, [a])]
groupOn k = foldr (uncurry (insertWith (++))) [] . map (\x -> (k x, [x]))

The problem here is that it’s not lazy enough. insertWith is strict in its last argument, which means that using foldr doesn’t gain us anything laziness-wise.

There is some extra information we can use to drive the result: we know that the result will have keys that are in the same order as they appear in the list, with duplicates removed:

groupOn :: Eq k => (a -> k) -> [a] -> [(k, [a])]
groupOn k xs = map _ ks
  where
    ks = map k xs

From here, we can get what the values should be from each key by filtering the original list:

groupOn :: Eq k => (a -> k) -> [a] -> [(k,[a])]
groupOn key xs = map (\k -> (k, filter ((k==) . key) xs)) (nub (map key xs))

Using a kind of Schwartzian transform yields the following slight improvement:

groupOn :: Eq k => (a -> k) -> [a] -> [(k,[a])]
groupOn key xs = map (\k -> (k , map snd (filter ((k==) . fst) ks))) (nub (map fst ks))
  where
    ks = map (\x -> (key x, x)) xs

But this traverses the same list multiple times unnecessarily. The problem is that we’re repeating a lot of work between nub and the rest of the algorithm.

The following is much better:

groupOn :: Eq k => (a -> k) -> [a] -> [(k,[a])]
groupOn key = go . map (\x -> (key x, x))
  where
    go [] = []
    go ((k,x):xs) = (k,x:map snd y) : go ys
      where
        (y,ys) = partition ((k==).fst) xs

First, we perform the Schwartzian transform optimisation. The work of the algorithm is done in the go helper. The idea is to filter out duplicates as we encounter them: when we encounter (k,x) we can keep it immediately, but then we split the rest of the list into the components that have the same key as this element, and the ones that differ. The ones that have the same key can form the collection for this key, and those that differ are what we recurse on.

This partitioning also avoids re-traversing elements we know to be already accounted for in a previous group. I think that this is the most efficient (modulo some inlining and strictness improvements) algorithm that can do groupOn with just an Eq constraint.

A Faster Version

The reason that the groupOn above is slow is that every element returned has to traverse the entire rest of the list to remove duplicates. This is a classic pattern of quadratic behaviour: we can improve it by using the same trick as quick sort, by partitioning the list into lesser and greater elements on every call.

groupOnOrd :: Ord k => (a -> k) -> [a] -> [(k,[a])]
groupOnOrd key = go . map (\x -> (key x, x))
  where
    go [] = []
    go ((k,x):xs) = (k,x:e) : go lt ++ go gt
      where
        (e,lt,gt) = foldr split ([],[],[]) xs
        split ky@(k',y) ~(e,lt,gt) = case compare k' k of
          LT -> (e, ky:lt, gt)
          EQ -> (y:e, lt, gt)
          GT -> (e, lt, ky:gt)

While this is $\mathcal{O}(n \log n)$, and it does group elements, it also reorders the underlying list. Let’s fix that by tagging the incoming elements with their positions, and then using those positions to order them back into their original configuration:

groupOnOrd :: Ord k => (a -> k) -> [a] -> [(k,[a])]
groupOnOrd k = map (\(_,k,xs) -> (k,xs)) . go . zipWith (\i x -> (i, k x, x)) [0..]
  where
    go [] = []
    go ((i, k, x):xs) = (i, k, x : e) : merge (go l) (go g)
      where
        (e, l, g) = foldr split ([],[],[]) xs

        split ky@(_,k',y) ~(e, l, g) = case compare k' k of
          LT -> (e  , ky : l,      g)
          EQ -> (y:e,      l,      g)
          GT -> (e  ,      l, ky : g)

    merge [] gt = gt
    merge lt [] = lt
    merge (l@(i,_,_):lt) (g@(j,_,_):gt)
      | i <= j    = l : merge lt (g:gt)
      | otherwise = g : merge (l:lt) gt

This is close, but still not right. This isn’t yet lazy. The merge function is strict in both arguments.

However, we have all the information we need to unshuffle the lists without having to inspect them. In split, we know which direction we put each element: we can store that info without using indices.

groupOnOrd :: Ord k => (a -> k) -> [a] -> [(k,[a])]
groupOnOrd k = catMaybes . go . map (\x -> (k x, x))
  where
    go [] = []
    go ((k,x):xs) = Just (k, x : e) : merge m (go l) (go g)
      where
        (e, m, l, g) = foldr split ([],[],[],[]) xs

        split ky@(k',y) ~(e, m, l, g) = case compare k' k of
          LT -> (  e, LT : m, ky : l,      g)
          EQ -> (y:e, EQ : m,      l,      g)
          GT -> (  e, GT : m,      l, ky : g)

    merge []        lt     gt     = []
    merge (EQ : xs) lt     gt     = Nothing : merge xs lt gt
    merge (LT : xs) (l:lt) gt     = l       : merge xs lt gt
    merge (GT : xs) lt     (g:gt) = g       : merge xs lt gt

What we generate here is a [Ordering]: this list tells us the result of all the compare operations on the input list. Then, in merge, we invert the action of split, rebuilding the original list without inspecting either lt or gt.

And this solution works! It’s $\mathcal{O}(n \log n)$, and fully lazy.

>>> map fst . groupOnOrd id $ [1..]
[1..]

>>> groupOnOrd id $ cycle [1,2,3]
(1,repeat 1):(2,repeat 2):(3,repeat 3):⊥

>>> groupOnOrd (`rem` 3) [1..]
(1,[1,4..]):(2,[2,5..]):(0,[3,6..]):⊥

The finished version of these two functions, along with some benchmarks, is available here.

I haven’t written much on this blog recently: since starting a PhD all of my writing output has gone towards paper drafts and similar things. Recently, though, I’ve been thinking about streams, monoids, and comonads and I haven’t managed to wrangle those thoughts into something coherent enough for a paper. This blog post is a collection of those (pretty disorganised) thoughts. The hope is that writing them down will force me to clarify things, but consider this a warning that the rest of this post may well be muddled and confusing.

Streams

The first thing I want to talk about is streams.

record Stream (A : Type) : Type where
  coinductive
  field head : A
        tail : Stream A

This representation is coinductive: the type above contains infinite values. Agda, unlike Haskell, treats inductive and coinductive types differently (this is why we need the coinductive keyword in the definition). One of the differences is that it doesn’t check termination for construction of these values:

alternating : Stream Bool
alternating .head       = true
alternating .tail .head = false
alternating .tail .tail = alternating

alternating :: [Bool]
alternating = True : False : alternating

We have the equivalent in Haskell on the right. We’re also using some fancy syntax for the Agda code: copatterns (Abel and Pientka 2013).

Note that this type is only definable in a language with some notion of laziness. If we tried to define a value like alternating above in OCaml we would loop. Haskell has no problem, and Agda—through its coinduction mechanism—can handle it as well.

Update 4-5-22: thanks to Arnaud Spiwack (@aspiwack) for correcting me on this, it turns out the definition of alternating above can be written in Ocaml, even without laziness. Apparently Ocaml has a facility for strict cyclic data structures. Also, I should be a little more precise with what I’m saying above: even without the extra facility for strict cycles, you can of course write a lazy list with some kind of lazy wrapper type.

There is, however, an isomorphic type that can be defined without coinduction:

ℕ-Stream : Type → Type
ℕ-Stream A = ℕ → A

ℕ-alternating : ℕ-Stream Bool
ℕ-alternating zero          = true
ℕ-alternating (suc zero)    = false
ℕ-alternating (suc (suc n)) = ℕ-alternating n

(notice that, in this form, the function ℕ-alternating is the same function as even : ℕ → Bool)

In fact, we can convert from the coinductive representation to the inductive one. This conversion function is more familiarly recognisable as the indexing function:

_[_] : Stream A → ℕ-Stream A
xs [ zero  ] = xs .head
xs [ suc n ] = xs .tail [ n ]

I’m not just handwaving when I say the two representations are isomorphic: we can prove this isomorphism, and, in Cubical Agda, we can use this to transport programs on one representation to the other.

tabulate : ℕ-Stream A → Stream A
tabulate xs .head = xs zero
tabulate xs .tail = tabulate (xs ∘ suc)

stream-rinv : (xs : Stream A) → tabulate (xs [_]) ≡ xs
stream-rinv xs i .head = xs .head
stream-rinv xs i .tail = stream-rinv (xs .tail) i

stream-linv : (xs : ℕ-Stream A) (n : ℕ) → tabulate xs [ n ] ≡ xs n
stream-linv xs zero    = refl
stream-linv xs (suc n) = stream-linv (xs ∘ suc) n

stream-reps : ℕ-Stream A ⇔ Stream A
stream-reps .fun = tabulate
stream-reps .inv = _[_]
stream-reps .rightInv = stream-rinv
stream-reps .leftInv xs = funExt (stream-linv xs)

One final observation about streams: another way to define a stream is as the cofree comonad of the identity functor.

record Cofree (F : Type → Type) (A : Type) : Type where
  coinductive
  field root : A
        step : F (Cofree F A)

𝒞-Stream : Type → Type
𝒞-Stream = Cofree id

Concretely, the Cofree F A type is a possibly infinite tree, with branches shaped like F, and internal nodes labelled with A. It has the following characteristic function:

{-# NON_TERMINATING #-}
trace : ⦃ _ : Functor 𝐹 ⦄ → (A → B) → (A → 𝐹 A) → A → Cofree 𝐹 B
trace ϕ ρ x .root = ϕ x
trace ϕ ρ x .step = map (trace ϕ ρ) (ρ x)

Like how the free monad turns any functor into a monad, the cofree comonad turns any functor into a comonad. Comonads are less popular and widely-used than monads, as there are less well-known examples of them. I have found it helpful to think about comonads through spatial analogies. A lot of comonads can represent a kind of walk through some space: the extract operation tells you “what is immediately here”, and the duplicate operation tells you “what can I see from each point”. For the stream, these two operations are inhabited by head and the following:

duplicate : Stream A → Stream (Stream A)
duplicate xs .head = xs
duplicate xs .tail = duplicate (xs .tail)

Generalising Streams

There were three key observations in the last section:

1. Streams are coinductive. This requires a different termination checker in Agda, and a different evaluation model in strict languages.
2. They have an isomorphic representation based on indexing. This isomorphic representation doesn’t need coinduction or laziness.
3. They are a special case of the cofree comonad.

Going forward, we’re going to look at generalisations of streams, and we’re going to see what these observations mean in the contexts of the new generalisations.

The thing we’ll be generalising is the index of the stream. Currently, streams are basically structures that assign a value to every ℕ: what does a stream of—for instance—rational numbers look like? To drive the intuition for this generalisation let’s first look at the comonad instance on the ℕ-Stream type:

ℕ-extract : ℕ-Stream A → A
ℕ-extract xs = xs zero

ℕ-duplicate : ℕ-Stream A → ℕ-Stream (ℕ-Stream A)
ℕ-duplicate xs zero    = xs
ℕ-duplicate xs (suc n) = ℕ-duplicate (xs ∘ suc) n

This is the same instance as is on the Stream type, transported along the isomorphism between the two types (we could have transported the instance automatically, using subst or transport; I have written it out here manually in full for illustration purposes).

The ℕ-duplicate method here can changed a little to reveal something interesting:

ℕ-duplicate₂ : ℕ-Stream A → ℕ-Stream (ℕ-Stream A)
ℕ-duplicate₂ xs zero    m = xs m
ℕ-duplicate₂ xs (suc n) m = ℕ-duplicate₂ (xs ∘ suc) n m

ℕ-duplicate₃ : ℕ-Stream A → ℕ-Stream (ℕ-Stream A)
ℕ-duplicate₃ xs n m = xs (go n m)
  where
  go : ℕ → ℕ → ℕ
  go zero    m = m
  go (suc n) m = suc (go n m)

ℕ-duplicate₄ : ℕ-Stream A → ℕ-Stream (ℕ-Stream A)
ℕ-duplicate₄ xs n m = xs (n + m)

In other words, duplicate basically adds indices.

There is something distinctly monoidal about what’s going on here: taking the (ℕ, +, 0) monoid as focus, the extract method above corresponds to the monoidal empty element, and the duplicate method corresponds to the binary operator on monoids. In actual fact, there is a comonad for any function from a monoid, often called the Traced comonad.

Traced : Type → Type → Type
Traced E A = E → A

extractᵀ : ⦃ _ : Monoid E ⦄ → Traced E A → A
extractᵀ xs = xs ε

duplicateᵀ : ⦃ _ : Monoid E ⦄ → Traced E A → Traced E (Traced E A)
duplicateᵀ xs e₁ e₂ = xs (e₁ ∙ e₂)

Reifying Traced

The second observation we made about streams was that they had an isomorphic representation which didn’t need coinduction. What we can see above, with Traced, is a representation that also doesn’t need coinduction. So what is the corresponding coinductive representation? What does a generalised reified stream look like?

So the first approach to reifying a function to a data structure is to simply represent the function as a list of pairs.

C-Traced : Type → Type → Type
C-Traced E A = Stream (E × A)

This representation obviously isn’t ideal: it isn’t possible to construct an isomorphism between C-Traced and Traced. We can—kind of—go in one direction, but even that function isn’t terminating:

{-# NON_TERMINATING #-}
lookup-env : ⦃ _ : IsDiscrete E ⦄ → C-Traced E A → Traced E A
lookup-env xs x = if does (x ≟ xs .head .fst)
                     then xs .head .snd
                     else lookup-env (xs .tail) x

I’m not too concerned with being fast and loose with termination and isomorphisms for the time being, though. At the moment, I’m just interested in exploring the relationship between streams and the indexing functions.

As a result, let’s try and push on this representation a little and see if it’s possible to get something interesting and almost isomorphic.

Segmented Streams

To get a slightly nicer representation we can exploit the monoid a little bit. We can do this by storing offsets instead of the absolute indices for each entry. The data structure I have in mind here looks a little like this:

┏━━━━━━━━━━┳━━━━━━┳━━━━━━┉
┃x         ┃y     ┃z     ┉
┡━━━━━━━━━━╇━━━━━━╇━━━━━━┉
╵⇤a╌╌╌╌╌╌╌⇥╵⇤b╌╌╌⇥╵⇤c╌╌╌╌┈

Above is a stream containing the values x, y, and z. Instead of each value corresponding to a single entry in the stream, however, they each correspond to a segment. The value x, for instance, labels the first segment in the stream, which has a length given by a. y labels the second segment, with length b, z with length c, and so on.

The Traced version of the above structure might be something like this:

str :: Traced m a
str i | i < a         = x
      | i < a + b     = y
      | i < a + b + c = z
      | ...

So the index-value mapping is also segmented. The stream, in this way, is kind of like a ruler, where different values mark out different quantities along the ruler, and the index function takes in a quantity and tells you which entry in the ruler that quantity corresponds to.

In code, we might represent the above data structure with the following type:

record Segments (E : Type) (A : Type) : Type where
  field
    length : E
    label  : A
    next   : Segments E A

open Segments

The question is, then, how do we convert this structure to an Traced representation?

Monuses

We need some extra operations on the monoid in the segments in order to enable this conversion to the Traced representation. The extra operations are encapsulated by the monus algebra: I wrote about this in the paper I submitted with Nicolas Wu to ICFP last year (2021). It’s a simple algebra on monoids which basically encapsulates monoids which are ordered in a sensible way.

The basic idea is that we construct an order on monoids which says “x is smaller than y if there is some z that we can add to x to get to y”.

_≼_ : ⦃ _ : Monoid A ⦄ → A → A → Type _
x ≼ y = ∃ z × (y ≡ x ∙ z)

A monus is a monoid where we can extract that z, when it exists. On the monoid (ℕ, +, 0), for instance, this order corresponds to the normal ordering on ℕ.

Extracting the z above corresponds to a kind of difference operator:

_∸_ : ℕ → ℕ → ℕ
x     ∸ zero  = x
suc x ∸ suc y = x ∸ y
_     ∸ _     = zero

This operator is sometimes called the monus. It is a kind of partial, or truncating, subtraction:

_ :  ∸  ≡ 
_ = refl

_ :  ∸  ≡ 
_ = refl

And, indeed, this operator “extracts” the z, when it exists.

∸‿is-monus : ∀ x y → (x≼y : x ≼ y) → y ∸ x ≡ fst x≼y
∸‿is-monus zero    _       (z , y≡0+z) = y≡0+z
∸‿is-monus (suc x) (suc y) (z , y≡x+z) = ∸‿is-monus x y (z , suc-inj y≡x+z)
∸‿is-monus (suc x) zero    (z , 0≡x+z) = ⊥-elim (zero≢suc 0≡x+z)

Our definition of a monus is simple: a monus is anything where the order ≼, sometimes called the “algebraic preorder”, is total and antisymmetric. This is precisely what lets us write a function which takes the Segments type and converts it back to the Traced type.

{-# NON_TERMINATING #-}
Segments→Traced : ⦃ _ : Monus E ⦄ → Segments E A → Traced E A
Segments→Traced xs i with xs .length ≤? i
... | yes (j , i≡xsₗ∙j) = Segments→Traced (xs .next) j
... | no  _             = xs .label

This function takes an index, and checks if that length is greater than or equal to the first segment in the stream of segments. If it is, then it continues searching through the rest of the segments with the index reduced by the size of that first segment. If not, then it returns the label of the first segment.

Taking the old example, we are basically converting to ∸ from +:

str :: Traced m a
str i | i         < a = x
      | i ∸ a     < b = y
      | i ∸ a ∸ b < c = z
      | ...

The first issue here is that this definition is not terminating. That might seem an insurmountable problem at first—we are searching through an infinite stream, after all—but notice that there is one parameter which is decreasing on each recursive call: the index. Well, it only decreases if the segment is non-zero: this can be enforced by changing the definition of the segments type:

record ℱ-Segments (E : Type) ⦃ _ : Monus E ⦄ (A : Type) : Type where
  coinductive
  field
    label    : A
    length   : E
    length≢ε : length ≢ ε
    next     : ℱ-Segments E A

open ℱ-Segments

This type allows us to write the following definition:

module _ ⦃ _ : Monus E ⦄ (wf : WellFounded _≺_) where
  wf-index : ℱ-Segments E A → (i : E) → Acc _≺_ i → A
  wf-index xs i a with xs .length ≤? i
  ... | no _ = xs .label
  wf-index xs i (acc wf) | yes (j , i≡xsₗ∙j) =
    wf-index (xs .next) j (wf j (xs .length , i≡xsₗ∙j ; comm _ _ , xs .length≢ε))

  ℱ-Segments→Traced : ℱ-Segments E A → Traced E A
  ℱ-Segments→Traced xs i = wf-index xs i (wf i)

Trying to build an isomorphism

So the ℱ-Segments type is interesting, but it only really gives one side of the isomorphism. There is no way to write a function Traced E A → ℱ-Segments E A.

The problem is that there’s no way to get the “next” segment from a function E → A. We can find the label of the first segment, by applying the function to ε, but there’s no real way to figure out the size of this segment. We can change Traced little to provide this size, though.

Ind : ∀ E → ⦃ _ : Monus E ⦄ → Type → Type
Ind E A = E → A × Σ[ length ⦂ E ] × (length ≢ ε)

This new type will return a tuple consisting of the value indicated by the supplied index, along with the distance to the next segment. For instance, on the example stream given in the diagram earlier, supplying an index i that is bigger than a but smaller than a + b, this function should return y along with some j such that i + j ≡ a + b. Diagrammatically:

╷⇤i╌╌╌╌╌╌╌╌⇥╷⇤j╌╌⇥╷
┢━━━━━━━━┳━━┷━━━━━╈━━━━━━┉
┃x       ┃y       ┃z     ┉
┡━━━━━━━━╇━━━━━━━━╇━━━━━━┉
╵⇤a╌╌╌╌╌⇥╵⇤b╌╌╌╌╌⇥╵⇤c╌╌╌╌┈

This can be implemented in code like so:

module _ ⦃ _ : Monus E ⦄ where
  wf-ind : ℱ-Segments E A → (i : E) → Acc _≺_ i → A × ∃ length × (length ≢ ε)
  wf-ind xs i _ with xs .length ≤? i
  ... | no xsₗ≰i =
    let j , _ , j≢ε = <⇒≺ i (xs .length) xsₗ≰i
    in xs .label , j , j≢ε
  wf-ind xs i (acc wf) | yes (j , i≡xsₗ∙j) =
    wf-ind (xs .next) j (wf j (xs .length , i≡xsₗ∙j ; comm _ _ , xs .length≢ε))

  ℱ-Segments→Ind : WellFounded _≺_ → ℱ-Segments E A → Ind E A
  ℱ-Segments→Ind wf xs i = wf-ind xs i (wf i)

Again, if the monus has finite descending chains, this function is terminating. And the nice thing about this is that it’s possible to write a function in the other direction:

Ind→ℱ-Segments : ⦃ _ : Monus E ⦄ → Ind E A → ℱ-Segments E A
Ind→ℱ-Segments ind =
  let x , s , s≢ε = ind ε
  in λ where .label    → x
             .length   → s
             .length≢ε → s≢ε
             .next     → Ind→ℱ-Segments (ind ∘ (s ∙_))

The problem here is that this isomorphism is only half correct. We can prove that converting to Ind and back is the identity, but not the other direction. There are too many functions in Ind.

Nonetheless, it’s still interesting!

State Comonad

There is a comonad on state (Waern 2018; Kmett 2018) that is different from store. Notice that above the Ind type has the same type (almost) as State E A.

This is interesting in two ways: first, it gives some concrete, spatial intuition for what’s going on with the state comonad.

Second, it gives a kind of interesting monad instance on the stream. If we apply the Ind→ℱ-Segments function to the implementation of join on state, we should get a join on ℱ-Segments. And we do!

First, we need to redefine Ind to the following:

𝒜-Ind : ∀ E → ⦃ _ : Monus E ⦄ → Type → Type
𝒜-Ind E A = (i : E) → A × Σ[ length ⦂ E ] × (i ≺ length)

This is actually isomorphic to the previous definition, but we return the absolute value of the next segment, rather than the distance to the next segment.

𝒜-iso : ⦃ _ : Monus E ⦄ → 𝒜-Ind E A ⇔ Ind E A
𝒜-iso .fun xs i =
  let x , s , k , s≡i∙k , k≢ε = xs i
  in  x , k , k≢ε
𝒜-iso .inv xs i =
  let x , s , s≢ε = xs i
  in  x , i ∙ s , s , refl , s≢ε
𝒜-iso .rightInv _ = refl
𝒜-iso .leftInv  xs p i =
  let x , s           , k , s≡i∙k                   , k≢ε = xs i
  in  x , s≡i∙k (~ p) , k , (λ q → s≡i∙k (~ p ∨ q)) , k≢ε

The implementation of join on this type is the following:

𝒜-join : ⦃ _ : Monus E ⦄ → 𝒜-Ind E (𝒜-Ind E A) → 𝒜-Ind E A
𝒜-join xs i =
  let x , j , i<j = xs i
      y , k , k<j = x j
  in  y , k , ≺-trans i<j k<j

This is the same definition of join as for State, modulo the < fiddling.

On a stream, this operation corresponds to taking a stream of streams and collapsing it to a single stream. It does this by taking a prefix of each internal stream equal in size to the segment of the outer entry. Diagrammatically:

┏━━━━━━━━━━┳━━━━━━┳━━━━━━┉
┃xs        ┃ys    ┃zs    ┉
┡━━━━━━━━━━╇━━━━━━╇━━━━━━┉
╵⇤a╌╌╌╌╌╌╌⇥╵⇤b╌╌╌⇥╵⇤c╌╌╌╌┈
          ╱        ╲
         ╱          ╲
        ╱            ╲
       ╱              ╲
      ╱                ╲
     ╷⇤b╌╌╌╌╌╌╌╌╌╌╌╌╌╌╌⇥╷
     ┢━━━━━━━┳━━━━━━┳━━━┷┉
ys = ┃xʸ     ┃yʸ    ┃zʸ  ┉
     ┡━━━━━━━╇━━━━━━╇━━━━┉
     ╵⇤aʸ╌╌╌⇥╵⇤bʸ╌╌⇥╵⇤cʸ╌┈

Here we start with a stream consisting of the streams xs, ys, and zs, followed by some other streams. Zooming in on ys, we see that it is in a segment of length b, and consists of three values xʸ, yʸ, and zʸ, with segment lengths aʸ, bʸ, and cʸ, respectively.

Calling join on this stream will give us the following stream:

┏━┉━┳━━━━┳━━━━┳━━━━━┳━━━━┉
┃ ┉ ┃xʸ  ┃yʸ  ┃zʸ   ┃    ┉
┡━┉━╇━━━━╇━━━━╇━━━━━╇━━━━┉
│   │⇤aʸ⇥╵⇤bʸ⇥╵⇤╌╌┈⇥│
╵⇤a⇥╵⇤b╌╌╌╌╌╌╌╌╌╌╌╌⇥╵⇤c╌╌┈

Again, we’re focusing on the ys section here, which occupies the segment from a to a ∙ b. After join, this segment is occupied by three elements, xʸ, yʸ, and zʸ.

Notice that this isn’t quite the normal join on streams. That join takes a stream of streams, and turns the ith entry into the ith entry in the underlying stream. It’s a diagonalisation, in other words.

This one is kind of similar, but it takes chunks of the outer stream.

Theory

All of this so far is very hand-wavy. We have an almost isomorphism (a split surjection, to be precise), but not much in the way of concrete theoretical insights, just some vague gesturing towards spatial metaphors and so on.

Thankfully, there are two separate areas of more serious research that seem related to the stuff I’ve talked about here. The first is update monads and directed containers, and the second is graded comonads. I think I understand graded comonads and the related work better out of the two, but update monads and directed containers seems more closely related to what I’m doing here.

Update Monads and Directed Containers

There are a few papers on this topic: Ahman, Chapman, and Uustalu (2012), Ahman and Uustalu (2013; Ahman and Uustalu 2014; Ahman and Uustalu 2016).

The first of these, “When Is a Container a Comonad?” constructs, as the title suggests, a class for containers which are comonads in a standard way.

Here’s the definition of a container:

Container : Type₁
Container = Σ[ Shape ⦂ Type ] × (Shape → Type)

⟦_⟧ : Container → Type → Type
⟦ S , P ⟧ X = Σ[ s ⦂ S ] × (P s → X)

Containers are a generic way to describe a class of well-behaved functors. Any container is a pair of a shape and position. Lists, for instance, are containers, where their shape is described by the natural numbers (the shape here is the length of the list). The positions in such a list are the numbers smaller than the length, in dependently-typed programming we usually use the Fin type for this:

Fin : ℕ → Type
Fin n = ∃ m × (m <ℕ n)

The container version of lists, then, is the following:

ℒ𝒾𝓈𝓉 : Type → Type
ℒ𝒾𝓈𝓉 = ⟦ ℕ , Fin ⟧

Here’s the same list represented in the standard way, and as a container:

someBools : List Bool
someBools = true ∷ true ∷ false ∷ true ∷ []

someBools′ : ℒ𝒾𝓈𝓉 Bool
someBools′ .fst = 
someBools′ .snd ( , _) = true
someBools′ .snd ( , _) = true
someBools′ .snd ( , _) = false
someBools′ .snd ( , _) = true

The benefit of using containers is that it gives a standard, generic, and composable way to construct functors that have some nice properties (like strict positivity). They’re pretty annoying to use in practice, though, which is a shame.

Directed containers are container that have three extra operations.

• A tail-like operation, where a position can be converted into the shape of containers that the suffix from that position.
• A head-like operation, where you can always return the root position.
• A +-like operation, where you take a position on some tail and translate it into a position on the original container, by adding it.

As the paper observes, these are very similar to a “dependently-typed” version of the monoid methods. This seems to me to be very similar to the indexing stuff we were doing earlier on.

The real interesting part is in the paper “Update Monads: Cointerpreting Directed Containers” (Ahman and Uustalu 2014). This paper presents a variant on state monads, called “update monads”.

These are monads that use a monoid action:

record RightAction (𝑃 : Type) (𝑆 : Type) : Type where
  infixl  _↓_
  field
    ⦃ monoid⟨𝑃⟩ ⦄ : Monoid 𝑃
    _↓_ : 𝑆 → 𝑃 → 𝑆
    ↓-assoc : ∀ x y z → (x ↓ y) ↓ z ≡ x ↓ (y ∙ z)
    ↓-ε : ∀ x → x ↓ ε ≡ x

A (right) monoid action is a monoid along with a function ↓ that “acts” on some other set, in a way that coheres with the monoid methods. The definition is given above. One way to think about it is that if a monoid 𝑃 has an action on 𝑆 it means that elements of 𝑃 can kind of be transformed into elements of 𝑆 → 𝑆.

Upd : (𝑃 𝑆 : Type) ⦃ _ : RightAction 𝑃 𝑆 ⦄ → Type → Type
Upd 𝑃 𝑆 X = 𝑆 → 𝑃 × X

η : ⦃ _ : RightAction 𝑃 𝑆 ⦄ → A → Upd 𝑃 𝑆 A
η x s = ε , x

μ : ⦃ _ : RightAction 𝑃 𝑆 ⦄ → Upd 𝑃 𝑆 (Upd 𝑃 𝑆 A) → Upd 𝑃 𝑆 A
μ xs s = let p , x = xs s
             q , y = x (s ↓ p)
         in  (p ∙ q , y)

It turns out that the dependently-typed version of this gives directed containers.

Grading and the Cofree Comonad

I’m still in the early stages of understanding all of this material, but at the moment graded comonads and transformers are concepts that I’m much more familiar and comfortable with.

The idea behind graded monads and comonads is similar to the idea behind any indexed monad: we’re adding an extra type parameter to the monad or type, which can constrain the operations involved. The graded monads and comonads use a monoid as that index. This works particularly nicely, in my opinion: just allowing any index at all sometimes feels a little unstructured. The grading construction seems to constrain things to the right degree: the use of the monoid, as well, works really well with comonads.

That preamble out of the way, here’s the definition of a graded comonad:

record GradedComonad (𝑆 : Type) ⦃ _ : Monoid 𝑆 ⦄ (𝐶 : 𝑆 → Type → Type) : Type₁ where
  field
    extract : 𝐶 ε A → A
    extend  : (𝐶 y A → B) → 𝐶 (x ∙ y) A → 𝐶 x B

  _=<=_ : (𝐶 x B → C) → (𝐶 y A → B) → 𝐶 (x ∙ y) A → C
  (g =<= f) x = g (extend f x)

  field
    idˡ : (f : 𝐶 x A₀ → B₀) → PathP (λ i → 𝐶 (ε∙ x i) A₀ → B₀) (extract =<= f) f
    idʳ : (f : 𝐶 x A₀ → B₀) → PathP (λ i → 𝐶 (∙ε x i) A₀ → B₀) (f =<= extract) f
    c-assoc : (f : 𝐶 x C₀ → D₀) (g : 𝐶 y B₀ → C₀) (h : 𝐶 z A₀ → B₀) →
          PathP (λ i → 𝐶 (assoc x y z i) A₀ → D₀) ((f =<= g) =<= h) (f =<= (g =<= h))

This seems to clearly be related to the stream constructions. Grading is all about the monoidal information about a comonad: the streams above are a comonad which indexes its entries with a monoid.

There are now two constructions I want to show that suggest a link between the stream constructions and graded comonads. First of these is the Cofree degrading comonad:

record G-CofreeF (𝐹 : Type → Type) (𝐶 : 𝑆 → Type → Type) (A : Type) : Type where
  coinductive; constructor _◃_
  field here : A
        step : 𝐹 (∃ w × 𝐶 w (G-CofreeF 𝐹 𝐶 A))
open G-CofreeF

G-Cofree : ⦃ _ : Monoid 𝑆 ⦄ → (Type → Type) → (𝑆 → Type → Type) → Type → Type
G-Cofree 𝐹 𝐶 A = 𝐶 ε (G-CofreeF 𝐹 𝐶 A)

This construction is similar to the cofree comonad transformer: it is based on the cofree comonad, but with an extra (graded) comonad wrapped around each level. For any functor 𝐹 and graded comonad 𝐶, G-Cofree 𝐹 𝐶 is a comonad. The implementation of extract is simple:

extract′ : ⦃ _ : Monoid 𝑆 ⦄ ⦃ _ : GradedComonad 𝑆 𝐶 ⦄ → G-Cofree 𝐹 𝐶 A → A
extract′ = here ∘ extract

extend is more complex. First, we need a version of extend which takes a proof that the grade is of the right form:

module _ { 𝐶 : 𝑆 → Type → Type } where
  extend[_] : ⦃ _ : Monoid 𝑆 ⦄ ⦃ _ : GradedComonad 𝑆 𝐶 ⦄ →
              x ∙ y ≡ z → (𝐶 y A → B) → 𝐶 z A → 𝐶 x B
  extend[ p ] k = subst (λ z → 𝐶 z _ → _) p (extend k)

Then we can implement the characteristic function on the free comonad: traceT. On graded comonads it has the following form:

module Trace ⦃ _ : Monoid 𝑆 ⦄ ⦃ _ : GradedComonad 𝑆 𝐶 ⦄ ⦃ _ : Functor 𝐹 ⦄ where
  module _ {A B} where
    {-# NON_TERMINATING #-}
    traceT : (𝐶 ε A → B) → (𝐶 ε A → 𝐹 (∃ w × 𝐶 w A)) → 𝐶 ε A → G-Cofree 𝐹 𝐶 B
    traceT ϕ ρ = ψ
      where
      ψ : 𝐶 x A → 𝐶 x (G-CofreeF 𝐹 𝐶 B)
      ψ = extend[ ∙ε _ ] λ x → ϕ x ◃ map (map₂ ψ) (ρ x)

This function is basically the unfold for the free degrading comonad. If G-Cofree is a internally-labelled tree, then ϕ above is the labelling function, and ρ is the “next” function, returning the children for some root.

Using this, we can implement extend:

  extend′ : (G-Cofree 𝐹 𝐶 A → B) → G-Cofree 𝐹 𝐶 A → G-Cofree 𝐹 𝐶 B
  extend′ f = traceT f (step ∘ extract)

The relation between this and the stream is that the stream can be defined in terms of this: Stream W = G-Cofree id (GC-Id W).

Finally, the last construction I want to introduce is the following:

module _ ⦃ _ : Monus 𝑆 ⦄ where
  data Prefix-F⊙ (𝐹 : Type → Type) (𝐶 : 𝑆 → Type → Type) (i j : 𝑆) (A : Type) : Type where
    prefix : ((i≤j : i ≤ j) → A × 𝐹 (∃ k × 𝐶 k (Prefix-F⊙ 𝐹 𝐶 k (fst i≤j) A))) → Prefix-F⊙ 𝐹 𝐶 i j A

  Prefix⊙ : (𝐹 : Type → Type) (𝐶 : 𝑆 → Type → Type) (j : 𝑆) (A : Type) → Type
  Prefix⊙ 𝐹 𝐶 j A = 𝐶 ε (Prefix-F⊙ 𝐹 𝐶 ε j A)

  Prefix : (𝐹 : Type → Type) (𝐶 : 𝑆 → Type → Type) (A : Type) → Type
  Prefix 𝐹 𝐶 A = ∀ {i} → Prefix⊙ 𝐹 𝐶 i A

This type is designed to mimic sized type definitions. It has an implicit parameter which can be set, by the user of the type, to some arbitrary depth. Basically the parameter means “explore to this depth”; by using the ∀ we say that it is defined up to any arbitrary depth.

When the ≺ relation on the monus is well founded it is possible to implement traceT:

  module _ ⦃ _ : GradedComonad 𝑆 𝐶 ⦄ ⦃ _ : Functor 𝐹 ⦄ (wf : WellFounded _≺_) {A B : Type} where
    traceT : (𝐶 ε A → B) → (𝐶 ε A → 𝐹 (∃ w × (w ≢ ε) × 𝐶 w A)) → 𝐶 ε A → Prefix 𝐹 𝐶 B
    traceT ϕ ρ xs = extend[ ∙ε _ ] (λ xs′ → prefix λ _ → ϕ xs′ ,  map (map₂ (ψ (wf _))) (ρ xs)) xs
      where
      ψ : Acc _≺_ y → (x ≢ ε) × 𝐶 x A → 𝐶 x (Prefix-F⊙ 𝐹 𝐶 x y B)
      ψ (acc wf) (x≢ε , xs) =
        extend[ ∙ε _ ]
          (λ x → prefix
            λ { (k , y≡x∙k) →
              ϕ x , map
                (λ { (w , w≢ε , xs) →
                  w , ψ (wf k (_ , y≡x∙k ; comm _ _ , x≢ε)) (w≢ε , xs)}) (ρ x)})
          xs

Conclusion

Comonads are much less widely used than monads in Haskell and similar languages. Part of the reason, I think, is that they’re too powerful in a non-linear language. Monads are often used to model sublanguages where it’s possible to introduce “special” variables which interact with the monadic context.

pyth = do
  x <- [1..10]
  y <- [1..10]
  z <- [1..10]
  guard (x*x + y*y == z*z)
  return (x,y,z)

The x variable here semantically spans over the range [1..10]. In the following two examples we see the semantics of state and maybe:

sum :: [Int] -> Int
sum xs = flip evalState 0 $ do
  put 0
  for_ xs $ \x -> do
    n <- get
    put (n + x)
  m <- get
  return m

data E = Lit Int | E :+: E | E :/: E

eval :: E -> Maybe Int
eval (Lit n) = n
eval (xs :+: ys) = do x <- eval xs
                      y <- eval ys
                      return (x + y)
eval (xs :/: ys) = do x <- eval xs
                      y <- eval ys
                      guard (y /= 0)
                      return (x / y)

The variables n and m introduced in the state example are “special” because their values depend on the computations that came before. In the maybe example the variables introduced could be Nothing.

You can’t do the same thing with comonads because you’re always able to extract the “special” variable with extract :: m a -> a. Instead of having special variable introduction, comonads let you have special variable elimination. But, since Haskell isn’t linear, you can always just discard a variable so this isn’t much use.

Looking at the maybe example, we have a function eval :: E -> Maybe Int that introduces an Int variable with a “catch”: it is wrapped in a Maybe. We want to use the eval function as if it were a normal function E -> Int, with all of the bookkeeping managed for us: that’s what monads and do notation (kind of) allow us to do.

An analagous example with comonads might be having a function consume :: m V -> String. This “handles” a V value, but the “catch” is that it needs an m context to do so. If we want to treat the consume function as if it were a normal function V -> String then comonads (and codo notation Orchard and Mycroft 2013) would be a perfect fit.

The reason that this analagous case doesn’t arise very often is that we don’t have many handlers that look like m V -> String in Haskell. Why? Because if we want to “handle” a V we can just discard it: as a non-linear language, you do not need to perform any ceremony to discard a variable in Haskell.

Graded comonads, though, seem to be much more useful than normal comonads. I think it is becuase they basically get rid of the m a -> a function, changing it into a much more restricted form. In this way, they give a kind of small linear language, but just for the monoidal type parameter.

And there are a lot of uses for the graded comonads. Above we’ve used them for termination checking. A recursive function might have the form a -> b, where a is the thing being recursed on. If we’re using well-founded recursion to show that it’s terminating, though, we add an extra parameter, an Acc _<_ proof, turning this function into Acc _<_ w × a -> b. The Acc _<_ here is the graded comonad, and this recursive function is precisely the “handler”.

Other examples might be privacy or permissions: a function might be able to work on some value, but only if it has particular permission regarding that value. The permission here is the monoid.

There are other examples I’m sure, those are just the couple that I have been thinking about.

References

I have packaged up the more interesting bits from the Algebras for Weighted Search paper and put it up on hackage.

You can see it here.

It contains the HeapT monad, the Monus class, and an implementation of Dijkstra’s algorithm, the Viterbi algorithm, and probabilistic parsing.

Check it out!

The paper “Algebras for Weighted Search” has just been accepted unconditionally to ICFP. I wrote it with my supervisor, Nicolas Wu, and it covers a lot of the topics I’ve written about on this blog (including hyperfunctions and breadth-first traversals).

The preprint is available here.

Check out this type:

newtype a -&> b = Hyp { invoke :: (b -&> a) -> b }

This a hyperfunction (J. Launchbury, Krstic, and Sauerwein 2013; 2000; 2000), and I think it’s one of the weirdest and most interesting newtypes you can write in Haskell.

The first thing to notice is that the recursion pattern is weird. For a type to refer to itself recursively on the left of a function arrow is pretty unusual, but on top of that the recursion is non-regular. That means that the recursive reference has different type parameters to its parent: a -&> b is on the left-hand-side of the equals sign, but on the right we refer to b -&> a.

Being weird is reason enough to write about them, but what’s really shocking about hyperfunctions is that they’re useful. Once I saw the definition I realised that a bunch of optimisation code I had written (to fuse away zips in particular) was actually using hyperfunctions (Ghani et al. 2005). After that, I saw them all over the place: in coroutine implementations, queues, breadth-first traversals, etc.

Anyways, since coming across hyperfunctions a few months ago I thought I’d do a writeup on them. I’m kind of surprised they’re not more well-known, to be honest: they’re like a slightly more enigmatic Cont monad, with a far cooler name. Let’s get into it!

What Are Hyperfunctions?

The newtype noise kind of hides what’s going on with hyperfunctions: expanding the definition out might make things slightly clearer.

type a -&> b = (b -&> a) -> b
             = ((a -&> b) -> a) -> b
             = (((b -&> a) -> b) -> a) -> b
             = ((((... -> b) -> a) -> b) -> a) -> b

So a value of type a -&> b is kind of an infinitely left-nested function type. One thing worth noticing is that all the as are in negative positions and all the bs in positive. This negative and positive business basically refers to the position of arguments in relation to a function arrow: to the left are negatives, and to the right are positives, but two negatives cancel out.

((((... -> b) -> a) -> b) -> a) -> b
           +     -     +     -     +

All the things in negative positions are kind of like the things a function “consumes”, and positive positions are the things “produced”. It’s worth fiddling around with very nested function types to get a feel for this notion. For hyperfunctions, though, it’s enough to know that a -&> b does indeed (kind of) take in a bunch of as, and it kind of produces bs.

By the way, one of the ways to get to grips with polarity in this sense is to play around with the Cont monad, codensity monad, or selection monad (Hedges 2015). If you do, you may notice one of the interesting parallels about hyperfunctions: the type a -&> a is in fact the fixpoint of the continuation monad (Fix (Cont a)). Suspicious!

Hyperfunctions Are Everywhere

Before diving further into the properties of the type itself, I’d like to give some examples of how it can show up in pretty standard optimisation code.

Zips

Let’s say you wanted to write zip with foldr (I have already described this particular algorithm in a previous post). Not foldr on the left argument, mind you, but foldr on both. If you proceed mechanically, replacing every recursive function with foldr, you can actually arrive at a definition:

zip :: [a] -> [b] -> [(a,b)]
zip xs ys = foldr xf xb xs (foldr yf yb ys)
  where
    xf x xk yk = yk x xk
    xb _ = []

    yf y yk x xk = (x,y) : xk yk
    yb _ _ = []

In an untyped language, or a language with recursive types, such a definition would be totally fine. In Haskell, though, the compiler will complain with the following:

• Occurs check: cannot construct the infinite type:
    t0 ~ a -> (t0 -> [(a, b)]) -> [(a, b)]

Seasoned Haskellers will know, though, that this is not a type error: no, this is a type recipe. The compiler is telling you what parameters it wants you to stick in the newtype:

newtype Zip a b =
  Zip { runZip :: a -> (Zip a b -> [(a,b)]) -> [(a,b)] }

zip :: forall a b. [a] -> [b] -> [(a,b)]
zip xs ys = xz yz
  where
    xz :: Zip a b -> [(a,b)]
    xz = foldr f b xs
      where
        f x xk yk = runZip yk x xk
        b _ = []

    yz :: Zip a b
    yz = foldr f b ys
      where
        f y yk = Zip (\x xk -> (x,y) : xk yk)
        b = Zip (\_ _ -> [])

And here we see the elusive hyperfunction: hidden behind a slight change of parameter order, Zip a b is in fact the same as [(a,b)] -&> (a -> [(a,b)]).

zip :: forall a b. [a] -> [b] -> [(a,b)]
zip xs ys = invoke xz yz
  where
    xz :: (a -> [(a,b)]) -&> [(a,b)]
    xz = foldr f b xs
      where
        f x xk = Hyp (\yk -> invoke yk xk x)
        b = Hyp (\_ -> [])

    yz :: [(a,b)] -&> (a -> [(a,b)])
    yz = foldr f b ys
      where
        f y yk = Hyp (\xk x -> (x,y) : invoke xk yk)
        b = Hyp (\_ _ -> [])

BFTs

In another previous post I derived the following function to do a breadth-first traversal of a tree:

data Tree a = a :& [Tree a]

newtype Q a = Q { q :: (Q a -> [a]) -> [a] }

bfe :: Tree a -> [a]
bfe t = q (f t b) e
  where
    f :: Tree a -> Q a -> Q a
    f (x :& xs) fw = Q (\bw -> x : q fw (bw . flip (foldr f) xs))

    b :: Q a
    b = Q (\k -> k b)

    e :: Q a -> [a]
    e (Q q) = q e

That Q type there is another hyperfunction.

bfe :: Tree a -> [a]
bfe t = invoke (f t e) e
  where
    f :: Tree a -> ([a] -&> [a]) -> ([a] -&> [a])
    f (x :& xs) fw = Hyp (\bw -> x : invoke fw (Hyp (invoke bw . flip (foldr f) xs)))

    e :: [a] -&> [a]
    e = Hyp (\k -> invoke k e)

One of the problems I had with the above function was that it didn’t terminate: it could enumerate all the elements of the tree but it didn’t know when to stop. A similar program (Allison 2006; described and translated to Haskell in Smith 2009) manages to solve the problem with a counter. Will it shock you to find out this solution can also be encoded with a hyperfunction?

bfe t = invoke (f t (Hyp b)) e 1
  where
    f :: Tree a -> (Int -> [a]) -&> (Int -> [a])
                -> (Int -> [a]) -&> (Int -> [a])
    f (x :& xs) fw =
      Hyp (\bw n -> x : invoke fw (Hyp (\k m -> invoke bw (foldr f k xs) (m+1))) n)

    e :: (Int -> [a]) -&> (Int -> [a])
    e = Hyp (\k -> invoke k e)

    b x 0 = []
    b x n = invoke x (Hyp b) (n-1)

(my version here is actually a good bit different from the one in Smith 2009, but the basic idea is the same)

Coroutines

Hyperfunctions seem to me to be quite deeply related to coroutines. At the very least several of the types involved in coroutine implementations are actual hyperfunctions. The ProdPar and ConsPar types from Pieters and Schrijvers (2019) are good examples:

newtype ProdPar a b = ProdPar (ConsPar a b -> b)
newtype ConsPar a b = ConsPar (a -> ProdPar a b -> b)

ProdPar a b is isomorphic to (a -> b) -&> b, and ConsPar a b to b -&> (a -> b), as witnessed by the following functions:

fromP :: ProdPar a b -> (a -> b) -&> b
fromP (ProdPar x) = Hyp (x . toC)

toC ::  b -&> (a -> b) -> ConsPar a b
toC (Hyp h) = ConsPar (\x p -> h (fromP p) x)

toP :: (a -> b) -&> b -> ProdPar a b
toP (Hyp x) = ProdPar (x . fromC)

fromC :: ConsPar a b -> b -&> (a -> b)
fromC (ConsPar p) = Hyp (\h x -> p x (toP h))

In fact this reveals a little about what was happening in the zip function: we convert the left-hand list to a ProdPar (producer), and the right-hand to a consumer, and apply them to each other.

Hyperfunctions Are Weird

Aside from just being kind of weird intuitively, hyperfunctions are weird in theory. Set-theoretically, for instance, you cannot form the set of a -&> b: if you tried, you’d run into those pesky size restrictions which stop us from making things like “the set of all sets”. Haskell types, however, are not sets, precisely because we can define things like a -&> b.

For slightly different reasons to the set theory restrictions, we can’t define the type of hyperfunctions in Agda. The following will get an error:

record _↬_ (A : Type a) (B : Type b) : Type (a ℓ⊔ b) where
  inductive; constructor hyp
  field invoke : (B ↬ A) → B

And for good reason! Agda doesn’t allow recursive types where the recursive call is in a negative position. If we turn off the positivity checker, we can write Curry’s paradox (example proof taken from here):

yes? : ⊥ ↬ ⊥
yes? .invoke h = h .invoke h

no! : (⊥ ↬ ⊥) → ⊥
no! h = h .invoke h

boom : ⊥
boom = no! yes?

Note that this isn’t an issue with the termination checker: the above example passes all the normal termination conditions without issue (yes, even if ↬ is marked as coinductive). It’s directly because the type itself is not positive.

Interestingly, there is a slightly different, and nearly equivalent, definition of hyperfunctions which doesn’t allow us to write the above proof:

record _↬_ (A : Type a) (B : Type b) : Type (a ℓ⊔ b) where
  inductive; constructor hyp
  field invoke : ((A ↬ B) → A) → B

This is basically a slightly expanded out version of the hyperfunction type, and importantly it’s positive. Not strictly positive however, since the recursive call does occur to the left of a function arrow: it’s just positive, in that it’s to the left of an even number of function arrows.

I found in a blog post by Sjöberg (2015) some interesting discussion regarding the question of this extra strictness: in Coq, allowing certain positive but not strictly positive types does indeed introduce an inconsistency (Coquand and Paulin 1990). However this inconsistency relies on an impredicative universe, which Agda doesn’t have. As far as I understand it, it would likely be safe to allow types like ↬ above in Agda (Coquand 2013), although I’m not certain that with all of Agda’s newer features that’s still the case.

The connection between non-strictly-positive types and breadth-first traversals has been noticed before: Berger, Matthes, and Setzer (2019) make the argument for their inclusion in Agda and Coq using a breadth-first traversal algorithm by Hofmann (1993), which uses the following type:

data Rou
  = Over
  | Next ((Rou -> [Int]) -> [Int])

Now this type isn’t a hyperfunction (but it’s close); we’ll see soon what kind of thing it is.

Hyperfunctions Are a Category

So we’ve seen that hyperfunctions show up kind of incidentally through certain optimisations, and we’ve seen that they occupy a strange space in terms of their theoretical interpretation: we haven’t yet seen much about the type itself in isolation. Luckily Ed Kmett has already written the hyperfunctions package (2015), where a laundry list of instances are provided, which can tell us a little more about what hyperfunctions can actually do on their own.

The Category instance gives us the following:

instance Category (-&>) where
  id = Hyp (\k -> invoke k id)
  f . g = Hyp (\k -> invoke f (g . k))

We’ve actually seen the identity function a few times: we used it as the base case for recursion in the breadth-first traversal algorithms.

Composition we actually have used as well but it’s more obscured. An analogy to help clear things up is to think of hyperfunctions as a kind of stack. id is the empty stack, and we can use the following function to push items onto the stack:

push :: (a -> b) -> a -&> b -> a -&> b
push f q = Hyp (\k -> f (invoke k q))

Understood in this sense, composition acts like a zipping operation on stacks, since we have the following law:

push f p . push g q ≡ push (f . g) (p . q)

While we can’t really pop elements off the top of the stack directly, we can get close with invoke, since it satisfies the following law:

invoke (push f p) q ≡ f (invoke q p)

Along with the id implementation we have, this will let us run a hyperfunction, basically folding over the contents of the stack:

run :: a -&> a -> a
run f = invoke f id

This analogy helps us understand how the breadth-first traversals worked: the hyperfunctions are kind of like stacks with $\mathcal{O}(1)$ push and zip, which is precisely what you need for an efficient breadth-first traversal.

bfe :: Tree a -> [a]
bfe = run . f
  where
    f (x :& xs) = push (x:) (zips (map f xs))

    zips = foldr (.) id

Finally, hyperfunctions are of course monads:

instance Monad ((-&>) a) where
  m >>= f = Hyp (\k -> invoke (f (invoke m (Hyp (invoke k . (>>=f))))) k)

I won’t pretend to understand what’s going on here, but it looks a little like a nested reader monad. Perhaps there’s some intuition to be gained from noticing that a -&> a ~ Fix (Cont a).

Hyper Arrows Are…?

As I said in the introduction I’m kind of surprised there’s not more research out there on hyperfunctions. Aside from the excellent papers by J. Launchbury, Krstic, and Sauerwein (2013) there’s just not much out there. Maybe it’s that there’s not that much theoretical depth to them, but all the same there are some clear questions worth looking into.

For example: is there a hyperfunction monad transformer? Or, failing that, can you thread a monad through the type at any point, and do you get anything interesting out?

I have made a little headway on this question, while fiddling with one of the bfe definitions above. Basically I wanted to remove the Int counter for the terminating bfe, and I wanted to use a Maybe somewhere instead. I ended up generalising from Maybe to any m, yielding the following type:

newtype HypM m a b = HypM { invokeM :: m ((HypM m a b -> a) -> b) }

This does the job for the breadth-first traversal:

bfe t = r (f t e)
  where
    f :: Tree a -> HypM Maybe [a] [a] -> HypM Maybe [a] [a]
    f (x :& xs) fw = HypM (Just (\bw -> x : fromMaybe (\k -> k e) (invokeM fw) (bw . flip (foldr f) xs)))

    e :: HypM Maybe [a] [a]
    e = HypM Nothing

    r :: HypM Maybe [a] [a] -> [a]
    r = maybe [] (\k -> k r) . invokeM

(In fact, when m is specialised to Maybe we have the same type as Rou)

This type has a very practical use, as it happens, which is related to the church-encoded list monad transformer:

newtype ListT m a = ListT { runListT :: forall b. (a -> m b -> m b) -> m b -> m b }

Just like -&> allowed us to write zip on folds (i.e. using foldr), HypM will allow us to write zipM on ListT:

zipM :: Monad m => ListT m a -> ListT m b -> ListT m (a,b)
zipM xs ys = ListT (\c n ->
  let
    xf x xk = pure (\yk -> yk (HypM xk) x)
    xb = pure (\_ -> n)

    yf y yk = pure (\xk x -> c (x, y) (join (invokeM xk <*> yk)))
    yb = pure (\_ _ -> n)
  in join (runListT xs xf xb <*> runListT ys yf yb))

I actually think this function could be used to seriously improve the running time of several of the functions on LogicT: my reading of them suggests that interleave is $\mathcal{O}(n^2)$ (or worse), but the zip above could be trivially repurposed to give a $\mathcal{O}(n)$ interleave. This would also have knock-on effects on, for instance, >>- and so on.

Another question is regarding the arrows of the hyperfunction. We’ve seen that a hyperfunction kind of adds “stacking” to functions, can it do the same for other arrows? Basically, does the following type do anything useful?

newtype HypP p a b = HypP { invokeP :: p (HypP p b a) b }

Along a similar vein, many of the breadth-first enumeration algorithms seem to use “hyperfunctions over the endomorphism monoid”. Basically, they all produce hyperfunctions of the type [a] -&> [a], and use them quite similarly to how we would use difference lists. But we know that there are Cayley transforms in other monoidal categories, for instance in the applicative monoidal category: can we construct the “hyperfunction” version of those?

References

The final version of my master’s thesis got approved recently so I thought I’d post it here for people who might be interested.

Here’s the university record.

Here’s the pdf.

And all of the theorems in the thesis have been formalised in Agda. The code is organised to follow the structure of the pdf here.

The title of the thesis is “Finiteness in Cubical Type Theory”: basically it’s all about formalising the notion of “this type is finite” in CuTT. I also wanted to write something that could serve as a kind of introduction to some components of modern dependent type theory which didn’t go the standard length-indexed vector route.

The Cayley monoid is well-known in Haskell (difference lists, for instance, are a specific instance of the Cayley monoid), because it gives us $O(1)$ <>. What’s less well known is that it’s also important in dependently typed programming, because it gives us definitional associativity. In other words, the type x . (y . z) is definitionally equal to (x . y) . z in the Cayley monoid.

data Nat = Z | S Nat

type family (+) (n :: Nat) (m :: Nat) :: Nat where
  Z   + m = m
  S n + m = S (n + m)

I used a form of the type-level Cayley monoid in a previous post to type vector reverse without proofs. I figured out the other day another way to use it to type tree flattening.

Say we have a size-indexed tree and vector:

data Tree (a :: Type) (n :: Nat) :: Type where
  Leaf  :: a -> Tree a (S Z)
  (:*:) :: Tree a n -> Tree a m -> Tree a (n + m)

data Vec (a :: Type) (n :: Nat) :: Type where
  Nil  :: Vec a Z
  (:-) :: a -> Vec a n -> Vec a (S n)

And we want to flatten it to a list in $O(n)$ time:

treeToList :: Tree a n -> Vec a n
treeToList xs = go xs Nil
  where
    go :: Tree a n -> Vec a m -> Vec a (n + m)
    go (Leaf    x) ks = x :- ks
    go (xs :*: ys) ks = go xs (go ys ks)

Haskell would complain specifically that you hadn’t proven the monoid laws:

• Couldn't match type ‘n’ with ‘n + 'Z’
• Could not deduce: (n2 + (m1 + m)) ~ ((n2 + m1) + m)

But it seems difficult at first to figure out how we can apply the same trick as we used for vector reverse: there’s no real way for the Tree type to hold a function from Nat to Nat.

To solve this problem we can borrow a trick that Haskellers had to use in the good old days before type families to represent type-level functions: types (or more usually classes) with multiple parameters.

data Tree' (a :: Type) (n :: Nat) (m :: Nat) :: Type where
  Leaf  :: a -> Tree' a n (S n)
  (:*:) :: Tree' a n2 n3
        -> Tree' a n1 n2
        -> Tree' a n1 n3

The Tree' type here has three parameters: we’re interested in the last two. The first of these is actually an argument to a function in disguise; the second is its result. To make it back into a normal size-indexed tree, we apply that function to zero:

type Tree a = Tree' a Z

three :: Tree Int (S (S (S Z)))
three = (Leaf 1 :*: Leaf 2) :*: Leaf 3

This makes the treeToList function typecheck without complaint:

treeToList :: Tree a n -> Vec a n
treeToList xs = go xs Nil
  where
    go :: Tree' a x y -> Vec a x -> Vec a y
    go (Leaf    x) ks = x :- ks
    go (xs :*: ys) ks = go xs (go ys ks)

Consider the following puzzle:

Given a list of $n$ labels, list all the trees with those labels in order.

For instance, given the labels [1,2,3,4], the answer (for binary trees) is the following:

┌1     ┌1      ┌1     ┌1     ┌1
┤      ┤      ┌┤     ┌┤     ┌┤
│┌2    │ ┌2   ││┌2   │└2    │└2
└┤     │┌┤    │└┤    ┤     ┌┤
 │┌3   ││└3   │ └3   │┌3   │└3
 └┤    └┤     ┤      └┤    ┤
  └4    └4    └4      └4   └4

This problem (the “enumeration” problem) turns out to be quite fascinating and deep, with connections to parsing and monoids. It’s also just a classic algorithmic problem which is fun to try and solve.

The most general version of the algorithm is on forests of rose trees:

data Rose a = a :& Forest a
type Forest a = [Rose a]

It’s worth having a go at attempting it yourself, but if you’d just like to see the slick solutions the following is one I’m especially proud of:

enumForests :: [a] -> [Forest a]
enumForests = foldrM f []
  where
    f x xs = zipWith ((:) . (:&) x) (inits xs) (tails xs)

In the rest of this post I’ll go through the intuition behind solutions like the one above and I’ll try to elucidate some of the connections to other areas of computer science.

A First Approach: Trying to Enumerate Directly

I first came across the enumeration problem when I was writing my master’s thesis: I needed to prove (in Agda) that there were finitely many binary trees of a given size, and that I could list them (this proof was part of a larger verified solver for the countdown problem). My first few attempts were unsuccessful: the algorithm presented in the countdown paper (Hutton 2002) was not structurally recursive, and did not seem amenable to Agda-style proofs.

Instead, I looked for a type which was isomorphic to binary trees, and which might be easier to reason about. One such type is Dyck words.

Dyck Words

A “Dyck word” is a string of balanced parentheses.

()()
(()())()
(())()

It’s (apparently) well-known that these strings are isomorphic to binary trees (although the imperative descriptions of algorithms which actually computed this isomorphism addled my brain), but what made them interesting for me was that they are a flat type, structured like a linked list, and as such should be reasonably straightforward to prove to be finite.

Our first task, then, is to write down a type for Dyck words. The following is a first possibility:

data Paren = LParen | RParen
type Dyck = [Paren]

But this type isn’t correct. It includes many values which don’t represent balanced parentheses, i.e. the expressions [LParen,RParen] :: Dyck are well-typed. To describe dyck words properly we’ll need to reach for the GADTs:

data DyckSuff (n :: Nat) :: Type where
  Done :: DyckSuff Z
  Open :: DyckSuff (S n) -> DyckSuff n
  Clos :: DyckSuff n     -> DyckSuff (S n)

type Dyck = DyckSuff Z

The first type here represents suffixes of Dyck words; a value of type DyckSuff n represents a string of parentheses which is balanced except for n extraneous closing parentheses. DyckSuff Z, then, has no extraneous closing parens, and as such is a proper Dyck word.

>>> Open $ Clos $ Open $ Clos $ Done :: Dyck
()()

>>> Clos $ Open $ Clos $ Done :: DyckSuff (S Z)
)()

>>> Open $ Open $ Clos $ Open $ Clos $ Clos $ Open $ Clos $ Done :: Dyck
(()())()

>>> Open $ Open $ Clos $ Clos $ Open $ Clos $ Done :: Dyck
(())()

The next task is to actually enumerate these words. Here’s an $O(n)$ algorithm which does just that:

enumDyck :: Int -> [Dyck]
enumDyck sz = go Zy sz Done []
  where
    go, zero, left, right :: Natty n -> Int -> DyckSuff n -> [Dyck] -> [Dyck]

    go n m k = zero n m k . left n m k . right n m k

    zero Zy 0 k = (k:)
    zero _  _ _ = id

    left (Sy n) m k = go n m (Open k)
    left Zy     _ _ = id

    right _ 0 _ = id
    right n m k = go (Sy n) (m-1) (Clos k)

>>> mapM_ print (enumDyck 3)
"()()()"
"(())()"
"()(())"
"(()())"
"((()))"

A variant of this function was what I needed in my thesis: I also needed to prove that it produced every possible value of the type Dyck, which was not too difficult.

The difficult part is still ahead, though: now we need to convert between this type and a binary tree.

Conversion

First, for the conversion algorithms we’ll actually need another GADT:

infixr 5 :-
data Stack (a :: Type) (n :: Nat) :: Type where
  Nil  :: Stack a Z
  (:-) :: a -> Stack a n -> Stack a (S n)

The familiar length-indexed vector will be extremely useful for the next few bits of code: it will act as a stack in our stack-based algorithms. Here’s one of those algorithms now:

dyckToTree :: Dyck -> Tree
dyckToTree dy = go dy (Leaf :- Nil)
  where
    go :: DyckSuff n -> Stack Tree (S n) -> Tree
    go (Open d) ts               = go d (Leaf :- ts)
    go (Clos d) (t1 :- t2 :- ts) = go d (t2 :*: t1 :- ts)
    go Done     (t  :- Nil)      = t

This might be familiar: it’s actually shift-reduce parsing dressed up with some types. The nice thing about it is that it’s completely total: all pattern-matches are accounted for here, and when written in Agda it’s clearly structurally terminating.

The function in the other direction is similarly simple:

treeToDyck :: Tree -> Dyck
treeToDyck t = go t Done
  where
    go :: Tree -> DyckSuff n -> DyckSuff n
    go Leaf        = id
    go (xs :*: ys) = go xs . Open . go ys . Clos

A Compiler

Much of this stuff has been on my mind recently because of this (2020) video on the computerphile channel, in which Graham Hutton goes through using QuickCheck to test an interesting compiler. The compiler itself is explored more in depth in Bahr and Hutton (2015), where the algorithms developed are really quite similar to those that we have here.

The advantage of the code above is that it’s all total: we will never pop items off the stack that aren’t there. This is a nice addition, and it’s surprisingly simple to add: let’s see if we can add it to the compiler presented in the paper.

The first thing we need to change is we need to add a payload to our tree type: the one above is just the shape of a binary tree, but the language presented in the paper contains values.

data Expr (a :: Type) where
  Val   :: a -> Expr a
  (:+:) :: Expr a -> Expr a -> Expr a

We’ll need to change the definition of Dyck similarly:

data Code (n :: Nat) (a :: Type) :: Type where
  HALT :: Code (S Z) a
  PUSH :: a -> Code (S n) a -> Code n a
  ADD  :: Code (S n) a -> Code (S (S n)) a

After making it so that these data structures can now store contents, there are two other changes worth pointing out:

• The names have been changed, to match those in the paper. It’s a little clearer now that the Dyck word is a bit like code for a simple stack machine.
• The numbering on Code has changed. Now, the HALT constructor has a parameter of 1 (well, S Z), where its corresponding constructor in Dyck (Done) had 0. Why is this? I am not entirely sure! To get this stuff to all work out nicely took a huge amount of trial and error, I would love to see a more principled reason why the numbering changed here.

With these definitions we can actually transcribe the exec and comp functions almost verbatim (from page 11 and 12 of 2015).

exec :: Code n Int -> Stack Int (n + m) -> Stack Int (S m)
exec HALT         st              = st
exec (PUSH v is)  st              = exec is (v :- st)
exec (ADD    is) (t1 :- t2 :- st) = exec is (t2 + t1 :- st)

comp :: Expr a -> Code Z a
comp e = comp' e HALT
  where
    comp' :: Expr a -> Code (S n) a -> Code n a
    comp' (Val     x) = PUSH x
    comp' (xs :+: ys) = comp' xs . comp' ys . ADD