blogstrapping

Fast Die Roll Permutations Still Elusive

In How Can I Find the Number of Permutations per Sum?, I discussed a little programming problem I was trying to solve:

I started trying to automate the process of measuring the probabilities of achieving particular rolls. I wrote some code in Ruby to do so, and the code works -- as long as system resources hold up. The problem is that my code stores numbers in arrays and for any nontrivial example combination of dice the number of possible rolls quickly outstrips my computer's capacity to store all those numbers in arrays that rapidly expand to ridiculous sizes.

Since publishing that essay at blogstrapping, my little calculation problem has garnered some attention from others -- some of them people I know personally, some who I know somewhat impersonally, and others that I do not know at all. Fellow TechRepublic contributor Justin James recently wrote a program in IronRuby to take a whack at solving my problem, and wrote an article about the exercise for TechRepublic. His article, My first IronRuby application, includes source code for his solution. Unfortunately, something awful happened to the formatting, and I am not feeling the urge to try untangling it right now. Luckily, he sent me some code prior to publication that I think is probably pretty close to what ended up in the article:

def calculate (iteration, low, high, currentsum, output)
    if (iteration == 1)
        low.upto(high) do |value|
            newsum = currentsum + value
            output[newsum] += 1
        end
    else
        low.upto(high) do |value|
            calculate(iteration - 1, low, high, currentsum + value, output)
        end
    end

    return output
end

diceInput = ARGV[0].to_i
lowInput = ARGV[1].to_i
highInput = ARGV[2].to_i

if (diceInput < 1)
    puts "You must use at least one dice."
    return
end

initResults = Hash.new

(lowInput * diceInput).upto(highInput * diceInput) do |value|
    initResults[value] = 0
end

results = calculate(diceInput, lowInput, highInput, 0, initResults)

results.each do |result|
    puts "#{result}"
end

puts "Press Enter to quit..."
gets

For the record, I think his "Press Enter" thing at the end might be an artifact of the way he wrote the code for IronRuby; I just delete that crap when I run it on my FreeBSD laptop from the shell prompt.

Justin's result is certainly faster than my first, naively brute-force, attempt. It just is not anywhere near enough faster. It starts getting unusably slow around the point where you are counting permutations for fifteen dice or so. This is because, while parts of the program are better optimized than mine, it suffers the same problem: it takes a brute force recursive approach to calculating permutations, where at least one operation has to be performed for half the total possible permutations.

By way of a friend, someone suggested using Pascal's Triangle as the basis of a function that would generate the correct number of permutations for each possible result of any given die roll. I beat my head against that problem for a little while, then gave up for a while. In the process, though, I ended up writing two or three different versions of a Pascal's Triangle class for Ruby. I expect to polish it up a little and stick it in BitBucket soonish.

The problem with the Pascal's Triangle approach was not speed. It was blazing fast. It was not that I could not generate the correct numbers for many cases (such as three six-sided dice). Before explaining the problem, I will explain the approach.

Pascal's Triangle is an infinite set of numbers canonically generated by starting with a one, then a row containing two ones below it; for every successive row, add together every pair of numbers from the preceding row, and add a one to each end. Thus, a twelve-row Pascal's Triangle looks like this:

                            1

                         1    1

                       1    2    1

                    1    3    3    1

                  1    4    6    4    1

               1    5   10   10    5    1

             1    6   15   20   15    6    1

          1    7   21   35   35   21    7    1

        1    8   28   56   70   56   28    8    1

     1    9   36   84  126  126   84   36    9    1

   1   10   45  120  210  252  210  120   45   10    1

1   11   55  165  330  462  462  330  165   55   11    1

The basic insight that makes this look like the source of a solution to the dice roll permutations problem is the fact that, for any die type with N faces, a rightward-sloping diagonal starting on the row corresponding with the number of dice rolled contains a series of numbers that matches the number of permutations for the first N possible rolls. For example, if you roll 3d6, count down to the third row (which reads 1 2 1); starting with the leftmost 1, follow the rightward-sloping diagonal and you will find numbers that match the number of permutations that will produce the first six results (3, 4, 5, 6, 7, 8). Those numbers of permutations are 1, 3, 6, 10, 15, and 21.

Past that point, of course, things change. For a roll of 3d6, the permutations look like the following:

 3       1
 4       3
 5       6
 6      10
 7      15
 8      21
 9      25
10      27
11      27
12      25
13      21
14      15
15      10
16       6
17       3
18       1

