Tidying up a list comprehension in Haskell

Question

So I'm trying to generate a list of taxicab numbers in Haskell. Taxicab numbers are numbers that can be written as the sum of two distinct cubes in two different ways - the smallest is 1729 = 1^3 + 12^3 = 9^3 + 10^3.

For now, I'm just bothered about generating the four numbers that "make up" the taxicab number e.g. (1,12,9,10), and have been told to use a list comprehension (which I'm not too familiar with). This function will generate all 4-tuples where the largest number is at most n:

taxi n = [(a,b,c,d) | a <- [1..n], b <- [1..n], c <- [1..n], d <- [1..n], a^3 + b^3 == c^3 + d^3, a < b, a < c, c < d]

However, it is troublesome for a few reasons:

• The domain for a,b,c,d are all identical but I can't work out how to simplify this code so [1..n] is only written once.
• I want an infinite list, without an upper limit, so I can just leave the program running for as long as possible and terminate it when I please. Clearly if I just set a <- [1..] etc. then the program will never end up evaluating anything.
• The program is very slow: just taxi 50 takes 19 seconds.

Any speed optimisations would also be good, but if not the naive method I'm using is sufficient.

Answer 1

Your constraints imply a < c < d < b. So let b run outermost and let the others run in appropriate lower ranges:

taxi n = [ (a,b,c,d) | b <- [1..n],
                       d <- [1..b-1],
                       c <- [1..d-1],
                       a <- [1..c-1],
                       a^3 + b^3 == c^3 + d^3 ]

To go infinite, simply use b <- [1..].

---

A further big improvement is to compute one of the four variables from the other three:

taxi = [ (a,b,c,d) | b <- [1..],
                     c <- [1..b-1],
                     a <- [1..c-1],
                     let d3 = a^3 + b^3 - c^3,
                     let d = round(fromIntegral(d3)**(1/3)),
                     c < d,
                     d^3 == d3 ]

---

Benchmarking taxi 50 in GHCi with :set +s like you did:

Yours:          (16.49 secs, 17,672,510,704 bytes)
My first:        (0.65 secs,    658,537,184 bytes)
My second:       (0.09 secs,     66,229,376 bytes)  (modified to use b <- [1..n] again)
Daniel's first:  (1.94 secs,  2,016,810,312 bytes)
Daniel's second: (2.87 secs,  2,434,309,440 bytes)

Answer 2

The domain for a,b,c,d are all identical but I can't work out how to simplify this code so [1..n] is only written once.

Use [a,b,c,d] <- replicateM 4 [1..n].

The program is very slow: just taxi 50 takes 19 seconds.

One cheap improvement is to bake your a<b, a<c, and c<d conditions into the comprehension.

taxi n = [(a,b,c,d) | a <- [1..n], b <- [a+1..n], c <- [a+1..n], d <- [c+1..n], a^3 + b^3 == c^3 + d^3]

This makes things significantly faster on my machine.

Or, for better composability with the next (and previous) part of the answer, treat b, c, and d as offsets.

taxi n =
    [ (a,b,c,d)
    | a <- [1..n]
    , b_ <- [1..n], let b = a+b_, b<=n
    , c_ <- [1..n], let c = a+c_, c<=n
    , d_ <- [1..n], let d = c+d_, d<=n
    , a^3 + b^3 == c^3 + d^3
    ]

I want an infinite list, without an upper limit.

See my answer to Cartesian product of 2 lists in Haskell for a hint. tl;dr use choices.