carlo.m
(Monte Carlo simulation of 3 dice problem)


The recreational mathematician and chess master Sam Loyd was also one of the most famous authors of puzzles of his day (1841 - 1911). Among the thousands of his puzzle creations is a game he called the Carnival Dice Game which works as follows: You bet 1 dollar to play which allows you to roll 3 dice: The question is: Is this a good bet to make?

Actually, this problem is not that difficult to solve if you have any experience with counting problems. Consider that there are 216 different possible outcomes (63) and each of these outcomes is equally likely. So let's assume we play the game 216 times and during that time we see each of the possible outcomes once. Now we can calculate our winnings: So our total winnings for the 216 bets is $199 (4 + 45 + 150). But we had to pay 1 dollar to make each of these bets so we end up with a net loss of $17. So on average, we will lose about .0787037 dollars (17/217) or nearly 8 cents on each bet. (This calculation of the average loss per bet is repeated inside carlo.m in lines 65 thru 69 where the last of those lines also puts the result into the third trace of the plot (using "expected value" for its trace ID and depicted as a horizontal dotted white line).

But suppose we didn't know how to make this calculation (or more likely we were dealing with a much more complicated problem where the calculations are intractable). Faced with this problem, can we still get an approximate answer to the question? Well yes, we can start rolling 3 dice and after each roll, we add the winnings to our total (as well as subtracting the dollar we pay for each bet). At first, this will yield a poor approximation of the true answer but the more times we repeat the process the better our answer becomes.

This idea is probably old, but the first time that we know about where a computer was used to execute this type of simulation was during the Manhattan project when they programmed the ENIAC computer to estimate the average distance a neutron would travel before it collided with an atomic nucleus and how much energy would be produced from this collision. Until this time, the method did not have a name, so one of the physicists on the project suggested using the name of a casino (the Monte Carlo in Monaco) near where he grew up. This seemed appropriate since games of chance are central to the business of a casino.

Although we don't need the Monte Carlo method to solve Sam Loyd's dice game puzzle, knowing the answer before running the Monte Carlo simulation allows us to observe how well the method estimates the true answer.

There are three different ways you can start this application:

carlo      When called with no arguments, the simulation is set up but no bets occur until you click a button.
carlo(25)
carlo 25  		Sets up simulation & makes 25 bets and then stops.
(This is an example of course ... you can use any positive number)
carlo(0)
carlo 0
carlo('run')
carlo run 		Sets up the simulation & makes bets continuously until you click stop.

Suppose we start the application using the first way (i.e. no arguments). Then we can simulate making one bet at a time by clicking the "1 bet" button. Every time we click, it rolls the 3 dice and we see the faces of the dice in the green squares above the plot. It calculates how much you would win or lose based on the number of sixes and displays this value to the right of the 3 dice and then appends the total earnings so far to the "Accumulated earnings" (purple trace) on the left hand axis. It also computes "Earnings per bet" which it appends to the green trace on the right hand axis. Doing one bet at a time helps us see what is happening in the application as we continue to play the game, but this is too slow to make any significant progress in estimating the Earnings per bet expected value. For this, we need to use the "run" button to execute the bets repeatedly.

The speed slider


The speed slider is used to control how fast we make the bets once the run button is pressed. The negative numbers (when the slider is close to the left edge) are for very slow speeds. A setting of -20 (for example) means that a delay of 20/100 seconds (or two-tenths of a second) is added after each bet is made. After clicking run, you will see a report below the plot (to the left of the Xaxis label) that will show approximately 5 bets/sec. Moving the slider all the way to the left (-100) gives a delay of 1 second and the readout will show about 1 bet/sec. This is too slow to make much progress towards approximating the Earnings per bet, but it allows us to watch the dice being rolled and the earnings accumulate. If you set the slider to zero, the plot will be updated after each bet and then the next bet will be made with no delay. On my (2020 vintage) computer, it reports that about 50 bets are being made every second. (You may see a faster or slower rate depending on the speed of your computer.) This is still pretty slow, but we can go faster. For example, if we move the slider to 900 (all the way to the right, 900 bets will be made before updating the display. This is much faster since updating the display takes a lot longer than making a bet and accumulating the results. On my computer when the slider is set to 900 it reports that it is making about 45 thousand bets every second. Even though that is the top end of the slider, we can go even faster than that by typing in a number into the slider edit box. For example, typing in "5E5" results in a speed of about 10 million bets per second. That's about as fast as we can go, since when we are making 500,000 bets for every display update, the display updates are taking up a small fraction of the cpu resources. Typing in an even bigger number will slow the display update rate, but it won't significantly increase the number of bets made per second.



Here I have set the speed to a fairly slow setting (10 bets per display update) so I could watch the progression and run it for 10 seconds or so. Notice that the green line is starting to approach the dotted line (at the expected value of -.0787 dollars) but is still not that close. We can see the losses steadily building by the decline of the purple trace (accumulated earnings) into negative territory. But as expected for a random process this decline is somewhat jagged with many small peaks and valleys.



Now let's set a fairly high bet rate by typing the number 9000 into the speed slider, and run it for around another 20 seconds.

Note that the Xaxis label shows that the Xaxis units are now in millions of bets instead of just the number of bets in the previous plot. (The plt('metricp') function makes this easy.) Also, note that the Yaxis units are now in thousands of dollars instead of just dollars in the previous plot.

The green trace is now getting very close to the true expected value shown by the dotted white line. Note that the right hand axis has zoomed in by about a factor of 40 compared with the previous plot. This happens because some of the less accurate early data gets dropped off due to data decimation (which is discussed later). This allows us to see the difference between the estimated average earnings per bet and the exact computed value. Because we are now showing so much data, the purple trace (accumulated earnings) looks like a straight line. However, if we could zoom in to get a detailed view of any portion of this trace it would look just as jagged as it does in the previous figure.

This application is unique in that the plot is updated every time the number of points plotted for each trace increases by one. By doing this continuously and rapidly we would soon be plotting lines containing millions of points which would start to slow down the update rate as well as the speed of the simulation (bets per second). To combat that problem, whenever the program sees that the lines contain more than 6000 points, it decimates the data by a factor of 10 (i.e. nine points are removed between each successive point that we are keeping). Because we are starting with so many points you can't tell on the display that we are not plotting all the original points. (After all, there are only a fixed number of pixels in the plotting area, so no more than about 600 of the data points can actually be plotted.)

By running this program from demoplt we can see that the carlo application consists of 87 lines of code (not counting the comments). Perhaps you will agree that this is pretty short considering the complexity of the simulation and the plotting process as we have described above.

Copyright © 2024
Paul Mennen