Contents

2008 Code Jam Beta

Triangle Trilemma

We compute the (signed) area of a triangle using the cross-product to. If zero, then the triangle is degenerate. Otherwise, we compare squares of distances to determine if the triangle is isosceles or scalene, then use the cosine rule to see if there is an obtuse (negative cosine) or right (zero cosine) angle.

We never need the distances between points; only their squares, and we only need the sign of the cosine. Thus we can restrict ourselves to integers since we only need addition, subtraction, and multiplication; we escape the vagaries of floating-point numbers.

import Jam
import Data.Bool
import Data.List

main = jam $ getints >>= \ns@[x1, y1, x2, y2, x3, y3] -> let
  a = (x2 - x1) * (y3 - y1) - (x3 - x1) * (y2 - y1)
  norm [a0, b0, a1, b1] = (a1 - a0)^2 + (b1 - b0)^2
  sqs = map norm [drop 2 ns, take 4 ns, take 2 ns ++ drop 4 ns]
  t = minimum $ (\[a, b, c] -> a + b - c) <$> permutations sqs
  in pure $ (++ " triangle") $ if a == 0 then "not a" else
    bool "isosceles" "scalene" (length (nub sqs) == 3) ++ " " ++
      case compare t 0 of LT -> "obtuse"
                          EQ -> "right"
                          GT -> "acute"

When finding the largest angle by looking for the smallest cosine, we’ve wastefully applied the rule to all permutations of the three points. This is terser than what I originally wrote:

t = minimum $ (\(a:b:c:_) -> a + b - c) . (`drop` cycle sqs) <$> [0..2]

and of course the extra computations are negligible.

The Price is Wrong

A subproblem is to find all longest increasing subsequences in the list, which we accomplish with a well-known dynamic programming method.

For each item, we record the longest subsequence up to that item, along with the next item required to start each such subsequence. This is the association list es in our code.

For the example given in the problem description, we wind up associate "google" with 3, the length of the longest subsequences that go up this item, along with "foo" and "bar", because one of the subsequences starts with "google" and "foo", and the other starts with "google" and "bar".

Brute force suffices afterwards: we enumerate every longest subsequence and remember the absent items written in alphabetical order (namely the items we can reprice to fix the sequence), then print the lexicographically smallest of these.

import Jam
import Data.List
import Data.Ord
import qualified Data.Map as M

main = jam $ do
  ws <- words <$> gets
  ns <- getints
  let
    pr = ((M.fromList $ zip ws ns) M.!)
    merge best new@(x, (m, ys)) = case find ((== x) . fst) best of
      Nothing -> new:best
      Just old@(_, (m', ys'))
        | m' > m    -> best
        | m' == m   -> (x, (m, ys ++ ys')) : delete old best
        | otherwise -> new : delete old best
    lcs best w = foldl' merge best $ (w, (1, [])):[(w, (m + 1, [x])) |
      (x, (m, ys)) <- best, pr x < pr w]
    es = foldl' lcs [] ws
    longest = maximum (fst . snd <$> es)
    starts = fst <$> filter ((longest ==) . fst . snd) es
    m = (M.fromList es M.!)
    f xs = concat [let ys = snd $ m x in case ys of
      [] -> [[x]]
      _  -> [x:y | y <- f ys] | x <- xs]

  pure $ unwords $ minimum $ map (sort . (ws \\)) $ f starts

Random Route

We use Dijkstra’s algorithm to find all shortest paths from the given start city, which naturally also finds all reachable cities. For each reachable city, the map best stores the shortest distance to the city, along with the list of the previous city in each such path. This is similar to the bread crumbs we left behind in the previous problem when computing longest increasing subsequences.

For the small input we may now use brute force: we enumerate each path to each city, and add the appropriate probability to each edge in the path, namely the reciprocal of the product of the number of reachable cities excluding the start city and the number of paths to this city.

To avoid any floating-point nastiness, we use arbitrary precision rationals and convert at the last minute. It’s probably fine to convert earlier, but I’d rather avoid thinking about it!

import Jam
import Data.List
import qualified Data.Map as M
import Data.Ratio
import Text.Printf

