Tetrahedral Numbers

September 13, 2011

Before I started writing my own programming exercises, I enjoyed solving other programming exercises available on the internet, many of them mathematical in nature. Today’s exercise comes from http://open-cs.net/problems.php?id=18 and concerns triangular and tetrahedral numbers.

A triangular number tells the number of ways that balls can be stacked in a triangle. The first triangular number is 1, the second is 3 (a row of 1 plus a row of 2), the third triangular number is 6 (the first two rows plus a row of 3), the fourth triangular number is 10 (adding a row of 4), the fifth triangular number is 15 (think of the 15 balls on a pool table), and so on.

A tetrahedral number is the three-dimensional equivalent of a triangular number; think of cannonballs stacked in a three-sided pyramid. The top layer has one ball, the second layer has 3 balls (the second triangular number) so the second tetrahedral number is 1 + 3 = 4, the third layer has 6 balls giving a total of 10 balls in the tetrahedron, and so on.

Your task is to find the base of the tetrahedron that contains 169179692512835000 balls. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.

About these ads

Pages: 1 2

15 Responses to “Tetrahedral Numbers”

  1. Tante Hedwig said

    import Data.List

    tria = tail $ scanl (+) 0 [1..]
    tetra = 1:zipWith (+) tetra (tail tria)

    main = do
    let (Just n) = findIndex (==169179692512835000) tetra
    putStrLn $ show n

  2. My Haskell solution (see http://bonsaicode.wordpress.com/2011/09/13/programming-praxis-tetrahedral-numbers/ for a version with comments):

    import Data.List
    
    triangular :: [Integer]
    triangular = scanl1 (+) [1..]
    
    tetrahedral :: [Integer]
    tetrahedral = scanl1 (+) triangular
    
    main :: IO ()
    main = print . maybe 0 succ $
           findIndex (== 169179692512835000) tetrahedral
    
  3. [...] today’s Programming Praxis, our goal is to find the base of the three-sided pyramid that has [...]

  4. Graham said

    For my Python solution, I wrote a linear search as well as a binary one:

    def tetra(n):
    return n * (n + 1) * (n + 2) / 6

    def prob18_linear(target):
    n = 1
    while tetra(n) < target:
    n += 1
    return n

    def bin_search(target, func):
    lo, hi = 1, 10
    while func(hi) < t:
    hi *= 10
    mid = (lo + hi) / 2
    fmid = func(mid)
    while fmid != target:
    if fmid < target:
    lo, mid = mid, (mid + hi) / 2
    else:
    hi, mid = mid, (lo + mid) / 2
    fmid = func(mid)
    return mid

    def prob18_log(t):
    return bin_search(t, tetra)

  5. Graham said

    Sorry for the formatting error, I forgot a second quotation mark. Trying again:

    def tetra(n):
        return n * (n + 1) * (n + 2) / 6
    
    def prob18_linear(target):
        n = 1
        while tetra(n) < target:
            n += 1
        return n
    
    def bin_search(target, func):
        lo, hi = 1, 10
        while func(hi) < t:
            hi *= 10
        mid = (lo + hi) / 2
        fmid = func(mid)
        while fmid != target:
            if fmid < target:
                lo, mid = mid, (mid + hi) / 2
            else:
                hi, mid = mid, (lo + mid) / 2
            fmid = func(mid)
        return mid
    
    def prob18_log(t):
        return bin_search(t, tetra)
    
  6. Lautaro Pecile said
    def triangular():
        t = 0
        count = 1
        while True:
            t += count
            yield t
            count += 1
    
    def tetrahedral():
        tr = triangular()
        a = 0
        b = tr.next()
        while True:
            a += b
            yield a
            b = tr.next()
    
    t = tetrahedral()
    for i, n in enumerate(t):
        if n == 169179692512835000:
            print i
            break
    
  7. Mike said

    Python 3 added accumulate() to the itertools library, which makes implementing
    a linear search a snap.

    from itertools import accumulate, count
    
    tgt = 169179692512835000
    
    triangles = accumulate(count(1))
    tetrahedrals = accumulate(triangles)
    print(next(n for n,t in enumerate(tetrahedrals, 1) if t >= tgt))
    
    

    Inspired by Graham’s binary search, I tried to find a way to use ‘bisect’
    from the standard library. The Tetrahedrals class fakes being a list of
    tetrahedral numbers by calculating the numbers on the fly based on the index
    passed to __getitem__(). This is enough to get bisect() to work when using the
    lo and hi parameters.

    
    # binary search
    from bisect import bisect
    
    class Tetrahedrals:
        def __getitem__(self, n):
            return n * (n+1) * (n+2) // 6
    
    tetras = Tetrahedrals()
    
    tgt = 169179692512835000
    n,nn = 0,1
    while True:
        t = bisect(tetras, tgt, n, nn)
        if t < nn:
            print(t - 1)
            break
    
        n, nn = nn, 2*nn
    
    
  8. slabounty said

    In Ruby …

    def linear(target, f)
        n = 1
        while ( f.call(n) != target)
            n = n + 1
        end
        n
    end
    
    def binary(target, f)
        low, high = 1, 2
        while (f.call(high) < target) do high = high*2 end
        mid = (high + low) / 2
        while (fmid = f.call(mid)) != target do
            fmid < target ?  (low, mid = mid, (mid + high) / 2) : (high, mid = mid, (low + mid) / 2)
        end
        mid
    end
    
    tetrahedral = lambda { |n| n * (n + 1) * (n + 2) / 6 }
    
    1.upto(10) { |i| puts tetrahedral.call(i) }
    
    puts linear(169179692512835000, tetrahedral)
    puts binary(169179692512835000, tetrahedral)
    

    We use a lambda to create a tetrahedral function we can pass around specifically to the linear (very slow in this case) and binary (much faster).

  9. Axio said

    ;; With binary search. Quite naive.

    (define *number-of-balls* 169179692512835000)
    
    (define (tetra n)
      (/ (* n (1+ n) (+ 2 n)) 6))
    
    ;; bounds
    (define (prax balls)
      (let loop ((n 1) (last 1))
        (let ((tn (tetra n)))
          (cond
            ((= tn balls) n)
            ((< tn balls)
             (loop (* n 2) n))
            (else
              (dicho balls last n))))))
    
    ;; dichotomic search
    (define (dicho balls mi ma)
      (let* ((mid (/ (+ ma mi) 2))
             (mid-tn (tetra mid)))
        (cond
          ((= balls mid-tn) mid)
          ((> balls mid-tn) (dicho balls mid ma))
          (else (dicho balls mi mid)))) )
    
    (define (main)
      (pretty-print (prax *number-of-balls*)))
    
    (main)
    
  10. Axio said

    Optimised/more clean version:

    (define *number-of-balls* 169179692512835000)
    
    (define (tetra n)
      (/ (* n (1+ n) (+ 2 n)) 6))
    
    ;; dichotomic search
    (define (dicho balls mi ma)
      (pp (list mi ma))
      (let* ((mid (/ (+ ma mi) 2))
             (mid-tn (tetra mid)))
        (cond
          ((< (tetra ma) balls) (dicho balls ma (* 2 ma)))
          ((= balls mid-tn) mid)
          ((> balls mid-tn) (dicho balls mid ma))
          (else (dicho balls mi mid)))) )
    
    (define (main)
      (pretty-print (dicho *number-of-balls* 1 1)))
    
  11. Graham said

    @Mike: very elegant. I was looking for a way to let bisect do the search for me, but didn’t think beyond having a very large list.

  12. ardnew said

    all brute force solutions?

    here’s one using newton’s method:

    use strict;
    use warnings;
    
    my $TOLERANCE = 0.1;
    
    #
    # the n'th tetrahedral number: 
    #
    # T(n) = n(n + 1)(n + 2) / 6
    #      = (n^3 + 3n^2 + 2n) / 6
    #
    sub t
    {
      my $n = shift;
      return ($n * $n * $n + 3.0 * $n * $n + 2.0 * $n) / 6.0;
    }
    
    
    #
    # we are looking for the solution to:
    #
    # f(n) = T(n) = 169179692512835000
    #      = T(n) - 169179692512835000 = 0
    #
    sub f
    {
      my $n = shift;
      return t($n) - 169179692512835000.0;
                                                   
    }
    
    #
    # the derivative of f(n) with respect to n:
    #
    # f'(n) = 3n^2 + 6n + 2
    #
    sub df
    {
      my $n = shift;
      return 3.0 * $n * $n + 6.0 * $n + 2.0;
    }
    
    
    # 
    # perform newton's method to find the (unique!) root of f(n):
    #
    my $xn   = 1000000.0; # initial guess
    my $xnp1 = 0.0;
    
    do
    {
      $xnp1 = $xn;
      $xn = $xn - f($xn) / df($xn);  
      #print "$xn\n";
    }
    until (abs($xn - $xnp1) <= $TOLERANCE);
    
    
    #
    # naively round the solution
    #
    $xn = int($xn + 0.5);
    
    print "base number of balls = $xn\n";
    
  13. ardnew said

    and since the theme seems to be keeping it as concise as possible, good luck reading this one

    do{$b=$a;$a=$a-($a**3+3*$a**2+2*$a-1015078155077010000)/(3*$a**2+6*$a+2)}until(abs($a-$b)<=.1);
    
    print "number of balls = ".int($a+.5)."\n";
    
  14. THE G MAN said

    (x^3+3*x^2+2*x)/6 = 169179692512835000; one real root: 1005000

    Constant time. Thanks for playing.

  15. GDean said

    quick and dirty solution:

    n<-500 #initial length of vector
    #b<-20958500
    b<-169179692512835000 #number of balls in tetrahedron of interest

    for (y in 1:20) {

    tri<-1
    for (a in 2:n) {
    tri<-c(tri,max(tri)+a)
    }

    dtri<-1

    for (a in 2:n) {
    dtri<-c(dtri,max(dtri)+tri[a])
    }

    c<-1
    for (a in 1:n) {
    if (dtri[a]==b) {
    print(a)
    c<-a
    }
    }

    if (c==1) {
    n<-n*10
    }
    else {
    break
    }
    }

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

Follow

Get every new post delivered to your Inbox.

Join 612 other followers

%d bloggers like this: