What is _wrong_ with Python? a.k.a. PoR probabilities

[Oedipussy Rex]

So. What's the difference between a d4 and a d6 that hasn't rolled a 5 or 6?

I decided I wanted to figure out the probability of getting N successes on M dice to beat TN using PoR mechanics ( standard comparative dice pool, ala White Wolf, with varying dice sizes where dice <= TN in size can be added together). When the dice are larger than TN, it's no big deal. It's not Binomial Distribution easy (all dice are the same size), but easy. I wrote a python script that calculated all probabilities for target numbers from 1 to 11 all possible dice pools from {2d4} to {11d12} where every die is greater than the target number in less than 15 minutes. It's not an especially efficient script, but I took measures to eliminate a lot of overhead. For instance, a pool is a set so order doesn't matter. Because order doesn't matter, we're only dealing with combinations, not permutations. Also, because order doesn't matter, there's no reason to not have the dice in order by size. Here's my code for die sequence:

Code:

if pool[0:MAX] == [ 12 ]*MAX:    #pool is a list of dice indicated by size (4, 6, 8, 10, 12), Max is the number of dice in the pool
       break    #All combinations have been used, go to the next TN
else:
       for idx in xrange (0,MAX):
         if idx == MAX-1:    #pool[0] through pool[idx-2] are unchanged
           pool[idx] += 2    #the next larger die
         elif pool[(idx+1):MAX] == [pool[idx]]*(MAX-idx-1):    #all dice from idx on are the same size
           pool[idx] += 2
           for j in xrange(idx+1,MAX):
             pool[j] = 4
           break   #idx    #only one die changes, so stop looking

The pools iterated thusly:

Code:

{4, 4, 4}
{6, 4, 4}
{6, 6, 4}
{6, 6, 6}
{8, 4, 4}
{8, 6, 4}
{8, 6, 6}
{8, 8, 4}
{8, 8, 6}
{8, 8, 8}
{10, 4, 4}
etc.

Now, if the target number is 7, then {8, 4, 4}, {8, 6, 4}, and {8, 6, 6} are all equivalent to {8} and each were calculate as such. Duplicate processing. Wasteful and inefficient. I should have set it up so that the dice ranged not from 4, but from MIN = max([(TN/2+1)*2, 4]). (For TN=7, MIN is 8, TN=8, MIN is 10, just how we want it.) But like I said earlier, working this non-binomial distribution is easy. You just need to know the odds for each Bernoulli Trial (single die) and figure out the 2^n combinations of successes. Easy. Done in less than 15 minutes.

Adding dice to exceed a target number, on the other hand....

Not only do you need to iterate through every pool that contains dice smaller than TN, you then need to iterate the dice rolls, then you need to see how many successes you can get from each roll of the pool. Now, to help speed things along, if the sum of all the dice is less than the target number, there aren't any successes. If the sum is less than twice the target number, there can only be one success. But if the number of dice is greater than 3 and the sum of the dice is greater than 2TN, well, now you need to do some calculating. Here's my recursive function:

Code:

def cycle( rolls, TN ):     # rolls is a user Class that contains a list where the index indicates\
                           # the result of a die and the value is the number of dice that rolled it.
   rollsum = rolls.rollsum()
   if rollsum < TN:
     return 0
   elif rollsum < 2 * TN:
     return 1
   else:
     large = rolls.newlarge()
     if rolls.rolls[ TN - large ]:     # Here TN is actually matched not exceeded. Makes for cleaner code to get the same result.
       rolls.rolls[ TN - large ] += -1
       rolls.rolls[large] += -1
       return 1 + cycle( rolls, TN )
     else:
       sum = large
       rolls.rolls[large] += -1
       while sum < TN:
         small = rolls.newsmall()
         sum += small
         rolls.rolls[small] += -1
       return 1 + cycle( rolls, TN )

This is a greedy algorithm that starts with the largest number and, if its complement doesn't exist, adds the smallest numbers until TN is matched or exceeded. I don't have proof that this method always gives the correct result, but anecdotally, I couldn't create a situation where it didn't.

And yes, I'm aware of the dangers of recursive functions in terms of time and memory. However, the number of times the function is called is at most ceiling(MAX/2). For a pool of 11 dice that's 6 times. So, the recursive function isn't the problem. The problem is if you have a pool of 8 d10s, you're looking at 100,000,000 different rolls. You know that fifteen minutes I mentioned earlier? It takes longer than than to calculate 100,000,000 die rolls to beat 10. Then we get to do it again for 11. And 12.

Which brings us back to my original question: What's the difference between a d4 and a d6 that hasn't rolled a 4 or 5? The answer: nothing. let's look at the roll combinations between the two for {2d4} and {d4,d6}.

Code:

d4     d4           d4     d6 that hasn't rolled 5 or 6
0      0            0      0    dice that roll 1 cannot be used, thus 0.
0      2            0      2
0      3            0      3
0      4            0      4
2      0            2      0
2      2            2      2
2      3            2      3
2      4            2      4
3      0            3      0
3      2            3      2
3      3            3      3
3      4            3      4
4      0            4      0
4      2            4      2
4      3            4      3
4      4            4      4

Looks fairly identical to me. So after calculating {2d4}, I only need to calculate

Code:

5    0
5    2
5    3
5    4
6    0
6    2
6    3
6    4

Remember that 100,000,000 calculations? Now it's only 20,000,000, the other 80,000,000 having already been done. Now that I've finally figured that out all I have to do is actually code it.