main = jam $ do
  [_n, start] <- words <$> gets
  let n = read _n
  es <- map ((\[v, w, d] -> (v, (w, read d))) . words) <$> getsn n
  let
    em = M.fromListWith (++) $ zipWith (\(v, _) i -> (v, [i])) es [0..]
    shorts (d1, is1) (d0, is0)
      | d0 < d1   = (d0, is0)
      | d0 > d1   = (d1, is1)
      | otherwise = (d0, is0 ++ is1)
    dijk acc done todo v
      | null vs = acc'
      | (v':todo') <- vs = dijk acc' (v:done) todo' v'
      where
        step m out = let (_, (w, dw)) = es!!out
          in M.insertWith shorts w (dv + dw, [out]) m
        acc' = foldl' step acc outs
        outs = M.findWithDefault [] v em
        (dv, _) = acc M.! v
        vs = (todo `union` map (fst . snd . (es!!)) outs) \\ done

    -- best :: Map city (best distance, [previous city])
    best = dijk (M.singleton start (0, [])) [] [] start

    allTo w
      | w == start = [[]]
      | otherwise  = concat [map (i:) $ allTo v |
        i <- snd $ best M.! w, let (v, _) = es!!i]

    count vs = zip (concat vs) $ repeat $
      1 % fromIntegral (length vs * (length (M.keys best) - 1))

    ds = M.fromListWith (+) $ concatMap (count . allTo) $ M.keys best

  pure $ unwords $ printf "%.7f" <$>
    [maybe 0 fromRational $ M.lookup i ds :: Double | i <- [0..n - 1]]

For the large input, we need a few more tricks. FIrstly, we use top-down dynamic programming to compute the number of shortest paths reaching each city, that is, we memoize the following recursion. If ways w is the number of ways to reach the city w by a shortest path, then:

ways w = sum [ways v | (v, w) <- inEdges w]

Next, we want to compute what we’ll call the weight of each city. Let dsts be the reachable cities. We define the weight of a city w to be:

weight w = sum [1 % ways v | v <- dsts, p <- bestPathsTo v, w `elem` p]

Since the number of shortest paths to a city might be exponential in the number of roads, computing the weight directly is infeasible. Instead, consider a reachable city w that is furthest away from the start city.

Trivially, w has a weight of 1, because a path goes through w if and only if it ends at w, so by definition we have a contribution of 1 from considering all the paths to w, and 0 from all paths to other cities. In fact, for any city v, the contribution from considering all paths to v is 1; the difficulty lies in computing the contributions from other paths.

With this in mind, we initialize the weight of each city to 1. (Actually, to 1 % 1 since we’re using arbitrary-precision rationals as much as possible.) We know the weight of w is correct, while we may need to increase the weights of other cities.

We can compute the contributions of the paths to w to the weights of adjacent cities:

M.fromListWith (+) [(v, weight w * (ways v % ways w)) | (v, w) <- inEdges w]

Now we remove w and recurse. In other words, in each step we take the city furthest away from the start, update the weight map with the above map, then remove the city. At each step, since the city we choose is furthest away, we know there are no paths left that can contribute to its weight, so the above multiplication by its weight uses the correct final value.

Part of this can be viewed as a topological sort on the subgraph consisting of reachable cities and shortest paths, except we’ve reversed all the directions of the roads because we’re starting from the furthest cities and working our way in.

In fact, in our code we break the above into a topological sort and an update phase. Our topological sort is inefficient because we stick to Haskell lists, and is also slightly confusing: again because of Haskell lists, we record the results backwards then reverse them all at the last minute.

At last we can figure out the probabilities for the edges. Let k be the number of reachable cities excluding the start city. For each edge from v to w that is part of a shortest path, the probability we travel on it is:

weight w * (ways v % ways w) / k

This has much in common with an expression in the weight update map, which suggests we could probably simplify somewhere. However, our solution works well enough.

I committed a series of blunders while writing this solution, most of which I found with the help of my brute force code.

  • During the topological sort, I considered all edges instead of only those that are part of shorest paths (the used list).

  • I skimmed too quickly when reading the guararantee about the input. I originally thought all input cities were reachable.

  • I also skimmed too quickly when reading about the output format. They want exactly 7 decimal places; it looks like over the years Code Jam has become more lenient in this regard, as after all, formatting issues like this are inconsequential details.

  • I considered unused roads when computing ways.

  • I made the opposite mistake and skipped the unused roads shen printing probabilities, instead of printing zero for them.

