Baby's First Type

If we're going to use van Laarhoven's trees, we have to create a type and define a meet operation on it that satisfies the semilattice properties, so it can take its rightful place in the Semilattice type class. Then we can use fromList to create a tree and search it.

(By the way, we've gone a long way with just a few built-in types. Initially we were using floating point, but that was in fact a bug because of the size of the integers we have to deal with. That leaves Int (fixed-size integer), Integer (arbitrary-precision integer aka "bignum"), and lists. Aside from Lisp, how many of the languages atop the TIOBE Index could one do this in purely with built-in types?)

So, how do you create a type in Haskell?

In our earlier code fragment, we supposed we could use a (Int, Integer)to hold a d-digit Y value and its ordinal position in the list of all d-digit Y values (though not in that order!). Perhaps we can use type, which just lets you refer to a type by a name of your choice. That name becomes a "synonym" for the type:

type TaggedY = (Int, Integer)

Will that do? Not for us, alas. The whole point is to create an instance of Semilattice, and the fine print of instance says "...furthermore, [the type you're declaring to be an instance of a type class] must not be a type synonym" (Haskell 2010 Report, 4.3.2, "Instance Declarations"). Shucks!

What now? We have two choices: newtype or data.

newtype is what the Haskell Report refers to as a data type "renaming". Unlike type, it gives an actual constructor that you have to use to create something of the new type, and I suspect that is what lets you say a newtype type is an instance of a type class.

newtype TaggedY = MakeTaggedY (Int, Integer)

We're really just renaming (Int, Integer), so we'll take this route. We don't need the added power of data.

The ordinal position is just along for the ride, so for us, meet is just the maximum of the Ys... but wait. If you take a look, van Laarhoven already has

newtype Max a = Max { getMax :: a } deriving (Show)
instance Ord a => Semilattice (Max a)
    where meet (Max a) (Max b) = Max (max a b)

newtype Ge a = Ge a deriving (Show)
instance Ord a => Satisfy (Ge a) (Max a)
    where satisfy (Ge q) = (>= q) . getMax

Huh? He's setting up the general case of what we want to do.

That first newtype doesn't look quite like ours:

• It's polymorphic, i.e. it works for any type in Ord. This is serious DRY, and Haskell makes it easy. van Laarhoven does here what should be my goal as I advance in Haskell.
• What's that {getMax :: a}? We've seen something like the insides before; it says getMax has type a, but what does it mean in this context? You have to read a bit further through the report than I expected; in this context it's a deconstructor. If x has type a, then Max x stuffs it into a Max a which it returns; if y has type Max a, getMax y gives you back the value of type a it contains.
• deriving (Show)? Show is the class of types that have printable representations. Remember that quick-and-dirty backwards, read . reverse . show? The ticket for a's admission to Show is a show function that you can hand an a and get back a String that represents it (and, if a is in Read, ideally read . show should be a (kind of slow) equivalent of id :: a, where id is the identity function, as you probably figured out). Haskell can figure out how to show some types for you, so you don't have to write your own, and this is one case.

The instance lines are what we're really after. The first one says "if a is in the type class Ord (the one for things that have a total order), then Max a is a Semilattice with max (working on the a values therein) as the meet function".

The second sets you up for queries looking for values greater than or equal to a given value. The Max a values in the tree get fed to getMax to extract the a, and then get compared with m.

But enough talk. I should have the courage of my convictions, and actually make the move to what should make the search through generated Ys take time on the O(logBase 2 (numNDigitYs n)) instead of O(numNDigitYs n). Back when the smoke clears.

Feeling like David Hilbert?

Remember the big brouhaha back in the late 19th and early 20th century about the foundations of mathematics? There was one radical school of thought, the Intuitionists, who insisted on constructive proofs. They held the method of "indirect proof", i.e. suppose what you wish to prove is false, and then show that a contradiction results, to be invalid, and disagreed with the "law of the excluded middle", that any given proposition is either true or false.

This was pretty controversial stuff, and David Hilbert wrote the following famous line about it: taking the principle of the excluded middle from the mathematician is the same as denying the telescope to the astronomer, or the boxer the use of his fists.

As a programmer brought up on imperative languages, I have to wonder whether those of similar training feel about pure applicative languages the way Hilbert did about intuitionism.

A must-read from Jeremy Bowers

"Why I'm Interested in Haskell"--not the writing of a fanboy, but of someone who I think fairly assessed the state of the programming language world at the time of writing (mid-April 2011).

A should-read from Evan Laforge

"What Haskell Doesn't Have"...thank goodness.

Advice on Haskell from that great philosopher, Bo Diddley

Don't let your instance write a check that your type can't cash.

Haskell's type checking is a glorious thing and saves one--let's be honest, saves me--from lots of stupid mistakes. However, there's one thing that it can't do for you. When you say a type is an instance of a type class, it has to trust that the operations you define on the type to gain admission to the type class behave the way they're supposed to. So be careful.

By golly...

...I think van Laarhoven's data structure is the way to go.

You'll recall that his "Search trees without sorting" stores the upper bounds of trees to speed up queries of the form "what's the first item with a key >= some value?" If we label the sorted d-digit Ys with their ordinal position:

addOrd :: [a] -> [(Int, a)]

addOrd xs = addOrd' xs 1
    where  addOrd' []   _ = []
           addOrd' x:xs n = (n, x) : addOrd' xs (n + 1)

(Is there a slick way to do that with a fold?  I'll have to look for one. UPDATE: no, but there's still a better way; see below)

then we put the results into one of van Laarhoven's semilattice trees and have

ysInRange m n d
   | nval == n = nPos - mPos + 1
   | otherwise = nPos - mPos
   where (mPos, mVal) = -- first entry w/ d-digit Y >= m
         (nPos, nVal) = -- first entry w/ d-digit Y >= n

OK, not quite, since n may be greater than the last d-digit Y, so we'd have to take care of that case--the easiest way would be to add a sentinel fake Y with value 10^d - 1 (i.e. tweak addOrd to add the sentinel value as a parameter that would be used in the base case).

ERRATUM: make that 10^d; we don't want anything to actually match the sentinel.

And BTW, on the way home from work it occurred to me that while Haskell's Maybe monad is the canonical way to indicate failure, you have to deal with the case in which you do find an actual Y that is greater than the upper bound, so I think a sentinel is better than Nothing.

UPDATE: Darn it, that sentinel will interfere with the lovely functioning of van Laarhoven's semilattice tree, because its maximum possible value swamps everything else. This needs a little more thought.

UPDATE: And indeed, the thing to do is create a data type that has a search tree and a separate sentinel. Then our search function returns

  case findFirst Ge n yTree of
    Just x  -> x
    Nothing -> sentinel

or something to that effect.

UPDATE: Not fold, zip:

addOrd = ([1..] `zip`)

zip takes two lists and gives you back a list of ordered pairs of corresponding elements, going on as long as there are corresponding elements. For example:

zip [1, 2, 4] [-1, 1, -1, 42] results in [(1, -1), (2, 1), (4, -1)]

The funky parentheses and `` are the way to construct a section. If you have a function f,

f :: a -> b -> c

then you can create two sections from f, one with an a value and one with a b value:

• (aValue `f`) bValue results in f aValue bValue.
• so does (`f` bValue) aValue.

 More concrete examples:

• (2 *) is a function that returns twice the value passed it
• (^ 2) is a function that returns the square of the value you pass it

The `` are just how you can use a normally prefix function as infix. We've been doing it all along with choose.

Great Leap Forward

So, we made the change. Note that we are in fact generating the lists for the combinations!

{-
  Here's our new approach: rather than generate the upper half separately and
  then have to reverse it a digit at a time, we generate both halves at the
  same time, when we have information we need to do so.

  So, out with backwards and {odd,even}DigitsPal, and twoTwos goes back to its
  simple self.
-}

twoTwos :: Int -> [Integer]

twoTwos n
    | n == 1    = []
    | even n    = [twoTwosNoOnes]
    | otherwise = [twoTwosNoOnes, twoTwosNoOnes + 10 ^ (n `div` 2)]               
    where twoTwosNoOnes = 2 * tenToThe (n - 1) + 2

oneTwos :: Int -> [Integer]

oneTwos n
    | n == 1     = [2]
    | even n     = []
    | otherwise  = [p + 2 * tenToThe halfN
                      | p <- justOnes n (min 1 (halfN - 1))]
    where halfN = n `div` 2

noTwos :: Int -> [Integer]

noTwos n
    | n == 1     = [1,3]
    | even n     = base
    | otherwise  = concat [[p, p + tenToThe halfN] | p <- base]
    where halfN = n `div` 2
          base = justOnes n (min 3 (halfN - 1)) 

choices :: Int -> Int-> [[Int]]

m `choices` n
    | n == 0           = [[]]
    | m == n           = [[m,m-1..1]]
    | otherwise        = [m : c | c <- (m - 1) `choices` (n - 1)]
                         ++ ((m - 1) `choices` n)

{-
  justOnes -- here's where the real action happens. Here we generate palindromes
  with a specified number of digits all of whose non-zero digits are 1s. 1 must
  thus be the most/least significant digit; in addition, between 0 and a
  specified number of 1s are scattered through the rest of each half. If the
  number of digits requested is odd, we will leave a 0 in the middle digit,
  which the caller may or may not replace with another digit.
-}

justOnes :: Int -> Int -> [Integer]

justOnes n o =
    let halfN = n `div` 2
        shell = tenToThe (n - 1) + 1
        innards = concat [(halfN - 1) `choices` k | k <- [0..o]]
        pair i = tenToThe i + tenToThe (n - (i + 1))
    in [shell + sum (map pair c) | c <- innards]

This gives a definite improvement. time output again varies a fair amount; the best we've seen so far is

real    0m0.509s
user    0m0.488s
sys     0m0.016s

and the worst is

real    0m0.581s
user    0m0.564s
sys     0m0.008s

I wish I knew the source of the variation; perhaps other processes running at the time? If we can believe the best result, we're down to 1.27 times the run time of the C++ solution we tested... and we haven't touched the real time sink. Check this out:

        total alloc = 345,155,656 bytes  (excludes profiling overheads)

COST CENTRE            MODULE  %time %alloc

numWithin              Main     44.7   11.3
solve                  Main     24.7   55.5
digitsIn.digitsIn'     Main      6.1    7.8
nDigitYs               Main      5.6    7.8
tenToThe               Main      5.3    0.0
justOnes               Main      3.0    6.4
justOnes.pair          Main      2.2    3.1
floorSqrt.floorSqrt'.y Main      2.1    1.9
noTwos                 Main      1.0    1.0
choices                Main      1.0    2.9

So we've cut down total allocation by over 200MB (not to be confused with the memory usage, which takes a final spike of about a quarter megabyte above the previous version). We could take a look at avoiding actually generating the list of chosen exponents for the half-Ys, but we have to face facts. numWithin is now eating nearly half the time, and if the above is to be believed, solve, of all things, is responsible for over half the allocation! (Maybe I don't understand just what the profiling output is trying to tell me.) We're within spitting distance of the C++ version; a move from O(n) to O(log n) for the range check ought to bring us ahead.

P.S. Note that I've deviated from the define-before-use ordering that I've stuck to in the past, and the results compiled without complaint. I infer from this that Haskell must be waiting until it's read the whole file to do at least some of the type checking, so I should be able to get away with putting definitions in top-down order.