Comparing the Probability of Parallel Multiple Dice Rolls

I've been working on a battle game with my son. We've got an idea for how combat could work, but I'm struggling to work out a way to calculate probability.

It's similar to how combat works in "Risk":

• The Attacker has x * 6 sided die.
• The Defender has y * 6 sided die.

Both roll their die & the dice are compared in parallel from highest roll to lowest roll. Example:

Attacker (rolls 4 dice): 6 4 3 2 
Defender (rolls 3 dice): 5 5 3

In the above example:

1. The Defender wins one dice pair.
2. The Attacker wins one dice pair.
3. The remaining dice pair is a draw.

With our proposed rules the Attacker & Defender lose one dice for each pair they lose & roll again:

Attacker (now rolls 3 dice): 6 3 1  
Defender (now rolls 2 dice): 4 2 

In this example the Defender loses both pairs, so they have been defeated.

I'm trying to work out how to compare the probabilities of different dice numbers - i.e. 4 dice vs 3 dice or 5 dice vs 2 dice.

Being a software developer I've thought of programmatically running various Attacker/ Defender ratios thousands (or millions) of times & recording the odds from the outcome.

However, this doesn't feel remotely elegant & I'm wondering if there's a formula I could use instead? I tend to use Google Sheets for these sorts of calculations.

Would it be possible to calculate the odds for this approach via a formula, or is it too complex?

An efficient pool evaluation algorithm

In general, if the function that scores a pool (or in this case, a pair of pools) is allowed to depend on the entire ordered rolled sequence, then the function would have to be evaluated on all possible such sequences, the number of which is exponential in the number of dice. Even if we restrict ourselves to functions that only depend on unordered rolled sequences (or equivalently, the sequence after sorting each pool), the number of such sequences grows as a high order polynomial (this can be derived from e.g. comments on this answer to another question).

However, it turns out that if we can phrase the pool scoring function as a state transition over tuples (outcome, number of dice in each pool that rolled that outcome) with not too many states, we can reduce the time complexity to a polynomial of reasonable order. Here's a fuller explanation of the algorithm. This RISK-like scoring can indeed be expressed like this.

Expressing a single stage

I implemented this approach as part of my hdroller Python library. You can try this computation in your browser using this JupyterLite notebook.

The first step is to express a single stage of a RISK-like pool evaluation as a state transition function as defined above. Given the number of dice on each side, this computes the distribution of the number of remaining dice on each side after a single stage.

class EvalRiskAttrition(hdroller.EvalPool):
    def next_state(self, state, outcome, a, b):
        if state is None:
            score_a, score_b, advantage = 0, 0, 0
        else:
            score_a, score_b, advantage = state
        # Advantage is the number of unpaired dice that rolled a previous (higher) number.
        # If positive, it favors side A, otherwise it favors side B.
        # We pair them off with newly-rolled dice of the disadvantaged side.
        if advantage > 0:
            score_a += min(b, advantage)
        elif advantage < 0:
            score_b += min(a, -advantage)
        advantage += a - b
        return score_a, score_b, advantage
    
    def final_outcome(self, final_state, pool_a, pool_b):
        score_a, score_b, advantage = final_state
        if score_a == 0 and score_b == 0:
            # No change. Eliminate this outcome to prevent infinite looping.
            # This is equivalent to rerolling the contest until at least one die is removed.
            return hdroller.Reroll
        # Each side loses dice equal to the other's hits.
        # The result is the number of remaining dice on each side.
        return pool_a.num_dice() - score_b, pool_b.num_dice() - score_a
    
    def direction(self, *_):
        # See outcomes in descending order.
        return -1

eval_risk = EvalRiskAttrition().bind_dice(hdroller.d6, hdroller.d6)

Finding the fixed point

Since we reroll each stage until at least one die is eliminated, the number of stages is bounded by the total number of dice. Therefore we can then repeatedly apply this attrition until everything reaches an absorbing state. We could compute the maximum number of stages explicitly, or just look for a fixed point.

def risk_attrition(a, b):
    if a == 0 or b == 0:
        # If one side has run out of dice, no more rolling is necessary.
        return a, b
    else:
        # Otherwise, run the contest.
        return eval_risk(a, b)

# 4 dice vs. 3 dice.
a = 4
b = 3
# Construct a die that always rolls the tuple (4, 3).
# Then, apply the risk_attrition function recursively until reaching a fixed point.
result = hdroller.Die((a, b)).sub(risk_attrition, max_depth=None, denominator_method='reduce')
# The result is how many dice are remaining at the end.
# The loser has 0 dice, and the winner has 1 or more dice.
print(result)

Denominator: 350190883979307611136000

If you just want to know the winner:

print(result.sub(lambda a, b: 'a' if a > 0 else 'b').reduce())

Denominator: 2244741412001587200

This can compute a dozen dice on each side without much issue, though the denominator can get quite massive.