import Data.List
import qualified Data.Map as M
import Data.MemoTrie
import Data.Ratio
import Text.Printf

main = jam $ do
  [_n, start] <- words <$> gets
  let n = read _n
  es <- map ((\[v, w, d] -> (v, (w, read d))) . words) <$> getsn n
  let
    em = M.fromListWith (++) $ zipWith (\(v, _) i -> (v, [i])) es [0..]
    shorts (d1, is1) (d0, is0)
      | d0 < d1   = (d0, is0)
      | d0 > d1   = (d1, is1)
      | otherwise = (d0, is0 ++ is1)
    dijk acc done todo v
      | null vs = acc'
      | (v':todo') <- vs = dijk acc' (v:done) todo' v'
      where
        step m out = let (_, (w, dw)) = es!!out
          in M.insertWith shorts w (dv + dw, [out]) m
        acc' = foldl' step acc outs
        outs = M.findWithDefault [] v em
        (dv, _) = acc M.! v
        vs = (todo `union` map (fst . snd . (es!!)) outs) \\ done

    -- best :: Map city (best distance, [previous city])
    best = dijk (M.singleton start (0, [])) [] [] start
    waysRecurse v
      | v == start      = 1
      | M.member v best = sum $ map (ways . fst . (es!!)) $ snd $ best M.! v
      | otherwise       = 0
    ways = memo waysRecurse

    dsts = M.keys best
    used = concat $ snd <$> best

    topoSort acc = maybe (reverse acc) (topoSort . (:acc)) $
      find (\c -> all (\(v, (w, _)) -> v /= c || w `elem` acc) $
        map (es!!) used) (dsts \\ acc)
    ts = topoSort []

    update m w = foldl' (\m v -> M.insertWith (+) v
      (m M.! w * (ways v % ways w)) m) m $
        map (fst . (es!!)) (snd $ best M.! w)
    weights = foldl' update (M.fromList $ zip dsts $ repeat $ 1 % 1) ts

  pure $ unwords [printf "%.7f" $ if i `notElem` used then 0 :: Double else
    let (v, (w, _)) = es!!i in fromRational $
       weights M.! w * (ways v % ways w) / (fromIntegral (length dsts) - 1)
    | i <- [0..n - 1]]

In our analysis, for clarity, we have written expressions such as (v, w) for edges. In reality, we actually refer to an edge by its index in the es list. This is so we can print the probabilities in the same order as the roads given in the input.

Hexagon Game

Hexagonal boards can be treated as checkerboards where we can move orthogonally, and also along one diagonal but not the other. For example, for a board with S = 5, we can take the following subset of a checkerboard:

1  2  3
4  5  6  7
8  9 10 11 12
  13 14 15 16
     17 18 19

and stipulate that a piece can move orthogonally, or along the northwest-southeast axis.

To go from one square to another in the fewest steps, we move diagonally until the destination is in the northeast or southwest quadrants, at which point we must move orthogonally the rest of the way.

We multiply the number of moves required to reach a particular position by the value of each checker to compute the corresponding change to our score.

Because there are only three possible ending configurations, we may as well find the minimum solution for each of them, then print the smallest of the three. In any case, it’s unclear how computations for one configuration can help wtih computations for another.

For a given ending configuration, the obvious method is to compute the score for all permutations of checkers in the final row, but this is infeasible even for the small input, where S may be as high as 15:

3 * product [1..15] == 3923023104000

Instead we apply the Held-Karp dynamic programming trick that speeds up the Traveling Salesman Problem. Roughly speaking, this reduces the running time by replacing a factorial with a power of 2. Yet again, Data.MemoTrie comes to the rescue once we’ve jotted down the basic recursion.

import Jam
import Data.Bits
import Data.List
import Data.MemoTrie

