Sententia cdsmithus

July 5, 2007

Find the Bug

Filed under: Haskell — cdsmith @ 2:47 pm

Here’s a pretty tricky bug I hunted down yesterday, with the help of oerjan on #haskell.  If you weren’t there to see, it’s a lot of fun.  Here’s the code:

foo n = product [24 - n .. 23] `div` product [1 .. n]
bar n = length (takeWhile ((< 1000000) . foo) [1 .. n])
baz n = if bar n == n then n else bar n

And here’s the result:

cdsmith@devtools:~/haskell/projects$ ghci Strange.hs
GHCi, version 6.7.20070703: http://www.haskell.org/ghc/  :? for help
[1 of 1] Compiling Main             ( Strange.hs, interpreted )
Ok, modules loaded: Main.
*Main> bar 10
9
*Main> baz 10
10

The question: Why does bar 10 evaluate to something different from baz 10?  By looking at the definition of baz, you might convince yourself that baz and bar do the same thing.  Apparently not.

Scroll down for the answer.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

The problem is numeric overflow.  The standard library length function has a result type of Int.  The intermediate result in foo is too large to fit in the data type Int.  The bar function has type forall a. Integral a => a -> Int.  When calling bar directly, this isn’t a problem because the argument defaults to type Integer and foo is evaluated entirely in arbitrary precision numbers.  However, when baz is evaluated, there is a comparison bar n == n, so the result of bar must have the same type as its argument.  The result is known to be Int because the return type of bar is not polymorphic.  So the argument n must also be of type Int.  Then foo's inferred type is Int -> Int, which causes the overflow.

All clear now?

If one wanted to take a lesson from this, it would be one of the following:

  • Be careful with types.
  • General purpose functions like length shouldn’t use dangerous types like Int.
  • Int should be made less dangerous by throwing an exception on overflow.

Let’s take a poll: which lesson do you draw from it?  Post in the comments.

June 25, 2007

Blogging in Haskell

Filed under: Haskell — cdsmith @ 8:12 pm

Yesterday, I finally decided to figure out syntax highlighting for Haskell in my blog.  I ended up finding two different ways to do it, so I’ll describe both of them, along with links to the (small bits of) code that I ended up writing in the process to help out.

Attempt 1: The Star-Light Approach

The first approach is to use Dean Edwards’ star-light JavaScript/CSS library.  This is a nice piece of client-side JavaScript and CSS that syntax highlights your source code for you, in a variety of languages.  All you have to do is enclose it in a tag like <pre class="javascript"> or <pre class="php">.  The CSS file (which you also need to reference from your page) takes care of invoking some JavaScript code that scours your page and does the syntax highlighting on the client.

Unfortunately, Dean Edward’s code doesn’t do Haskell by default.  So I wrote the code to do it, and you can download a modified version of star-light here.  It wasn’t terribly difficult, but neither is my code very elegant.  It mostly works, though.  Limitations that I know about are listed in the comments at the beginning of star-hs.js, and include poor handling of nested comments, “gaps” in string literals, and non-ASCII characters, and some overly excessive highlighting of Prelude symbols and the pseudo-keywords as, qualified, and hiding when they should be left alone.

Star-Light: Advantages

  1. Very little effort on your part.
  2. Support for several different languages.
  3. Client-side programming is all the rage!

Star-Light Disadvantages

  1. Uses non-standard CSS extensions, so limited portability.
  2. You have to trust my code. :)
  3. The parser is based on regular expressions, so it’s ultimately hopeless for Haskell.

In the end, though, neither of these is the reason that I didn’t end up using star-light for syntax highlighting.  I didn’t use it because the blogging host I’m using, WordPress, strips any “dangerous” stuff out of CSS and bans JavaScript from their blogs.  They claim it’s a security risk for them.  I’m a little worried, and you should be too, that a very popular blog platform believes that their security validly depends on client-side behavior.  I’d be even more worried if I thought they might be right.  Ah well; given the chance to choose again, I’d have avoided WordPress for a blog, but it’s a little late for that now.

Off to my second attempt…

Attempt 2: HsColour + Plugin

Having given up on client-side syntax highlighting, I turned to the obvious choice for server-side implementation: HsColour.  This nifty little project specializes in syntax highlighting Haskell, and can output HTML.  If only it were a little easier to use when writing a blog entry…

I write my blog entries using Windows Live Writer.  While I’m not at all addicted to it, it’s a nice alternative to the web-based editors in conventional blogs, and it doesn’t have the tendency to freeze up randomly that characterizes WordPress’s built-in editor.  (Did I mention that I’d never start a blog with WordPress again?  I thought so.)  It also has a plug-in interface, so I built a simple plug-in that asks for a block of code, runs it through HsColour, and sticks the result into a <pre> block in the blog.  Wanting to do the same for inline code fragments, I added a second plugin to ask for a line of code and put it into HTML <code> tags.

The plugin can be found here as a DLL file, which you can just drop into the plugins directory underneath Windows Live Writer’s installation directory to install it.  You’ll need HsColour installed and in your path.  If you want the source code, look here.  This plugin calls HsColour with the CSS option, so you’ll need to add a CSS stylesheet to define your syntax highlighting styles.  Alternatively, you could edit the plugin to use -html instead.  (WordPress, for example, charges extra for writing a CSS file even with their “security” limitations; did I mention I’d never start another blog with them?)  If you choose CSS, you’ll need styles for the selectors .keyword, .keyglyph, .layout, .comment, .conid, .varid, .conop, .varop, .str, .chr, and .num. The only confusing one is .layout, which I initially assumed had something to do with omitting squiggly braces. It turns out it’s just more reserved symbols and should probably be set to the same thing as .keyglyph.

HsColour + Plugin: Advantages

  1. Result works on all browsers with basic CSS support.
  2. The planet.haskell.org server has (boring) HsColour CSS entries, so it works there!
  3. HsColour is probably better at correctly parsing the language.