The curve is always symmetrical, so you only have to calculate the first half, then you can mirror that set of results. The problem is calculating the numbers in the first half (1, 3, 6, 10, 15, 21, 25, 27) that come after the first six: 25 and 27. These numbers do not match the next two numbers after 21 in the series taken from Pascal's Triangle: 28 and 36.

It turned out to be reasonably easy to figure out a function that will work for two or three dice. It works for two because, when you get through a number of results equal to the number of faces on a die, you're done with the first half of the "curve", which means you do not actually need to use anything other than the diagonal series from Pascal's Triangle. Unfortunately, it does not work for four. I screwed around looking for a generalized pattern that works for arbitrary numbers of dice, and failed to find one. At least I developed a better understanding of Pascal's Triangle, and maybe I will find a way to do this with Pascal's Triangle eventually. On the other hand, maybe there is no way to use Pascal's Triangle for this.

I received some feedback from various readers of blogstrapping via the contact form. Interesting suggestions included:

The first was an interesting approach that differed from both mine and Justin's, but looked like yet another (reduced) brute force approach. Maybe I misread it; I will look at it in more depth later. I am still burned out on this after fighting with the Pascal's Triangle approach for a while.

The second included two specific links in the sequence database:

I have not had the inspiration to look at the two pages from the sequence database in any depth yet, but I'm listing them here in part to provide myself an easy sort of "bookmark" to remind me to check into it later.

Then . . . someone Justin knows -- Fancisco Blanco-Silva -- produced an implementation in comments following his TechRepublic article. 100 100 in about a minute is the comment that holds the source. The claim is that it generates the numbers of permutations for all the possible outcomes of a 100d100 roll in about a minute -- later specifically described as seventy seconds.

Of course, my results showed about eleven or twelve minutes using MRI for Ruby 1.8.7, and about six minutes using YARV for Ruby 1.9.2 instead. Maybe IronRuby is faster at executing this specific algorithm than MRI and YARV, and maybe Francisco is using IronRuby; maybe Francisco is using a computer that totally outclasses my ThinkPad T60; maybe the stars just aren't aligning for me. All I know is that even seventy seconds is longer than I would like, and it will still end up taking longer periods of time very quickly as you increase the number of virtual dice.

Francisco wrote up an explanation of his approach and posted it to the Web: Unusual dice. He makes use of a pattern I had intuited while thinking about the problem early on, where listing the permutations for one die in a row, for a number of rows equal to the number of dice, each row shifted one column to the right from the placement of the row above it, provides the basis for simple addition to generate the needed permutation totals.

It is faster than my really naive approach; it is faster than Justin's; it is still too slow (and resource heavy) for practical use in a relatively unconstrained context -- like providing a Webpage that people can use to generate their own permutation counts. Permutations for all results in 100d100 can at least be calculated within a period of time that is by some lights reasonable, with Francisco's approach. If my only need was to produce a set of permutations for 100d100, I would have been done the first time I spent more than ten minutes waiting for the program to run. What I want, though, is a generally usable application, so it kinda falls short. Even one minute is a bit much for something like that.

I did a minimal rewrite of Francisco's program -- renaming variables, replacing a for loop and an Array#cycle(1) message with Array#each iterators, and monkeypatching the Array class for aesthetic purposes -- for purposes of making it a little easier for me to follow. Performance does not really change with my rewrite. Francisco asked me to share my version of his code (I think), but code sample display is not one of the strengths of TechRepublic's discussion software. I decided to post it here, instead of there:

#!/usr/bin/env ruby

class Array
  def sum
    self.inject(0) {|total,num| total + num }
  end

  def fprint
    print self.join(' '), "\n"
  end
end

def build_struct(width, height)
  Array.new(width) { Array.new(height, 0) }
end

def count_perms(faces, perms, dice)
  if dice == 1
    return perms
  else
    perm_struct = build_struct(faces, faces + perms.length - 1)
    (0...faces).each do |key|
      perm_struct[key][key...perms.length+key] = perms
    end

    mirror_struct = perm_struct.transpose
    count_perms(
      faces,
      Array.new(faces + perms.length - 1) {|x| mirror_struct[x].sum },
      dice - 1
    )
  end
end

def get_permutations(dice, faces)
  perms = Array.new faces, 1
  if dice == 1
    return perms
  else
    return count_perms(faces, perms, dice)
  end
end

if $0 == __FILE__
  dice  = ARGV[0].to_i
  faces = ARGV[1].to_i

  get_permutations(dice, faces).fprint
end

That's it for now.

To Be Continued . . .