### Solving XKCD’s Nerd Sniping problem

A while ago, during a period of free time, I implemented a solution in Haskell for the XKCD’s raptor problem. Now, I’ll try to solve another problem presented there, the one found in the Nerd Sniping comic. Of course, the implementation will still be in Haskell.

First, we have to define a data structure for keeping the network of resistors. I’ll use a structure used a while ago in a C algorithm for solving the network analysis problem but more simplified than what I used at the time. We only need to be interested in networks formed with resistors only, thus we can use only node elimination to solve this problem. Because of this, the data structures used are very simple:

2
3 type Resistance = Double
4 type Conductance = Double
5
6 type Network a = [(a, Node a)]
7 type Node a = [WireTo a]
8 type WireTo a = (a, Resistance)
9

As you see, we’ll use two association lists: one to keep the mapping between node ids (those can be almost anything) and one to keep the mapping between the node id and the resistance of the wire connected between that node and the actual node. While this structure prevents the know tying problem (homework: try solving this problem using techniques presented in the linked article), it adds duplicated information: for each wire between nodes a and b we store it in two places: once in node a and once in node b.

The previous paragraph hinted that there are two problems which need to be solved. Firstly, when comparing two wires and two nodes only the ids need to be taken care of, not everything in the pair. Thus, we need to implement our own equality functions.

10 (===) :: (Eq a) => (a, b) -> (a, b) -> Bool
11 (a, _) === (b, _) = a == b
12
13 (=/=) :: (Eq a) => (a, b) -> (a, b) -> Bool
14 (=/=) a = not . (a ===)

To solve the second problem, I wrote two functions to the user and will demand (expect) that the user will use them to construct any network. Thus, the first function will construct a part of the network: a simple wire (similar to the return function from monads).

66 -- builds a simple network: a simple wire
67 buildPart :: (Ord a) => a -> a -> Resistance -> Network a
68 buildPart a b r
69   | r < 0 = error "Negative resistance is not allowed"
70   | a == b = error "No wire can have the same ends"
71   | a > b = buildPart b a r
72   | otherwise = [(a, [(b, r)]), (b, [(a, r)])]

Using this function we can already construct the first test network, as seen in the following picture:

117 {-
118 Simple test: a network with only a wire.
119 -}
120 testSimple :: Network Int
121 testSimple = buildPart 0 1 5

Simple network with only a resistor

The second function will take two networks and construct a new one by joining them (similar to the mplus function from MonadPlus or mappend from Monoid).

74 -- joins two networks
75 joinParts :: (Eq a) => Network a -> Network a -> Network a
76 joinParts [] ns = ns
77 joinParts (n:ns) nodes
78   | other == Nothing = n : joinParts ns nodes
79   | otherwise = (fst n, snd n `combineNodes` m) : joinParts ns nodes'
80   where
81     other = lookup (fst n) nodes
82     m = fromJust other
83     nodes' = filter (=/= n) nodes

Care must be taken when joining two networks containing the same nodes. To combine them correctly, the following two functions are used:

