This Haskell cicada came about because sometimes, when I’m bored, I pick a random exercise from a random book and try to solve it. In this case the book was “Algorithm Design” by Jon Kleinberg and Éva Tardos.
"Algorithm Design takes a fresh approach to the algorithms course, introducing algorithmic ideas through the real-world problems that motivate them. In a clear, direct style, Jon Kleinberg and Eva Tardos teach students to analyze and define problems for themselves, and from this to recognize which ...
The exercise I randomly picked from the book is Exercise 8 on page 319:
The residents of the underground city of Zion defend themselves through a combination of kung fu, heavy artillery, and efficient algorithms. Recently they have become interested in automated methods that can help fend off attacks by swarms of robots.
Here is the exercise setup:
Robots attack at seconds $1, 2, \ldots n$. At second $i$ there are $x_i$ robots attacking. The people of Zion can deploy an electromagnetic pulse blaster which can destroy some number of robots as they arrive. This number is dependent on how many seconds have passed since it was last deployed (the blaster needs to recharge). We are given a charging function $f(t)$ which tells us how many robots can be destroyed given $t$ seconds of charging. So with the array of robot attacks $[x_1, x_2, \ldots, x_n]$, at second $k$ we can destroy $min(x_k, f(t))$ if $t$ seconds have passed since our last blast.
We are asked to choose the subset of seconds from $\{1, 2, \ldots n\}$ when we blast the electromagnetic pulse blaster such that we destroy the maximum number of robots.
The exercise is in the chapter on dynamic programming, so ahem: we will use dynamic programming.
There is a pretty obvious choice at second $n$: we either blast the EMP or we don’t. If we don’t, then the $x_n$ robots at second $n$ go through unharmed and don’t contribute to the total sum of destroyed robots. So if we don’t blast then we can defer to the smaller problem of $[x_1, x_2, \ldots, x_{n-1}]$ robots. If we do blast at second $n$, then we have $n$ choices to consider: one for each of the possible previous blasts (that is because the time of the previous blast influences how much the blaster can charge for this blast). There are $n$ possible previous blasts: at seconds $0, 1, 2, \ldots, n-1$ (second 0 means no previous blast).
Let $R(n)$ be the best solution if we have a blast at $n$. Then
$$ R(n) = max_{j = 0}^{n-1} (R(j) + min(x_n, f(n -j))) $$
If we populate the array $R$ then we can compute the best solution if we blast at $n$. $R(n)$ references all previous entries in the array, so we have to keep the whole array around.
We have to compare the blast choice with a solution where we don’t blast at $n$. Let’s denote $L(n)$ the overall solution. Then
$$ L(n) = max(R(n), L(n-1)) $$
And that’s it. We have to populate two arrays of size $n$ and the work at each step is linear, so the time needed is $O(n^2)$ and the space needed is $O(n)$.
Since this is a Haskell cicada, let’s implement it in Haskell. At this point it is a good idea to visit the Haskell Dynamic Programming page for a refresher.
The exercise doesn’t just ask us to sum up the maximum number of robots destroyed but also wants the seconds when we blasted. We need to keep track of indices in arrays, so a useful function will be maxInListWithIndex:
maxInListWithIndex :: [Int] -> (Int, Int)
maxInListWithIndex [] = (0, -1)
maxInListWithIndex [x] = (x, 0)
maxInListWithIndex (x:xs) = if x > y then (x, 0) else (y, (i+1))
where
(y, i) = maxInListWithIndex xs
The next function computes the array $R$. We are using Data.Vector and the generate vector creation from it:
blastsImpl :: [Int] -> (Int -> Int) -> (Vector (Int, [Int]))
blastsImpl xs f = r
where
r = generate ((length xs) + 1) (blasted r f xs)
This is very similar to the example in Haskell Dynamic Programming. The cool and surprising thing is that you can reference previous values of array $r$ in the function that you pass to generate. For us that function is blasted and here is its implementation:
blasted :: (Vector (Int, [Int])) -> (Int -> Int) -> [Int] -> Int -> (Int, [Int])
blasted _ _ [] _ = (0, [])
blasted _ _ _ 0 = (0, [])
blasted r f xs n = (mp, (n:iz))
where
x = xs !! (n - 1)
g :: (Vector (Int, [Int])) -> Int -> Int
g rr j = ((fst (rr!j))) + (min x (f (n - j)))
ps = [(g r j) | j <- [0..(n-1)]]
(mp, i) = maxInListWithIndex ps
iz = snd (r!i)
It is a direct translation of
$$ R(n) = max_{j = 0}^{n-1} (R(j) + min(x_n, f(n -j))) $$
into Haskell while also keeping track of the winning indices of the recursive blasts.
For the overall solution, we have to compare the blast vs non-blast totals and get:
overall :: (Vector (Int, [Int])) -> (Vector (Int, [Int])) -> Int -> (Int, [Int])
overall _ _ 0 = (0, [])
overall o r n = if b > c then (b, biz) else (c, oiz)
where
(b, biz) = r ! n
(c, oiz) = o ! (n - 1)
overallImpl :: [Int] -> (Int -> Int) -> (Int, [Int])
overallImpl xs f = o ! n
where
n = length xs
r = blastsImpl xs f
o = generate ((length xs) + 1) (overall o r)
The function overallImpl uses generate to compute the array $L$ by passing it the function overall which does the comparison of blast vs non-blast at specified second, referencing previous results.
You can grab the solution at the Haskell cicadas github repo.