main = jam $ do
  ps <- getints
  vs <- getints
  let
    s = length ps
    m = div (s + 1) 2

    offset r | r <= m    = 0
             | otherwise = r - m

    cols r   | r <= m    = m + r - 1
             | otherwise = s - (r - m)

    toRC r p | p <= cols r = (r, p + offset r)
             | otherwise   = toRC (r + 1) (p - cols r)

    dist (r0, c0) (r1, c1)
      | signum dr /= signum dc = ar + ac
      | otherwise = min ar ac + abs (ar - ac)
        where
          dr = r1 - r0
          dc = c1 - c0
          ar = abs dr
          ac = abs dc

    moves goal = map ((\x -> map (dist x) goal) . toRC 1) ps
    cost goal = zipWith (\v as -> map (v*) as) vs $ moves goal

    best goal = brute (2^s - 1 :: Int)
      where
        ts = transpose $ cost goal
        brute = memo bruteRaw
        bruteRaw bs
          | bs == 0   = 0
          | otherwise = minimum [cs!!b + brute (clearBit bs b) | b <- ones]
          where
            cs = ts!!(s - popCount bs)
            ones = filter (testBit bs) [0..s - 1]

  pure $ show $ minimum $ best <$>
    [ [(i, i) | i <- [1..s]]
    , [(m, i) | i <- [1..s]]
    , [(i, m) | i <- [1..s]]
    ]

For the large input, we notice that for a given ending configuration, we are trying to solve the assignment problem: we have an S x S matrix where the ith row and jth column is the cost of moving the ith piece to the jth position.

Thus we implement the Hungarian method:

import Jam
import Data.Array
import qualified Data.Bimap as B
import Data.List

main = jam $ do
  ps <- getints
  vs <- getints
  let
    s = length ps
    m = div (s + 1) 2

    offset r | r <= m    = 0
             | otherwise = r - m

    cols r   | r <= m    = m + r - 1
             | otherwise = s - (r - m)

    toRC r p | p <= cols r = (r, p + offset r)
             | otherwise   = toRC (r + 1) (p - cols r)

    dist (r0, c0) (r1, c1)
      | signum dr /= signum dc = ar + ac
      | otherwise = min ar ac + abs (ar - ac)
        where
          dr = r1 - r0
          dc = c1 - c0
          ar = abs dr
          ac = abs dc

    moves goal = map ((\x -> map (dist x) goal) . toRC 1) ps
    cost goal = zipWith (\v as -> map (v*) as) vs $ moves goal

    hungMax cs = f B.empty
      -- Initial potentials of nodes on either side.
      (listArray (1, s) $ map maximum cs)
      (listArray (1, s) $ repeat 0)
      -- Edge weights.
      (listArray ((1, 1), (s, s)) $ concat cs)

    -- Augment matching until perfect, adjusting potentials if needed.
    f ms lx ly w
      | B.size ms == s = sum $ map (w!) $ B.toList ms
      | otherwise      = g [] [root] [] lx ly
      where
        Just root = find (`B.notMember` ms) [1..s]
        -- Find a tight edge from ss to a node outside tt.
        g tr ss tt lx ly = case find (\(x, y) ->
            x `elem` ss && y `notElem` tt &&
            w!(x, y) == lx!x + ly!y) $ indices w of
          Nothing -> let
            a = foldl1' min [lx!x + ly!y - w!(x, y) |
              x <- ss, y <- [1..s] \\ tt]
            -- Adjust potentials to guarantee a tight edge exists.
            lx' = lx // [(v, lx!v - a) | v <- ss]
            ly' = ly // [(v, ly!v + a) | v <- tt]
            in g tr ss tt lx' ly'
          Just (x, y)
            -- If y is already matched, update ss and tt and keep searching.
            | B.memberR y ms -> g tr' (ms B.!> y:ss) (y:tt) lx ly
            -- Otherwise augment by path to y in the alternating tree.
            | otherwise      -> f (aug ms y tr') lx ly w
            where tr' = (y, x):tr
        aug ms y tr
          | x == root = ms'
          | otherwise = aug ms' (ms B.! x) tr
          where
            Just x = lookup y tr
            ms' = B.insert x y ms

    hungMin cs = -(hungMax $ map (0-) <$> cs)

  pure $ show $ minimum $ hungMin . cost <$>
    [ [(i, i) | i <- [1..s]]
    , [(m, i) | i <- [1..s]]
    , [(i, m) | i <- [1..s]]
    ]

The Hungarian algorithm can have cubic running time, but thankfully, the largest possible value for S is 75, so we can be sloppier. For example, we use a list to store the alternating tree.

We follow a description of the algorithm that seeks assignments of maximum weight. We want the mininum, so we simply negate all the weights first. An alternative would be flipping signs throughout the g function.

Our code handles matchings gracefully with Data.Bimap.

Contents
Ben Lynn blynn@cs.stanford.edu 💡