# [Haskell-cafe] Code Review: Sudoku solver

Sat Apr 8 14:28:29 EDT 2006

```I have finished my cleanup of the dancing links based solver for Sudoku.

I don't have time to compare performance with the other programs that have been
posted recently, or to do more profiling of my code.

For those who will say "It is ugly, imperative, and ugly!" please remember this
is a conversion of Knuth's c-code, which depended on five non-trivial goto jumps
because he did not have tail recursion.  And the whole point of the algorithm
are the imperative unlink & relink routines acting on a sparse binary matrix.

This does not use any clever logic, but it does pick the tightest constraint at
every step.  This means if there is only one obvious possibility for a
position/row/column/block then it will immediately act on it.

The literate source file is attached.

[ My clever logic solver may eventually be cleaned up as well. ]

--
Chris

-------------- next part --------------
A Sukodku solver by Chris Kuklewicz (haskell (at) list (dot) mightyreason (dot) com)

I compile on a powerbook G4 (Mac OS X, ghc 6.4.2) using
ghc -optc-O3 -funbox-strict-fields -O2 --make -fglasgow-exts

This is a translation of Knuth's GDANCE from dance.w / dance.c

http://www-cs-faculty.stanford.edu/~uno/preprints.html
http://www-cs-faculty.stanford.edu/~uno/programs.html

I have an older verison that uses lazy ST to return the solutions on
demand, which was more useful when trying to generate new puzzles to
solve.

> module Main where

> import Data.Array.IArray
> import Data.STRef.Strict
> import Data.Char(intToDigit,digitToInt)
> import Data.List(unfoldr,intersperse)

> new = newSTRef
> {-# INLINE new #-}
> write = writeSTRef
> {-# INLINE write #-}
> modify = modifySTRef
> {-# INLINE modify #-}

Data types to prevent mixing different index and value types

> type A = Int
> newtype R = R A deriving (Show,Read,Eq,Ord,Ix,Enum)
> newtype C = C A deriving (Show,Read,Eq,Ord,Ix,Enum)
> newtype V = V A deriving (Show,Read,Eq,Ord,Ix,Enum)
> newtype B = B A deriving (Show,Read,Eq,Ord,Ix,Enum)

Sudoku also has block constraints, so we want to look up a block
index in an array:

> lookupBlock :: Array (R,C) B
> lookupBlock = listArray bb [ toBlock ij | ij <- range bb ]
>     where ra :: Array Int B
>           ra = listArray (0,pred (rangeSize b)) [B (fst b) .. B (snd b)]
>           toBlock (R i,C j) = ra ! ( (div (index b j) 3)+3*(div (index b i) 3) )

The values for an unknown location is 'u'.
The bound and range are given by b and rng.  And bb is a 2D bound.

> u = V 0  -- unknown value
> b :: (Int,Int)
> b = (1,9) -- min and max bounds
> rng = enumFromTo (fst b)  (snd b)  -- list from '1' to '9'
> bb = ((R (fst b),C (fst b)),(R (snd b),C (snd b)))

A Spec can be turned into a parsed array with ease:

> type Hint = ((R,C),V)
> newtype Spec = Spec [Hint] deriving (Eq,Show)

> type PA = Array (R,C) V

> parse :: Spec -> PA
> parse (Spec parsed) = let acc old new = new
>                       in accumArray acc u bb parsed

The dancing links algorithm depends on a sparse 2D node structure.
Each column represents a constraint.  Each row represents a Hint.
The number of possible hints is 9x9x9 = 271

> type (MutInt st)  = (STRef st) Int

The pointer types:

> type (NodePtr st) = (STRef st) (Node st)

The structures is a 2D grid of nodes, with Col's on the top of
columns and a sparse collection of nodes.  Note that topNode of Head
is not a strict field.  This is because the topNode needs to refer to

> type HeadName = (Int,Int,Int) -- see below for meaning

>                     ,topNode:: (Node st) -- header node for this column
>                     ,len:: !(MutInt st)  -- number of nodes below this head
>                     }

> data Node st = Node {getHint:: !Hint
>                     ,up,down,left,right :: !(NodePtr st)  -- two doubly-linked lists
>                     }

> instance Eq (Head st) where

> instance Eq (Node st) where
>     a == b = up a == up b

To initialize the structures is a bit tedious.  Knuth's code reads in
the problem description from a data file and builds the structure
based on that.  Rather than short strings, I will use HeadName as the
identifier.

The columns are (0,4,5) for nodes that put some value in Row 4 Col 5
(1,2,3) for nodes that put Val 3 in Row 2 and some column
(2,7,4) for nodes that put Val 4 in Col 7 and some row
(3,1,8) for nodes that put Val 8 in some (row,column) in Block 1

The first head is (0,0,0) which is the root.  The non-root head data
will be put in an array with the HeadName as an index.

> headNames = let names = [0,1,2,3]
>             in (0,0,0):[ (l,i,j) | l<-names,i<-rng,j<-rng]

A "row" of left-right linked nodes is a move.  It is defined by a

Initial hints are enforced by making them the only legal move for
that location.  Blank entries with value 'u = V 0' have a move for
all possible values [V 1..V 9].

> parseSpec :: Spec -> [Move]
> parseSpec spec =
>   let rowsFrom :: Hint -> [Move]
>       rowsFrom (rc@(R r,C c),mv@(V v')) =
>           if mv == u then [ rsyms v | v <- rng ]
>           else [ rsyms v' ]
>         where (B b) = lookupBlock ! rc
>               rsyms :: A -> Move
>               rsyms v = map ( (,) (rc,V v) ) [(0,r,c),(1,r,v),(2,c,v),(3,b,v)]
>   in concatMap rowsFrom (assocs (parse spec))

constructor and the recursive "mdo" notation as documented at
http://www.cse.ogi.edu/PacSoft/projects/rmb/

For more fun with this, see the wiki page at

> mkDList :: (MonadFix m) => (b -> a -> b -> m b) -> [a] -> m b
> mkDList _ [] = error "must have at least one element"
> mkDList mkNode xs = mdo (first,last) <- go last xs first
>                         return first
>   where go prev []     next = return (next,prev)
>         go prev (x:xs) next = mdo this <- mkNode prev x rest
>                                   (rest,last) <- go this xs next
>                                   return (this,last)

toSimple takes a function and a header node and iterates (read . function)
until the header is reached again, but does not return the header
itself.

>     where loop y = if header/=y then liftM (y:) (read (step y) >>= loop)
>                                 else return []
>

forEach is an optimization of (toSimple step header >>= mapM_ act)

>      where loop y = if header/=y then (act y >> (read (step y)) >>= loop)
>                                  else return ()

Now make the root node and all the head nodes. This also exploits mdo:

>     where makeHead before name after = mdo
>             ~newTopNode <- liftM4 (Node ((R 0,C 0),V 0) newHead) (new newTopNode) (new newTopNode)
>                                                                 (new newTopNode) (new newTopNode)
>                                    (new 0) (new after) (new before)

The Head nodes will be places in an array for easy lookup while building moves:

> hBounds = ((0,1,1),(3,9,9))
> type Root st =  (Head st,HArray st)

The addMove function creates the (four) nodes that represent a move and adds
them to the data structure.  The HArray in Root makes for a fast

> addMove :: forall st. (Root st) -> Move -> (ST st) (Node st)
>     where addNode :: (Node st) -> (Hint,HeadName) -> (Node st) -> (ST st) (Node st)
>           addNode before (hint,name) after = do
>             let head = ha ! name
>             let below = topNode head
>             above <- read (up below)
>             newNode <- liftM4 (Node hint head) (new above) (new below)
>                                                (new before) (new after)
>             write (down above) newNode
>             write (up below) newNode
>             seq l (return newNode)

Create the column headers, including the fast lookup array.  These
will be resused between puzzles.

> initHA :: (ST st) (Root st)
> initHA = do
>   heads <- toSingle next root
>   return (root,ha)

Take the Root from initHA and a puzzle Spec and fill in all the Nodes.

> initRoot :: (Root st) -> Spec -> (ST st) ()
> initRoot root spec = do
>   let moves = parseSpec spec

Return the column headers to their condition after initHA

> resetRoot :: (Root st) -> (ST st) ()
> resetRoot (root,ha) = do
>   let heads@(first:_) = elems ha
>         let node = topNode head
>         write (down node) node
>         write (up node) node
>       reset (last:[]) = do
>         write (prev root) last
>         write (next root) first
>         reset xs
>         reset xs

getBest iterates over the unmet constraints (i.e. the Head that are
reachable from root). It locates the one with the lowest number of
possible moves that will solve it, aborting early if it finds 0 or 1
moves.

> getBest root = do
>   first <- read (next root)
>   if first == root then return Nothing
>     else do
>       let findMin m best head | head == root = return (Just best)
>                               | otherwise = do
>             if l <= 1 then return (Just head)
>       findMin 10 first first

The unlink and relink operations are from where Knuth got the name
"dancing links".  So long as "a" does not change in between, the
undo the changes of cover and unconverOthers will undo coverOthers.

> unlink :: (a->STRef st a) -> (a->STRef st a) -> a -> (ST st) ()
> unlink prev next a = do
>   before <- read (prev a)
>   after <- read (next a)
>   write (next before) after
>   write (prev after) before

> relink :: (a->STRef st a) -> (a->STRef st a) -> a -> (ST st) ()
> relink prev next a = do
>   before <- read (prev a)
>   after <- read (next a)
>   write (next before) a
>   write (prev after) a

> cover :: (Head st) -> (ST st) ()
>   let eachDown rr = forEach right rr eachRight
>       eachRight nn = do
>         modify (len \$ getHead nn) pred
>   forEach down (topNode head) eachDown

> uncover :: (Head st) -> (ST st) ()
>   let eachUp rr = forEach left rr eachLeft
>       eachLeft nn = do
>         modify (len \$ getHead nn) succ
>   forEach up (topNode head) eachUp

> coverOthers :: (Node st) -> (ST st) ()
> coverOthers node = forEach right node (cover . getHead)

> uncoverOthers :: (Node st) -> (ST st) ()
> uncoverOthers node = forEach left node (uncover . getHead)

A helper function for gdance:

> choicesToSpec :: [(Node st)] -> Spec
> choicesToSpec = Spec . (map getHint)

This is the heart of the algorithm.  I have altered it to return only
the first solution, or produce an error if none is found.

Knuth used several goto links to do what is done below with tail
recursion.

> gdance :: (Head st) -> (ST st) Spec -- [Spec]
> gdance root =
>     let
>         forward choices = do
>                 Nothing -> if null choices
>                              then error "No choices in forward" -- return [] -- for [Spec]
>                              else do -- nextSols <- recover choices -- for [Spec]
>                                      return \$ (choicesToSpec choices) -- :nextSols -- for [Spec]
>
>             let endOfRows = topNode (getHead newRow)
>             if (newRow == endOfRows)
>               then do uncover (getHead newRow)
>                       if (null oldChoices)
>                         then error "No choices in advace" -- return [] -- for [Spec]
>                         else recover oldChoices
>               else do coverOthers newRow
>                       forward choices
>
>         recover (oldRow:oldChoices) = do
>             uncoverOthers oldRow
>             newRow <- readSTRef (down oldRow)
>
>     in forward []

Convert a text board into a Spec

> parseBoard :: String -> Spec
> parseBoard s = Spec (zip rcs vs'check)
>   where rcs :: [(R,C)]
>         rcs = [ (R r,C c) | r <- rng, c <- rng ]
>         isUnset c = (c=='.') || (c==' ') || (c=='0')
>         isHint c = ('1'<=c) && (c<='9')
>         cs = take 81 \$ filter (\c -> isUnset c || isHint c) s
>         vs :: [V]
>         vs = map (\c -> if isUnset c then u else (V \$ digitToInt c)) cs
>         vs'check = if 81==length vs then vs else error ("parse of board failed\n"++s)

This is quite useful as a utility function which partitions the list into groups of n elements.
Used by showSpec.

> groupTake :: Int->[a]->[[a]]
> groupTake n b = unfoldr foo b
>     where foo [] = Nothing
>           foo b = Just (splitAt n b)

Make a nice 2D ascii board from the Spec (not used at the moment)

> showSpec :: Spec -> String
> showSpec spec = let pa = parse spec
>                     g = groupTake 9 (map (\(V v) -> if v == 0 then '.' else intToDigit v) \$ elems pa)
>                     addV line = concat \$ intersperse "|" (groupTake 3 line)
>                     addH list = concat \$ intersperse ["---+---+---"] (groupTake 3 list)

One line display

> showCompact spec = map (\(V v) -> intToDigit v) (elems (parse spec))

The main routine is designed to handle the input from http://www.csse.uwa.edu.au/~gordon/sudoku17

> main = do
>   all <- getContents
>   let puzzles = zip [1..] (map parseBoard (lines all))
>   root <- stToIO initHA
>   let act :: (Int,Spec) -> IO ()
>       act (i,spec) = do
>         answer <- stToIO (do initRoot root spec
>                              answer <- gdance (fst root)
>                              resetRoot root