HsColour + Plugin: Disadvantages

  1. Only Haskell works; other techniques needed for other languages.
  2. My plugin only works on this one piece of blog writing software.

Attempt 2.5: Another Loose End

Another annoying thing about WordPress is that their blog software converts normal old everyday quote characters to “smart quotes”.  This is merely annoying for regular text, but it’s absolutely fatal for source code.  (Did I mention I’d never start another blog with WordPress?)  One more quick change produces a half-fix for this.  You should only need this if you are using WordPress as well.  The idea is to hide quotes from the WordPress code-mangler by writing unnecessary HTML entities.  Adding the appropriate HTML escape codes (which WordPress won’t let me write, but are an ampersand, followed by #, followed by either 39 or 34, followed by a semicolon) at line 79 of HsColour’s HTML.hs does the trick.

This sort of fixes the problem.  Unfortunately, Windows Live Writer helpfully notices that these entity tags are not needed, and converts them back to spaces for me every once in a while.  For the time being, my strategy for solving this is to be vigilant.  If you notice a problem with smart quotes in source code in my blog, just say so and I’ll try to fix it.

That’s all I’ve got.  Hope it was helpful.

June 24, 2007

Learning Haskell and Number Theory: GCD and Higher Order Functions

Filed under: Haskell, math — cdsmith @ 9:24 pm

This is the fourth installment of my series “Learning Haskell and Number Theory”, in which I work through Gareth and Mary Jones’ text book, Elementary Number Theory, and translate its ideas and algorithms into Haskell.  In the last installment, I wrote a function to compute the greatest common divisor of two numbers.  In this installment, we will play with higher order functions, and culminate with the extension of greatest common divisors to an arbitrary number of inputs.

The second of the two GCD implementations in the previous installment looked like this:

gcd2 0 0 = undefined
gcd2 a b = let next (x,y)    = (y, x `rem` y)
               vals          = iterate next (abs a, abs b)
               Just (ans, _) = find ((== 0) . snd) vals
           in  ans

The expression iterate next (abs a, abs b) is interesting, in that its first parameter (next) is itself a function. This general idea can be extended to all sorts of tasks, and is one of the cornerstones of functional programming.  It goes by a name: higher order functions.  (The general intuition here is that a first order function is one that operates on basic pieces of data, such as integers or strings; a second order function is one that operates on first order functions, and so on.  I don’t like this intuition because it ignores the fact that what’s chosen as “simple” data in any given language is arbitrary; but I don’t get to pick these names.)  A higher order function can take other function as parameters, and it can evaluate to a function as its result.  If it does both of these, it can be seen as a sort of transformer of functions.

I intend to write a higher-order function, and then talk about some of them that are provided in the Haskell Prelude.  Before I can get there, though, I need to clean up something that many of you pointed to in comments earlier, and that’s the meaning of the greatest common divisor of zero and zero.  The Jones text calls it undefined, which makes very literal sense.  Everything divides zero, and the set of integers has no greatest element.  Therefore, there can be no greatest divisor of zero.  However, as many of you pointed out, defining gcd(0, 0) = 0 has some nice benefits.

One of those benefits is that it allows me to continue where I want to go.  (I’m still following the book, but in fact the book is wrong here.  As an exercise, if you’ve got the book, see if you can find the error in the book’s answer to exercise 1.9.  There must be an error, but the statement for which it wants a proof is not true!)

Therefore, we’ve got:

gcd1 :: forall a. Integral a => a -> a -> a
gcd1 a b = gcd1' (abs a) (abs b)
    where gcd1' a 0 = a
          gcd1' a b = gcd1' b (a `rem` b)

Suppose we want the greatest common divisor of more than two numbers.  Following example 1.9 of the Jones text, we could define:

gcd' :: forall a. Integral a => a -> a -> a -> a
gcd' a b c = gcd (gcd a b) c

gcd'' :: forall a. Integral a => a -> a -> a -> a -> a
gcd'' a b c d = gcd (gcd (gcd a b) c) d

gcd''' :: forall a. Integral a => a -> a -> a -> a -> a -> a
gcd''' a b c d e = gcd (gcd (gcd (gcd a b) c) d) e

Oops… looks like we’re going on forever here. We could define gcd'''' or even worse.  Better is to exploit this fact:

*NumTheory> :qc \a -> gcd1 0 a == abs a
+++ OK, passed 100 tests.

Great!  So that gives us a starting point, 0, into which we can add numbers one by one until we have the greatest common divisor of them all.  We need a function that converts a list of integral values into a single integral value representing the GCD.

gcdmany :: forall a. Integral a => [a] -> a
gcdmany [] = 0
gcdmany xs = let x = init xs
                 y = last xs
             in  gcd1 (gcdmany x) y

The new functions init and last do what they say. The init function returns the initial part of a list — more precisely, all but the last element. The last function returns the last element of the list. Note that x is again a list!  There is, of course, one more thing that may be tricky about this code: it is recursive, like most useful things in Haskell.  Some careful thought, though, should convince you that this does exactly what the book suggests: from left to right, it applies the greatest common divisor algorithm to each pair of numbers until it gets to the end.

It’s quite inefficient, though.  At turns out that in haskell, it’s much faster to pick the first element from a list than the last one, because lists are singly linked in the forward dimension.  A more complicated, but faster, implementation is:

gcdmany :: forall a. Integral a => [a] -> a
gcdmany xs = gcdmany' 0 xs
    where gcdmany' p []     = p
          gcdmany' p (x:xs) = gcdmany' (gcd1 p x) xs

If you’ve got the hang of Haskell, you can check that this works.  The only thing I haven’t mentioned yet is the list constructor, :, which is basically LISP’s cons function if you know other functional languages.  If the code remains a mystery to you (which would be normal if you’re new to the language) you can study it, or just trust me.  The problem is that this is a little more complicated than we’d like to write again and again, and it turns out that you want to do exactly this lots and lots in Haskell code!  A higher order function would be handy, and the standard library provides one, called foldl, that does exactly what we’re looking for.

gcdmany :: forall a. Integral a => [a] -> a
gcdmany xs = foldl gcd1 0 xs

The last thing we can do, just be fancy, is called eta reduction.  This is the observation that if f(x) = g(x) for all x, then f = g.  This property is sometimes known as extensional equivalence, and it’s a property of pure functions, by definition. The xs at the end of the left and right-hand sides of the Haskell equation match up just like the xs in the math equation, and we can remove them.

gcdmany :: forall a. Integral a => [a] -> a
gcdmany = foldl gcd1 0

And that’s pretty simple.

June 22, 2007

Remaining Questions about HAppS

Filed under: Haskell — cdsmith @ 3:40 pm

This is a list of questions I need answered about HAppS’ MACID monad before I can decide to use it for a large-scale application.  It’s a long list, but I don’t mean to imply that the answers are bad; indeed, I’m hopeful they will turn out okay.  I’m making the list in the hopes someone can point me in the right directions, as just reading the source code is slow going so far!

  1. Is it feasible (even by extending the HAppS framework, if needed) to load balance a HAppS application across multiple servers?  I’ve got performance problems already with the old version of this application; I’m hopeful that many of them can be solved by dropping this Java object-relational mapper - which generates SQL that I never see - with Haskell code for which I can use smart data structures and algorithms and fix the problems.  Still, the hope is that more people use this in the future, and I’d hate to permanently give up the option of clustering or load balancing.  Although it wouldn’t be immediately needed, it remains worthwhile to have the choice.
  2. Is it feasible (even by extending the HAppS framework, if needed) to have HAppS store only part of its data at a time in memory?  While servers are getting more and more powerful over time, applications are using more and more data as well.  The data backup files for the application I’ve currently porting to Haskell are 9 GB in size, after compression.  I’m hoping that this increases over time.  While most of this data could potentially be moved into auxiliary files outside of the HAppS MACID monad, it’s still a little questionable whether I’d want to limit the data to the size of the server’s RAM.  The web site says that everything has to fit in RAM; but I don’t know if that’s because of something fundamental about the way the system works, or just because that’s what is implemented today.  The talk of moving to S3 and EC2 seems to indicate the latter (more promising) state of affairs.
  3. Is it feasible (even by extending the HAppS framework, if needed) to tie HAppS’ MACID transactions to transactions in external information systems via a two-phase commit protocol?  One of the key goals of the software I’m developing is integrating with other systems seamlessly, and that likely means allowing changes to HAppS-controlled data to participate in a distributed transaction.
  4. When I change the data structures for a HAppS application, it fails to load with the old state.  Is there an easy technique for migrating data when one needs to change the data tracked by the application?
  5. What is the best way to manage (i.e., interactively query and update) live data from a HAppS application outside of the pre-planned application functions?  For example, if I’m asked by my boss how many people with a first name of Tiffany have used the application after midnight since 17 days ago, I can currently throw together some SQL and find out.  Now, can I throw together some Haskell and find out?  How would I do that?  The only thing that comes to mind is something like lambdabot, which uses hs-plugins to compile and then dynamically load some code and run it.  I can do this, but it’s a little odd and it seems like there’s be an easy answer to this.

June 21, 2007

A Hack for Approaching Haskell Libraries

Filed under: Haskell — cdsmith @ 10:55 am

One of the most difficult problems I find in learning new libraries and modules is where to start.  If you start in the obvious places, there are inevitably zillions of abstractions that just don’t make sense yet.  As a quick kludge, I write the following Haskell code; it examines a set of modules in a directory, and sorts them by the number of dependencies on other modules.  The result is that if module A depends on module B, then B will be listed before A in the result (except for circular dependencies).  If you read the source code in the order given here, there’s more chance of understanding it all.

So here’s the code to build such a list:

import System
import System.Directory
import IO
import Data.List
import Data.Char
import Data.Maybe
import Data.Ord

{-
    Read a Haskell *.hs file and get a list of all the modules
    that it imports.
-}
findDeps base pkg = do
        let hi = base ++ (map dotToSlash pkg) ++ ".hs"
        ex <- doesFileExist hi
        if not ex then return [] else do
            src <- readFile hi
            let imps = filter isImport (lines src)
            return $ catMaybes $ map getTargetModule imps

    where dotToSlash '.' = '/'
          dotToSlash x   = x

          isImport (' ':t)  = isImport t
          isImport ('\\t':t) = isImport t
          isImport t        = "import" `isPrefixOf` t

          getTargetModule s = let pre = takeWhile (/= '(') s in
                              find (isUpper . head) (words pre)

{-
    Find the transitive, reflexive closure of the relation defined
    by the findDeps function.  This returns a list of ALL modules
    that this one uses or depends on, directly or indirectly.
-}
allDeps base mod = allDeps' [mod] [mod] where
    allDeps' (m:ms) full = do d <- findDeps base m
                              let d' = filter (not . flip elem full) d
                              t <- allDeps' (d' ++ ms) (d' ++ full)
                              return (m : t)
    allDeps' [] _        = return []

{-
    Usage: OrderByComplexity  

        = directory where source code is found.  This MUST
                end in '/'
     = file that lists the modules you're interested in,
                one per line.  This is often taken from a .cabal
-}
main = do [ base, pkgFile ] <- getArgs
          pkgStr            <- readFile pkgFile
          let pkgs           = lines pkgStr
          mods              <- mapM (allDeps base) pkgs
          let deps           = zip pkgs mods
          let deps'          = sortBy (comparing fst) deps
          mapM_ print $ map fst $ sortBy (comparing $ length . snd) deps'

June 18, 2007

Is TIOBE Fatally Flawed?

Filed under: Haskell, Uncategorized — cdsmith @ 7:38 pm

update: As Bogdy mentions in the comments, my reasoning here was based on false assumptions. It still seems clear that ranking APL above Haskell, along with other anomalies, disqualifies TIOBE for any serious purpose, at least past the top ten or so languages. My rankings should be ignored, though.

During a debate at work about using Haskell for a project, a coworker pointed out that Haskell is ranked #41 on the TIOBE.  On further investigation, things look really fishy.  Common interpretations of TIOBE include the amount of “community”, “buzz”, or “excitement” around a language.  By none of these standards can APL reasonably edge out Haskell.  I dug further.

Summary of findings: the TIOBE is severely broken.  It is falling victim to the fact that search engines grossly overestimate their number of results.  For example, if I search Google for “haskell programming”, as TIOBE does, the resulting page proudly estimates 44,500 results.  However, if I click through the results, I hit the end of the list after only 652.  Nice for marketing Google, perhaps, but it seems the estimate was rather poor.  Similar things happen with other languages.

TIOBE, despite using several search engines, seems to correlate well with Googles estimated (i.e., phony) number of results.  It correlates very badly with the actual number of results.  Here’s my corrected TIOBE list, built only from the top 50 languages in the original list.  In order to comply with Google’s terms of service, I painstakingly did this by hand; so I didn’t go any further.

There are some things that are initially surprising; but some thought indicates they may be reasonably expected.  Languages near the top tend to be those that are somewhat old (more time to write about them) or commonly used – past or present – in business and/or the academic world.  That’s because these languages have a reason to have a lot of web pages written about them.  One example: Prolog clearly isn’t a commonly used language nor one with a lot of community, but it’s taught by just about every computer science department in the world’s “programming languages” intro courses, because they feel better including something besides imperative and functional languages.  Hence, it’s been written about a lot.  One can see the effect of the “big community” effect though, if only in languages that appear above where you’d expect to see them.

I also split Lisp/Scheme into Lisp and Scheme separately, and dropped Natural because Googling for “natural programming” turned up more irrelevant results than relevant ones.

Without further delay, the “Chris” update to the TIOBE list.

  1. Fortran
  2. COBOL
  3. C
  4. Logo
  5. JavaScript
  6. MATLAB
  7. Prolog
  8. RPG
  9. ML
  10. Pascal
  11. Lingo
  12. Scheme
  13. LISP
  14. REXX
  15. C++
  16. Forth
  17. Smalltalk
  18. Icon
  19. SAS
  20. ABAP
  21. Tcl
  22. IDL
  23. FoxPro
  24. Haskell
  25. Bash
  26. Java
  27. CL
  28. APL
  29. ColdFusion
  30. Delphi
  31. Perl
  32. BASIC
  33. Objective C
  34. Erlang
  35. Lua
  36. Ada
  37. Awk
  38. ActionScript
  39. VBScript
  40. Ocaml
  41. D
  42. Dylan
  43. C#
  44. Python
  45. Ruby
  46. Transact-SQL
  47. PHP
  48. LabView
  49. S-lang
  50. PL/SQL

June 6, 2007

Learning Number Theory and Haskell: Greatest Common Divisor

Filed under: Haskell, math — cdsmith @ 3:12 pm

This is part three of my series, “Learning Number Theory and Haskell”, in which I work through Gareth and Mary Jones’ Elementary Number Theory and translate the ideas and concepts into Haskell. In the previous two parts, I used QuickCheck to verify some theorems and translate some ideas into Haskell code.

It is widely speculated that the first real algorithm worthy of the name was due to Euclid, when he described how to calculate the greatest common divisor of two terms. What better place to start, then in writing a first serious algorithm in Haskell? (Besides, the textbook does this next.) In the process, I’ll use “equational reasoning” (a popular word in functional programming) to improve the code. Then I’ll try something completely different to finish off

First, let’s write a type signature. Greatest common divisors can be calculated with integral values only, and maps two numbers to a result. The signature is:

gcd1 :: forall a. Integral a => a -> a -> a

For a first version of the algorithm, I plan to transcribe the proof from the text as code. There are a few trivial cases involving zero gotten near the top of page 6. We can write those separately using pattern matching.

gcd1 0 0 = undefined
gcd1 a 0 = abs a
gcd1 0 b = abs b

There are a few more special cases for non-zero arguments.

gcd1 a b | a <  0    = gcd1 (-a) b
         | b <  0    = gcd1 a (-b)
         | a == b    = a
         | a <  b    = gcd1 b a
                                 -- Phew!
         | otherwise = gcd1 b (a `rem` b)

Let’s test against the built-in gcd function from the prelude. The only subtlety here is that since gcd 0 0 is undefined, we should avoid calling it.

*NumTheory> :qc \a b -> (a /= 0 || b /= 0) ==> gcd1 a b == gcd a b
+++ OK, passed 100 tests.

We’ve got a working gcd1. However, it is terribly unsatisfying. By following the textbook to the letter, we ended up with a function defined with eight separate cases. Yuck!

The problem is that our number theory textbook handled a lot of stuff as special cases to get it out of the way easily for the proof. Many of those special cases aren’t all that special after all, and we can clean this up by finding equivalent cases and combining them. This is really easy to do and get right in Haskell, once you’ve got the hang of it, because Haskell functions can be read and substituted just like algebra equations. This concept is often called “equational reasoning”, and is cited as an advantage of pure functional programming. The ease of doing this is a direct consequence of not having side effects.

The last two cases cover when a < b, and a > b, respectively. However, the last case subsumes the next to last one, because if a < b (with a and b both positive), then a `rem` b is equivalent to a. That eliminates one case. Similarly, the third from last case is unnecessary. If a = b, then gcd1 a (a `rem` b) evaluates to gcd1 a 0, which in turn evaluates to abs a. Since we know a > 0 if the evaluation gets this far, this is the same as a. This case is also unnecessary. A similar argument proves that gcd1 0 b = b is unnecessary. That leaves:

gcd1 0 0 = undefined
gcd1 a 0 = abs a
gcd1 a b | a <  0    = gcd1 (-a) b
         | b <  0    = gcd1 a (-b)
         | otherwise = gcd1 b (a `rem` b)

This is better, but it’s even shorter to promote the check for negative values of a and b and for gcd1 0 0 to an outer function. The where keyword can be used to define an inner function that is available only within a particular top-level function. The result looks like this:

gcd1 0 0 = undefined
gcd1 a b = gcd1' (abs a) (abs b)
    where gcd1' a 0 = a
          gcd1' a b = gcd1' b (a `rem` b)

Eliminating all of these special cases could be interpreted as an indication that we didn’t need them to begin with.  As an exercise, try rewriting the proof in the text without assuming anything except that a and b are non-negative and are not both zero.

Rerunning the test suite to verify:

*NumTheory> :qc \a b -> (a /= 0 || b /= 0) ==> gcd1 a b == gcd a b
+++ OK, passed 100 tests.

Another implementation of the greatest common divisor is interesting (though, ultimately, it doesn’t perform as well.) One way to look at the algorithm is that we are inductively defining a list of number pairs which all have the same divisors. We start with (|a|, |b|) as the first element in the list, and then define a function to get from there to the next element. That looks like this:

gcd2 0 0 = undefined
gcd2 a b = let next (x,y)    = (y, x `rem` y)
               vals          = iterate next (abs a, abs b)
               Just (ans, _) = find ((== 0) . snd) vals
           in  ans

We first define a next function that takes us from one pair of numbers to the next. Then we use a Haskell standard library function, iterate, to build a list of values resulting from successive applications of next to the pair of numbers. Finally, we look for the first entry from that list that has 0 for its second number. The result is the first number of that pair. You can use QuickCheck to verify that it gives the same answer as the other two algorithms above.

This introduces an explicit data structure to represent the structure of the algorithm. Sometimes, doing this can make algorithms clearer. In this case, it’s just an interesting alternative. It involves some new concepts: an infinitely long list, a higher-order function, and another interesting use of pattern matching. If you don’t understand it yet, take some time to poke around the standard library documentation. Don’t despair, though; I’ll talk about more of these ideas in future installments.

June 3, 2007

Learning Number Theory and Haskell: More QuickCheck and Divisors

Filed under: Haskell, math — cdsmith @ 9:44 am

This is part two of my series Learning Number Theory and Haskell. Last time, we implemented a division algorithm according to an exact specification from Jones Elementary Number Theory. We used the QuickCheck testing library to make ourselves fairly confident that the code was actually correct (in other words, behaved according to certain expectations).

More QuickCheck

Since QuickCheck may be feeling pretty new right now, let’s play around with it for a bit. Exercises 1.1 and 1.2 of the Jones text give some interesting theorems about the remainders of perfect squares. We will verify these theorems with a hundred or so random test cases each. (Note that I’ve now stopped using the jonesQuot and jonesRem functions defined earlier. It was nice to be able to reproduce that behavior, but it’s nicer to become skilled at using the standard library functions instead of inventing our own replacements. We’ll just be careful about negative numbers from here on.)

Prelude> :qc \\n -> (n^2 `rem` 4) `elem` [0, 1]
+++ OK, passed 100 tests.
Prelude> :qc \\n -> (n^2 `rem` 3) `elem` [0, 1]
+++ OK, passed 100 tests.
Prelude> :qc \\n -> (n^2 `rem` 5) `elem` [0, 1, 4]
+++ OK, passed 100 tests.
Prelude> :qc \\n -> (n^2 `rem` 6) `elem` [0, 1, 3, 4]
+++ OK, passed 100 tests.

Yep, they all check out. Of course, you should still write proofs because, as you’re probably quite familiar by now, QuickCheck is suitable for contradicting things or providing evidence in favor, but not for proving things.

Generalizing the Proof Method

I was interested in this concept, so I extended it to arbitrary divisors with a simple Haskell program that generalizes the method of the proofs in the text. If you’re following along to learn Haskell, this makes use of two shorthands for creating lists. The innermost list [0 .. d - 1] is an arithmetic sequence, and hopefully it’s fairly obvious what that does. The outer list is built using a list comprehension. The vertical bar should still be read as “such that”, and in this context, <- is read as “element of”. So the outermost list is the list of all (n^2) `rem` d such that n is an element of [0 .. d - 1]. These techniques make a lot of mathematical code in Haskell very concise, and keep it closer to math syntax. Finally, nub just removes the duplicates from a list. You can :edit NumTheory.hs again, and add the following:

import Data.List (nub) -- Add this at the top of the file, by the other import

squareRems d = nub $ [(n^2) `rem` d | n <- [0 .. d - 1] ]

We can check this, but we’ll need to be careful about negative values of the divisor since I used rem instead of jonesRem.

Prelude> :l NumTheory
[1 of 1] Compiling NumTheory        ( NumTheory.hs, interpreted )
Ok, modules loaded: NumTheory.
*NumTheory> :qc \\n d -> (d > 0) ==> ((n^2) `rem` d) `elem` squareRems d
+++ OK, passed 100 tests.

This algorithm can be supported by a proof that is not given in the text. The theorem is that the remainder of a2 when divided by d is equal to the remainder of some k2, where 0 <= k < d. This is proven by writing a as qd + r by the division algorithm. The a2 = q2d2 + 2qdr + r2. In turn, r2 can be written as pd + s by the division algorithm. Then a2 = q2d2 + 2qdr + pd + s, which in turn is equal to d(q2d + 2qr + p) + s. Thus, a2 gives a remainder of s, but so does r2, and we know that 0 <= r < d. It should, therefore, be no surprise that the QuickCheck test succeeds.

Divisors

The next set of theorems from the text deal with divisors.  If we’d like to verify these theorems using QuickCheck, we need a way to decide if one number divides another one or not. Following Jones, we’ll say (b | a), read as “b divides a,” iff there is some integer q such that a = qb. Recall that for any pair of numbers a and b, with b not equal to 0, there is a unique pair of integers (q, r) such that a = qb + r and 0 <= r < b. In the special case that b is not zero, it’s clear that b | a iff r = 0 for the unique values of q and r given by the division algorithm.

b `divides` a = a `rem` b == 0

This leaves out the possibility that b = 0. The function above will fail if that occurs, because rem a b is undefined when b = 0. So we go back to the definition, which is that b divides a iff there is a q so that a = qb. If b = 0, that equation becomes a = 0b, which is only true if a = 0. Adding that as a special case, the function looks like this:

0 `divides` a  = a == 0
b `divides` a  = a `rem` b == 0

(This is a definition by pattern matching, which is one of the nicest syntax aspects of the Haskell language. The order is significant here, as the second version of the function will match all possible parameters that aren’t matched beforehand.)

I’ve convinced myself that the function is correct, but it’s time for some good old QuickCheck. The text gives a few corner cases explicitly at the top of page 4, so we can test them.

*NumTheory> :qc \\n -> n `divides` 0
+++ OK, passed 100 tests.
*NumTheory> :qc \\n -> 1 `divides` n
+++ OK, passed 100 tests.
*NumTheory> :qc \\n -> n `divides` n
+++ OK, passed 100 tests.

Excellent! Now we’d like to test the general case so we can be really sure that divides is correct… but how? All we have is that b divides a iff for some q, a = qb. We need to know something more about q than it merely being an integer to make this work. A convenient property is given in exercise 1.3 (d) from the Jones text. If a satisfactory q exists, it must divide a itself, so it is between -a and a. The case not covered by this reasoning is when a = 0, but then q = 0 will definitely do the trick! Let’s write a test.

*NumTheory> :qc \\a b -> (b `divides` a) == any (\\q -> a == q * b) [(- abs a) .. abs a]
+++ OK, passed 100 tests.

Here, I build a list of all possible values for q, and then use the built-in any function to see if any of them work. Notice that the first parameter to any is another function, which I’ve written using the anonymous lambda syntax. This is far too slow to be used for an implementation, but it works fine as a test.

Testing a Property of Divisors

Let us now make use of our new function. Exercise 1.4 of the text asks whether it is the case that if a | b and c | d, then (a + c) | (b + d). If you’re following along with me in the text, which I hope you are, you should convince yourself of the answer ahead of time by more mundane means. QuickCheck easily gives the answer.

(Here’s a comprehension quiz. The mere fact that QuickCheck can provide a definite answer should tell you what the answer is. Why?)

*NumTheory> :qc \\a b c d -> ((a `divides` b) && (c `divides` d)) ==> (a + c) `divides` (b + d)
*** Failed! Falsifiable (after 5 tests and 7 shrinks):
1
1
1
0

So not only do we get an answer (no, the conjecture is false!), but we get a counter-example. If a = 1, b = 1, c = 1, and d = 0, then it is true that a | b because b = (-2)a and c | d because d = (-1)c, but it is not true that (a + c) | (b + d).  a + c is 2, and b + d is 1, and clearly there is no integer q such that 1 = 2q.

June 2, 2007

Learning Number Theory and Haskell: The Division Algorithm

Filed under: Haskell, math — cdsmith @ 1:06 pm

This series of blog posts is a chronicle of my working my way through Gareth and Mary Jones’ Elementary Number Theory and translating the ideas into the Haskell programming language. This is the first installment of the series, which I’m calling “Learning Number Theory and Haskell.”

(As a quick side note, I know many readers already know Haskell. I hope this will get interesting for you soon, but my stated goal is to write for beginners… so you may want to skim the first few entries, or treat this as an opportunity to comment on whether I’ve chosen the best way to present ideas for new Haskellers.)

The Jones text begins with a description of the “division algorithm.” (If you’re following along in the text, which I hope you are, I’m combining theorem 1.1 with corollary 1.2 here.) This is seemingly an odd name for a theorem, but here it is. For all integers a and b, where b is not 0, there is a unique q and r such that a = qb + r and 0 <= r < |b|. The text calls these numbers the quotient and remainder. I call them an excuse to go hunting through the Haskell standard library a bit.

There are a few candidates in the standard library for the divison of numbers into a quotient and remainder. Haskell has quot and rem functions, which look promising. It also has div and mod functions, though. Let’s play around a bit and see which (if any) do what we want.

Since I said I’d start from the beginning, here’s the beginning. There are two commands to run GHC: ghc and ghci. We will start out in GHCi, which allows you to type lines of Haskell code and have them evaluated immediately. In some social circles, this is known as an REPL, which stands for “Read-Eval-Print Loop”. It makes poking around a lot easier, since you don’t need to save, recompile, or write a main function before seeing results. Later on, we can use the ghc command to compile programs or libraries that can be run on their own.

$ ghci -fglasgow-exts

(I added an option there that I haven’t mentioned yet. That option enables some language extensions that I don’t intend to use. I added it because it causes GHCi to use a slightly different form of type signatures that I think are actually easier to understand than the default.)

Finally, it’s time to ask GHCi for basic information about those functions. Specifically, we will ask for their types. One of the major keys to Haskell is that if you know the type of something, you are three fourths of the way to understanding it.

Prelude> :t quot
quot :: forall a. (Integral a) => a -> a -> a
Prelude> :t rem
rem :: forall a. (Integral a) => a -> a -> a
Prelude> :t div
div :: forall a. (Integral a) => a -> a -> a
Prelude> :t mod
mod :: forall a. (Integral a) => a -> a -> a

First of all, Prelude> is the basic GHCi prompt. Anything after that on a line is something I typed. We’re more interested in the responses, which are all about the same. Those are type signatures, and types are fundamental to Haskell. I think it might be helpful, though, to work our way up from simpler types toward the more complicated types above.

Prelude> :t not
not :: Bool -> Bool

The response to the :t command uses the :: syntax, which provides type annotations in Haskell. You should read :: as “has type”, so that statement is: not has type Bool -> Bool. That, of course, raises the question of what Bool -> Bool means as a type. First, Bool is the Boolean data type, with values of True or False. Second, -> describes a function mapping from one type to another. The not function, then, takes either True or False as an argument, and returns either True or False as a result. Simple enough.

Prelude> :t length
length :: forall a. [a] -> Int

There are two extra elements here. The first, which we won’t be using yet, is that [Bool] denotes a list whose elements are type Bool. The important thing is that in Haskell, you can’t just have a list. You must have a list of Bool or a list of Integer, or some other type. The length function, though, wants to operate on all possible types of lists. Therefore, it uses polymorphism. The keyword forall introduces a generic name for a type. It’s basically a logic quantifier, which in other contexts is often written as an upside-down A. The type here says that for all possible types a, the length function can calculate an Int from a list of a.

I like the -fglasgow-exts command line argument to ghci because it makes this explicit. Without that argument, the type above would have just been written as [a] -> Int. This is a valid abbreviation for the type given above; if some name in a type signature begins with a lower-case letter, it is assumed to be a type variable, even if there is no explicit forall clause. On the other hand, types beginning with upper-case letters are specific type names, like Int or Bool.

Prelude> :t negate
negate :: forall a. (Num a) => a -> a

The new concept introduced here is that of type classes. There are many different types that all represent numbers in Haskell. It would be painful if things like negate (or +, for that matter) could only work on one of those many types. (As a side note, I suppose it’s not universally agreed that this would be painful, because ML does it that way.) On the other hand, it makes no sense for negate to operate on all possible types, either! How do you keep a programmer from negating a banana warehouse? (Hmm, I feel like there should be a punchline here.) The solution is type classes. A type class, which starts with an upper-case letter much like a type, defines a group of many types (called instances) for which some set of operations is possible. You can read the => operator as an if/then statement, so the statement is: negate has a type so that for all types a, if a belongs to the type class Num, then negate can take an argument of type a and return a result of type a. Whew, that’s a mouthful!  It says what we want, though.  Now back to the beginning.

Prelude> :t quot
quot :: forall a. (Integral a) => a -> a -> a
Prelude> :t rem
rem :: forall a. (Integral a) => a -> a -> a
Prelude> :t div
div :: forall a. (Integral a) => a -> a -> a
Prelude> :t mod
mod :: forall a. (Integral a) => a -> a -> a

The only new thing here is that these functions take two parameters of type a, and return a result of type a. The Integral type class is exactly what you’d expect: it describes types whose values are mathematical integers. The two interesting types in this type class are Integer, which is an arbitrary precision type (in other words, it never overflows), and Int, which stores at least 30-bit signed values, depending on the implementation.

We now know the types, so let’s see how these functions behave. In other words, let’s play! (Any function of two parameters can be used as an infix operator by enclosing it in back-ticks, which are probably above the tilde on your keyboard. That’s what I’m doing below.)

Prelude> 6 `div` 4
1
Prelude> 6 `mod` 4
2
Prelude> 6 `quot` 4
1
Prelude> 6 `rem` 4
2

So far, so good. It looks like both pairs of functions behave like the Jones text’s quotient and remainder.  However, a single test case does not make a correct function! We could try a few more test cases.  Let’s try 100 of them instead! A nice systematic way to try lots of test cases is with a nifty library called QuickCheck. This library will verify stuff with bunches of random test cases. It tells us whether a property is false, or likely to be true. The wording there is very deliberate. It will definitely not tell us that the statement is true. That would require a general proof! QuickCheck will just try some random inputs to look for counter-examples.

I get tired of typing long commands, so I’ve added the following line to my GHCi configuration file, ~/.ghci:

:def qc \\c -> return $ ":m + Test.QuickCheck.Property \\nTest.QuickCheck.quickCheck (" ++ c ++ ") \\n:m - Test.QuickCheck.Property"

To use QuickCheck, we need functions that take parameters (values for a and b) and evaluate to either True (success) or False (failure). That’s convenient, because the Jones text gives us two of them. Specifically, the theorem asks for: a = qb + r, and 0 <= r < |b|. We can define these as anonymous functions using the lambda syntax in Haskell. Here’s a simple example of lambda syntax:

Prelude> (\\ x -> x * x) 4
16

The anonymous function begins with a backslash (\) and a list of parameters. The arrow (->) introduces the body of the function. This example defines an anonymous function to calculate the square of a number, and then applies it to 4 to get 16. We’ll do the same thing for our test cases. We are given that b is not zero, and QuickCheck defines a handy operator ==> to express that. Let’s try this with quot and rem.

Prelude> :qc \\a b -> (b /= 0) ==> (a `quot` b) * b + (a `rem` b) == a
+++ OK, passed 100 tests.
Prelude> :qc \\a b -> (b /= 0) ==> (0 <= a `rem` b) && (a `rem` b < abs b)
*** Failed! Falsifiable (after 9 tests and 5 shrinks):
-5
-3

Oops! The first test case succeeded, but the second test case failed. Specifically, it failed with values of a = -5 and b = -3. QuickCheck generates random test cases, so you will probably see a different test case fail. Hopefully your test did fail, though! Clearly, then, the quot and rem functions don’t satisfy the requirements of the division algorithm, because the result of rem is out of range for this input! Looking closer at the test case, we find that using quot and rem gives q and r of 1 and -2, while we wanted 2 and 1.

Let’s try div and mod.

Prelude> :qc \\a b -> (b /= 0) ==> (a `div` b) * b + (a `mod` b) == a
+++ OK, passed 100 tests.
Prelude> :qc \\a b -> (b /= 0) ==> (0 <= a `mod` b) && (a `mod` b < abs b)
*** Failed! Falsifiable (after 7 tests and 3 shrinks):
1
-2

So neither of the Haskell standard library functions agree with the spec that we got from the text. Both pairs of functions satisfy the first equality from the text, a = qb + r, all the time. If you’re wondering how they differ, the answer is that when the signs of a and b are different, rem gives the remainder the same sign as a, while mod gives the remainder the same sign as b. Thus, quot and rem fail our purposes precisely when a is negative and r is non-zero, while div and mod fail precisely when b is negative and r is non-zero.

We now know that we have two different incorrect functions. What about a correct one? We’ll have to write it; but fortunately it’s given in the Jones text. We’re in GHCi, and I’d prefer not to write a lot of code at a command prompt, so use GHCi’s :edit command.

Prelude> :edit NumTheory.hs

That should open up a text editor. (If you prefer a different editor, then set your EDITOR environment variable.) Now enter the following code, save the file, and exit the editor.

module NumTheory where

import Data.Ratio

jonesQuot, jonesRem :: Integral a => a -> a -> a

a `jonesQuot` b | b < 0 = ceiling (a % b)
                | b > 0 = floor   (a % b)

a `jonesRem` b = a - b * (a `jonesQuot` b)

Here’s what’s going on. The module statement defines this as a module with the name NumTheory. The import statement lets us use the Rational type, and the % operator to build it. I chose the Rational type to avoid the rounding error that would occur with floating point numbers. Next, there’s a type signature. Haskell could infer the type signature, but there are two advantages to specifying it ourselves. First, the signature inferred by Haskell is too general and is confusing. Second, the explicit type signature makes it abundantly clear how to use the function without having to ask GHCi for help. Next, jonesQuot has two cases corresponding to the main theorem and the corollary in the Jones text. (The vertical bar can be read as “such that” or “for” and is called a guard clause.) In one case, a floor is needed, while the other case requires a ceiling (this is in the book, if you’re reading it.) The % operator constructs a Rational from two Integral numbers.

Use the GHCi :l command to load the module, and you have precisely the functions we went searching for. We have QuickCheck to help us be more certain.

Prelude> :l NumTheory
[1 of 1] Compiling NumTheory        ( NumTheory.hs, interpreted )
Ok, modules loaded: NumTheory.
*NumTheory> :qc \\a b -> (b /= 0) ==> (a `jonesQuot` b) * b + (a `jonesRem` b) == a
+++ OK, passed 100 tests.
*NumTheory> :qc \\a b -> (b /= 0) ==> (0 <= a `jonesRem` b) && (a `jonesRem` b < abs b)
+++ OK, passed 100 tests.

Ah, the feeling of success! Of course, QuickCheck didn’t prove anything except that failure probably happens less than one out of a hundred inputs. Proofs of the theorem and simple algebra to establish that the algorithm is correct are in the Jones text, and you’re encouraged to read them. Just remember that there is no (okay, only a little) competition between testing and proof. Knuth’s famous quote comes to mind: “Beware of bugs in the above code; I have only proved it correct, not tried it.” Proofs don’t always catch boneheaded coding mistakes.

It’s worth a brief mention that, in this case, the interaction between testing and proof is particularly appealing. The text proves that there is one and only one function for which those two QuickCheck tests would always succeed. In other words, we actually have a proof that our test cases test everything of importance! That’s pretty cool; we have proven the completeness of the test suite. (The passing tests still don’t prove the correctness of the code, though. Why not?)

Learning Number Theory and Haskell: The Plan

Filed under: Haskell, math — cdsmith @ 9:42 am

Now that I’ve been added to planet.haskell.org, I feel like I should write a lot about Haskell. I also need to learn number theory over the summer. I thought, why not combine the two? Thus, I have hatched the idea of this “series” of blog posts, which I am calling Learning Number Theory and Haskell. My plan for this summer is to structure a set of writings on Haskell around my reading of the textbook Elementary Number Theory by Gareth and Mary Jones.

My writing will be aimed at beginners in both Haskell and number theory. My goal will be to:

  • Talk about all of the basic ideas of Haskell, using an exploratory approach.
  • Try for as much “text coverage” of the number theory book as I can.
  • Respond to your feedback on what you find interesting and worthwhile.

If you share similar interests, I’d encourage you to follow along. It may help if you actually purchase a copy of the textbook; if not, you will just miss out on much of the number theory part, as I don’t intend to rewrite the book in the course of my blog.

By way of acknowledgement, I partially got this idea from Rydeheard and Burstall’s book Computational Category Theory, which is available online here. I should point out that I haven’t actually read more than the first chapter of that book, and it is really not at all like what I’m going to do.  However, Rydeheard and Burstall’s work does raise the intriguing idea that it may be possible to build a software library that, in a sense, embodies the theory of a mathematical field.  I find that concept fascinating.

I may or may not create a library at all. The first few entries, which I’ve already sketched out, don’t involve writing anything permanent. I will be somewhat limited by the fact that I’m working at a beginner’s level in Haskell. The series will start with some hacking around in GHCi, and grow from there. Again, I hope to get some feedback as I go along on which directions to pursue, and which to abandon.

So, I hope you’re with me! I’ll be posting the first installment later today.

« Previous PageNext Page »

Blog at WordPress.com.