85 {-
86 Combines two nodes (wires leaving the same node, declared in two parts of
87 network.
88 -}
89 combineNodes :: (Eq a) => Node a -> Node a -> Node a
90 combineNodes = foldr (\(i,r) n -> addWireTo n i r)
91
92 {-
93 Adds a wire to a node, doing the right thing if there is a wire there already.
94 -}
95 addWireTo :: (Eq a) => Node a -> a -> Resistance -> Node a
96 addWireTo w a r = parallel \$ (a,r) : w

When we encounter two parallel resistances, we will transform them immediately.

98 {-
99 Reduce a list of wires by reducing parallel resistors to a single one.
100 -}
101 parallel :: (Eq a) => [WireTo a] -> [WireTo a]
102 parallel [] = []
103 parallel (w:ws) = w' : parallel ws'
104   where
105     ws' = filter (=/= w) ws
106     ws'' = w : filter (=== w) ws
107     w' = if ws'' == [] then w else foldl1 parallel' ws''
108
109 {-
110 Reduce two wires to a single one, if they form parallel resistors.
111 -}
112 parallel' :: (Eq a) => WireTo a -> WireTo a -> WireTo a
113 parallel' (a, r1) (b, r2)
114   | a == b = (a, r1 * r2 / (r1 + r2))
115   | otherwise = error "Cannot reduce: not parallel resistors"

Right now, we can construct several more tests, presented in the following pictures

123 {-
124 Second test: a network with two parallel wires.
125 -}
126 testSimple' :: Network Int
127 testSimple' = buildPart 0 1 10 `joinParts` buildPart 0 1 6
128
129 {-
130 Third test: a network with two series resistors.
131 -}
132 testSimple'' = buildPart 0 1 40 `joinParts` buildPart 1 2 2
133
134 {-
135 Fourth test: tetrahedron
136 -}
137 testTetra :: Network Int
138 testTetra
139   = buildPart 0 1 1 `joinParts`
140     buildPart 0 2 1 `joinParts`
141     buildPart 1 2 1 `joinParts`
142     buildPart 0 3 1 `joinParts`
143     buildPart 1 3 1 `joinParts`
144     buildPart 2 3 1

A collection of networks

Now, we can even build 2D networks by giving Cartesian id’s to nodes. For example, the following is the simplest instance of the XKCD problem (reduced to a minimum configuration).

146 {-
147 Fifth test: 2 squares
148 -}
149 testSquares :: Network (Int, Int)
150 testSquares = foldl joinParts [] . map (\(x, y) -> buildPart x y 1) \$ list
151   where
152     list = [((0, 0), (0, 1)), ((0, 1), (0, 2)), ((0, 2), (1, 2)),
153       ((1, 2), (1, 1)), ((1, 1), (0, 1)), ((1, 1), (1, 0)), ((0, 0), (1, 0))]

Small 2D network

Before going into defining the XCKD problem and solving it, we need to implement the algorithm for solving any network of resistors. The following code is not optimized and I am really sure that it can be tweaked a little to allow for infinite structures due to the laziness of Haskell (to solve this is left as a homework). First, the code tests whether we compute a valid answer or not.

16 {-
17 Starts the solving phase testing if each node is defined.
18 -}
19 solve :: (Ord a) => Network a -> a -> a -> Resistance
20 solve n st en
21   | stn == Nothing = error "Wrong start node"
22   | enn == Nothing = error "Wrong end node"
23   | otherwise = solve' n st en
24   where
25     stn = lookup st n
26     stnd = fromJust stn
27     enn = lookup en n

After we are sure that the nodes in question exist in the network, we start the solving process

29 solve' :: (Ord a) => Network a -> a -> a -> Resistance
30 solve' n st en
31   | null candidates = getSolution n st en
32   | otherwise = solve' (removeNode n (head candidates)) st en
33   where
34     candidates = filter (\(x,_) -> x /= st && x /= en) n

When we have only the start and end nodes we return the solution

63 getSolution :: (Eq a) => Network a -> a -> a -> Resistance
64 getSolution n st en = fromJust . lookup en . fromJust . lookup st \$ n

Otherwise, when we have a candidate node to be removed, we remove it:

36 removeNode :: (Ord a) => Network a -> (a, Node a) -> Network a
37 removeNode net w@(tag, nod)
38   | length nod == 1 = filter (=/= w) net
39   | otherwise = filtered `joinParts` keep `joinParts` new
40   where
41     affectedTags = map fst nod
42     -- construct the unaffected nodes list
43     unaffectedNodes = filter (\(a,_) -> a `notElem` affectedTags) net
44     keep = filter (=/= w) unaffectedNodes
45     -- and the affected ones
46     change = filter (\(a,_)-> a `elem` affectedTags) net
47     -- remove the node from all the maps
48     filtered = map (purify tag) change
49     -- compute the sum of inverses
50     sumR = sum . map ((1/) . snd) \$ nod
51     pairs = [(x, y) | x <- affectedTags, y <- affectedTags, x < y]
52     new = foldl joinParts [] . map (buildFromTags nod sumR) \$ pairs
53
54 purify :: (Eq a) => a -> (a, Node a) -> (a, Node a)
55 purify t (tag, wires) = (tag, filter (\(a,_)->a /= t) wires)
56
57 buildFromTags :: (Ord a) => Node a -> Conductance -> (a, a) -> Network a
58 buildFromTags n s (x, y) = buildPart x y (s * xx * yy)
59   where
60     xx = fromJust . lookup x \$ n
61     yy = fromJust . lookup y \$ n

We can test each of the previously defined networks to see if the algorithm works.

*Main> solve testSimple 0 1
5.0
*Main> solve testSimple' 0 1
3.75
*Main> solve testSimple'' 0 1
40.0
*Main> solve testSimple'' 0 2
42.0
*Main> solve testTetra 0 3
0.5
*Main> solve testSquares (0, 0) (1, 2)
1.4000000000000001

Now, the solution to the XKCD’s problem. First, we define the network: we will take a finite case, between $(-n,-n)$ and $(n,n)$. By increasing n we will get closer and closer to the actual result (as you will see later, the convergence is pretty good).

155 build :: Int -> Network (Int, Int)
156 build n = foldl joinParts [] . map (\(x, y) -> buildPart x y 1) \$ list
157   where
158     list = [(x, y) | x <- points, y <- points, x `neigh` y, x < y]
159     points = fillBelow n
160
161 fillBelow :: Int -> [(Int, Int)]
162 fillBelow n = [(x, y) | x <- [-n .. n], y <- [-n ..n]]
163
164 neigh :: (Int, Int) -> (Int, Int) -> Bool
165 neigh (a, b) (c, d) = abs (a - c) + abs (b - d) == 1

An instance of the XKCD problem

Thus, to solve the problem for one iteration we have

167 solveXKCD :: Int -> Resistance
168 solveXKCD n = solve (build n) (0, 0) (1, 2)

And the main, used for printing the results

170 main = mapM (print . \x -> (x, solveXKCD x)) [3..]

And now the results. I let the program run until it reached a network of size 25 then I stopped it and plotted the results.

Resistance vs grid size, see the convergence speed

From the plot, it is easy to see that the valid result lies somewhere between $.774$ and $.772$ range, just like the valid answer is.

In the end, I’d like to relate to another post from that topic: please add more puzzles like this so that I can do something when I have free time instead of slacking off.