Dynamic programming is a weirdly named way to speed up (complexity-wise) recursive algorithms. If you are writing a recursive algorithm that ends up calling itself with the same arguments over and over, dynamic programming can solve the problem by “memoizing” the results in a table. Basically, it just means build an array big enough to hold all the possible inputs, then fill in the array as you go. If you ever need the answer to a sub-problem you have already solved, just look it up in the array. The poster-boy for dynamic programming is the knapsack problem: imagine you are a thief who has broken into a jewellery store. You want to steal as much as you can to maximise the total value of the items you’ve stolen, but there is a maximum weight you can carry in your knapsack. Which items do you steal?
The problem with doing this in Haskell is that you can’t modify arrays (well, not easily anyway). Once you create an array, all of its values are fixed. This means you can’t just build an array full of dummy values and go about inserting the actual values later. But, as with many things in Haskell, laziness can save us.
Let’s start with a simple solution to the unbounded knapsack problem (you can steal as many copies of each item as you like), taken from the Wikipedia article:
-- knapsack [(v0, w0), ..., (vn, wn)] W knapsack items wmax = m wmax where m 0 = 0 m w = maximum $ 0:[vi + m (w - wi) | (vi, wi) <- items, wi <= w]
You call this by passing a list of (v, w) pairs (the value and weight of each item), and the maximum weight the bag can hold, W. You can use the example given in the image at the top of the Wikipedia page:
knapsack [(4, 12), (2, 2), (2, 1), (1, 1), (10, 4)] 15
Which gives the correct answer, 36. It will recursively compute the m function, which tells us the maximum value we can fit into a bag that can carry weight w. However, it is insanely slow because it recomputes sub-problems many times each. To be more formal, it has worst-case time of O(nW), because it is going to test n possibilities to a depth of W. Even this simple example with a maximum weight of 15 takes about four seconds to compute on my laptop. I just bought a new computer yesterday and I’m waiting for it to arrive, which should cut the time down a bit, but in the mean time, we should optimise the algorithm 🙂
Dynamic programming with lazy arrays
It seems like it will be impossible to use dynamic programming in Haskell, because arrays cannot be mutated. But if you think about it, dynamic programming doesn’t actually require that we mutate any values — it merely requires that we have some values that are uncomputed which later become computed. That’s exactly what laziness is all about.
So here’s the plan: we build an array whose indices and elements are the inputs and outputs of the function m. Because of laziness, the array won’t actually contain the outputs of the function — it will just contain thunks which are willing to compute the output if we ask them to. When we recurse, we won’t call a function, we’ll instead pull a value out of the array. If that particular array element has never been requested before, this will cause it to be computed (which might trigger more recursion). If that particular array element has already been found, this will simply return the existing known value. Here we go:
knapsack items wmax = arr ! wmax where arr = array (0, wmax) [(w, m w) | w <- [0..wmax]] m 0 = 0 m w = maximum $ 0:[vi + arr ! (w - wi) | (vi, wi) <- items, wi <= w]
We first build the array arr as an array of thunks. If any element of the array is requested for the first time, it will call the function m. Then, we try to pull out element wmax, and since it hasn’t been requested, this calls m wmax, which in turn, starts pulling more values out of arr. You can see that this will only ever call m a maximum of wmax times, and furthermore, it won’t necessarily call m for all possible values of w. It now has a worst-case time of O(nW), since each recursion makes n tests, and it will recurse a maximum number of W times. Indeed, it is now very fast in the 15 case. If we change the 15 to, say, 100000, it takes just a few seconds.
Note that all we had to change (besides adding the definition of arr) was replacing all calls to m with a lookup on arr. We should be able to generalise this solution.
A general dynamic programming solution
We would like a solution that allows people to write completely naïve code, such as our first attempt, and automatically turn it into a dynamic programming solution. This is possible in Mercury (see Tabled evaluation) as a special language feature which allows you to say “the runtime should memoize all past outputs of this function”. But we can’t do it without a special language feature (such as in Haskell), because the naïve version recurses on itself. We need to change the recursion slightly to allow our generic tabling code to insert itself in between the recursive calls. So let us go back to our basic solution, but modify m to accept a function argument called “self“, and call that when it recurses. (This is a bit like the Y combinator approach to recursion in Lambda calculus.)
m _ 0 = 0 m self w = maximum $ 0:[vi + self (w - wi) | (vi, wi) <- items, wi <= w]
Note that m calls self on the recursive call instead of m. Now it is possible for whatever calls m in the first place to insert some additional code in self. But for now, let’s just make it work as simply as possible (without any tabling). We write a function called table_null:
table_null :: (Ix a) => (a, a) -> ((a -> b) -> a -> b) -> (a -> b) table_null _ f = let g x = f g x in g
You can pass m to table_null and it will give you back a simple function that computes knapsack values. The first argument to table_null is supposed to represent the domain of the function, but we will ignore it for now. Note that if you pass m as the second argument to table_null, it returns our original m function — when m calls self, table_null‘s g will just call m again. So, this is equivalent to our original, naïve (very slow) knapsack solution:
knapsack items wmax = table_null (0, wmax) m wmax where m _ 0 = 0 m self w = maximum $ 0:[vi + self (w - wi) | (vi, wi) <- items, wi <= w]
So what is the point of this? Now we have abstracted away what it means to do recursion. When we use table_null, it just does plain ordinary recursion. But we could make any other table function we like, which can change the way recursion is handled. For example, we could build a table_trace that shows us just how inefficient this approach is by printing out the input to every recursive call:
table_trace :: (Ix a, Show a) => (a, a) -> ((a -> b) -> a -> b) -> (a -> b) table_trace _ f = let g x = trace (show x) (f g x) in g
If you have knapsack call table_trace, you will see it printing out the values, because table_trace‘s g prints out its input before calling the supplied function.
And now, a general solution to dynamic programming:
table_array :: (Ix a) => (a, a) -> ((a -> b) -> a -> b) -> (a -> b) table_array dom f = g where arr = array dom [(x, f g x) | x <- range dom] g x = arr ! x
If you have knapsack call table_array, it will suddenly become efficient. That’s because the g in table_array doesn’t always call f. First, table_array builds an array of thunks as we did above. The thunks inside the array call f the first time they are required. The g function merely picks values out of the array. So anything you use table_array on will automatically memoize the inputs to the function. The only caveat is that you have to know in advance the domain of inputs the function is likely to accept, and it has to be contiguous.