I had an idea for a magic item for D&D not too long ago. It works something like this:

Pick some number of dice to roll, and if none of them are a 1, then you gain hit points equal to the sum of the dice. But if even one die is a 1, then you instead lose hit points equal to the sum of the dice.

So of course, being a fellow with some background in mathematics and a penchant for giving myself problems to solve, I decided to figure out how many dice you should roll in order to maximize your expected value.

This was not a trivial problem to solve. At first I just made a Python function—actually multiple functions because I had some helper functions in there too—that modeled the outcome of the player’s choice of how many dice to roll. The number of dice and the number of sides on each die were the two parameters, and I ran the function a few thousand times for each number of dice. That answered my question well enough, and maybe I could have tried finding a confidence interval or some such, but frankly I don’t remember that part of Math 407.

So I got to work trying to write a function to model the behavior of the item. My first attempt at modeling it came up with this:

(Forgive my sloppy syntactic style. This is pretty much verbatim how I had it in my notebook reproduced in LaTeX.)

In the above functions, the “1 epsilon” means “there is at least one 1 in the resulting rolls,” and likewise for “1 not-epsilon.” n is the number of dice rolled, and d is the number of faces on each die.

The second line (and thus also the sixth) ended up being horribly wrong, and ended up being a tricky problem to solve. I ended up making my first bit of solid progress by reducing the problem down to a simpler version of the same sort of issues: “If you flip N coins, and at least one of them is a heads, what is the expected value of the total number of heads?” This reduced the values down to just 0 or 1, which makes for simpler expressions than values 1 through 6.

So here, if n is the total number of heads, then the expected value of heads, E(n), is the following:

The summation is to account for all possible values of heads, with k being the number presently being considered. For each value of k (the number of heads being considered), there are n-choose-k ways of to get that, which I multiply by the number of heads to get the value of that count.

(The “k+0” should of course be just “k”, but the zero is just included in there to remind me that, while the tails have no value in this particular instance, in the original problem, there are no zero values. So the 0 is just a placeholder. Don’t mind that the function here is not simplified as much as it could be.)

So building on that simpler version, I can now get:

You can see that the zero got replaced by  above, since that is the average value of the dice that are not 1’s, multiplied by the number of dice I have that are not 1’s.

Finally this can be plugged back into the line five of the original attempt to get… drumroll please…

Or as it looks in my Python script I was using to test this:

def combo(n,r):
    return math.factorial(n)/(math.factorial(r)*math.factorial(n-r))
def expectedWith(d,n):
    sigma = [combo(n,k)*(d-1)**(n-k)*(k+(d+2)/float(2)*(n-k)) for k in range(1,n+1)]
    return sum(sigma) / float(d**n-(d-1)**n)
def probabilityWith(d,n):
    return 1-((d-1)/float(d))**n
def probabilityWithout(d,n):
    return ((d-1)/float(d))**n
def expectedWithout(d,n):
    return (d+2)/float(2)*n
def expectedTotal(d,n):
    return expectedWithout(d,n)*probabilityWithout(d,n) - expectedWith(d,n)*probabilityWith(d,n)

That’s even uglier to look at here than the syntax I used for my math. I promise it looks better in Sublime Text 2.

Anyway, I feel reasonably certain I didn’t fuck this up. And, with six-sided dice, the answer was two dice.