little math/programming problems April 2023
04-28-2023, 11:55 AM
Post: #1
 pier4r Senior Member Posts: 2,248 Joined: Nov 2014
little math/programming problems April 2023
I have a problem someone could probably simulate faster than I can work it out. I have a set of dice with values: 0, 0, 0, 0, 2, 3. You roll the dice, add up the total, then roll that number of dice. And repeat. If I start with five dice, how long will the game last?

Source: https://twitter.com/jamesgrime/status/16...46721?s=20 (contains solutions that spoil the challenge a bit)

Wikis are great, Contribute :)
04-28-2023, 10:16 PM
Post: #2
 David Hayden Senior Member Posts: 461 Joined: Dec 2013
RE: little math/programming problems April 2023
Simulating the experiment 100,000 times, I get an average of 5.49893 rolls with 5 dice. Here are the results with different numbers of starting dice:
Code:
 Initial # dice    Average rolls 1                 2.2266 2                 3.25986 3                 4.10384 4                 4.84957 5                 5.49893 6                 6.08402
Someone smarter than me will have to figure out the exact values.

Dave
04-29-2023, 06:50 PM
Post: #3
 Allen Member Posts: 224 Joined: Aug 2014
RE: little math/programming problems April 2023
Maybe possible with the generating function and some ugly recursions?

$$(z^3 + z^2 + 4)^n$$

where $$n$$ is the number of dice.
$$f(5) = z^{15} + 5z^{14} + 10z^{13} + 30z^{12} + 85z^{11} + 121z^{10} + 240z^9 + 500z^8 + 480z^7 + 800z^6 + 1280z^5 + 640z^4 + 1280z^3 + 1280z^2 + 1024$$

(Interesting it's not possible to roll a sum equal to 1 no matter how many dice you use!)

Since the number of dice changes based on the previous roll, you'd need to recursively add up the probabilities for every possible branch in order to calculate an (more?) exact answer - multiplying the conditional probability of each branch times the total chances of arriving at that. Using the initial state of $$f(5)$$ and creating a markov chain 24 layers deep you can cover 99% of the expected cases.

Among 10M games with 5 initial dice I got this estimation:
• 5.490012 expected rolls on average
• mode was 2 rolls
• the longest was 81 rolls!!
• 13% of the games ended on the first roll.
• 99% of the games ended in 24 rolls or less.

I also tried with 1000 to 1M dice at a time.. the average number of rolls asymptotically approaches something around $$4.32 * ln(n)$$ with 1000-dice games ending in 32 rolls on average.

17bii | 32s | 32sii | 41c | 41cv | 41cx | 42s | 48g | 48g+ | 48gx | 50g | 30b

04-30-2023, 01:32 AM (This post was last modified: 04-30-2023 01:41 AM by Allen.)
Post: #4
 Allen Member Posts: 224 Joined: Aug 2014
RE: little math/programming problems April 2023
I believe the first few terms are now calculated (hopefully correct). Starting with 5 dice, the chances of the game terminating on exactly x rolls is:

Code:
 x  approx                exact   0  0.0                    0 1  0.13168724279835392   32 / 3^5 2  0.1762589048590848    614576800 / 3^20 3  0.14345296021217352   1477716325360817158616745850400 / 3^65 4  0.11119048831439719   1127935854140837598232937023804660987401971125340384 / 3^109 ...

The total number of states for depth grows quickly to [1, 15, 375, 17655, 1469850, , ...] After that things get out of hand

Here's the difference between the Monte Carlo sampling and the exact calculation:

Code:
 halting count for 10,000,000 rounds sampling Monte Carlo starting with 5 dice rolls    count    Monte Carlo    Empirical    Abs Error    Percent error 1    1317393    0.1317393      0.131687243    -5.20572E-05    -0.039530938 2    1761181    0.1761181      0.176258905     0.000140805     0.079885246 3    1434023    0.1434023      0.14345296      5.06602E-05     0.03531486 4    1112492    0.1112492      0.111190488    -5.87117E-05    -0.052802795

10M rolls Total Expected Value 5.4931065

Code:
 rolls    count    Monte Carlo    PDF    cumulative 1    1317393    0.1317393    0.1317393    0.1317393 2    1761181    0.1761181    0.3522362    0.4839755 3    1434023    0.1434023    0.4302069    0.9141824 4    1112492    0.1112492    0.4449968    1.3591792 5    862577    0.0862577    0.4312885    1.7904677 6    673073    0.0673073    0.4038438    2.1943115 7    531444    0.0531444    0.3720108    2.5663223 8    424292    0.0424292    0.3394336    2.9057559 9    339510    0.033951    0.305559    3.2113149 10    274634    0.0274634    0.274634    3.4859489 11    223553    0.0223553    0.2459083    3.7318572 12    183459    0.0183459    0.2201508    3.952008 13    149193    0.0149193    0.1939509    4.1459589 14    123336    0.0123336    0.1726704    4.3186293 15    100750    0.010075    0.151125    4.4697543 16    83419    0.0083419    0.1334704    4.6032247 17    68726    0.0068726    0.1168342    4.7200589 18    57454    0.0057454    0.1034172    4.8234761 19    47370    0.004737    0.090003    4.9134791 20    39204    0.0039204    0.078408    4.9918871 21    32289    0.0032289    0.0678069    5.059694 22    26840    0.002684    0.059048    5.118742 23    22160    0.002216    0.050968    5.16971 24    18664    0.0018664    0.0447936    5.2145036 25    15493    0.0015493    0.0387325    5.2532361 26    12861    0.0012861    0.0334386    5.2866747 27    10782    0.0010782    0.0291114    5.3157861 28    9111    0.0009111    0.0255108    5.3412969 29    7672    0.0007672    0.0222488    5.3635457 30    6240    0.000624    0.01872    5.3822657 31    5192    0.0005192    0.0160952    5.3983609 32    4240    0.000424    0.013568    5.4119289 33    3544    0.0003544    0.0116952    5.4236241 34    2980    0.000298    0.010132    5.4337561 35    2461    0.0002461    0.0086135    5.4423696 36    2093    0.0002093    0.0075348    5.4499044 37    1686    0.0001686    0.0062382    5.4561426 38    1459    0.0001459    0.0055442    5.4616868 39    1204    0.0001204    0.0046956    5.4663824 40    1034    0.0001034    0.004136    5.4705184 41    809    0.0000809    0.0033169    5.4738353 42    715    0.0000715    0.003003    5.4768383 43    565    0.0000565    0.0024295    5.4792678 44    448    0.0000448    0.0019712    5.481239 45    396    0.0000396    0.001782    5.483021 46    326    0.0000326    0.0014996    5.4845206 47    311    0.0000311    0.0014617    5.4859823 48    235    0.0000235    0.001128    5.4871103 49    172    0.0000172    0.0008428    5.4879531 50    142    0.0000142    0.00071    5.4886631 51    135    0.0000135    0.0006885    5.4893516 52    102    0.0000102    0.0005304    5.489882 53    73    0.0000073    0.0003869    5.4902689 54    78    0.0000078    0.0004212    5.4906901 55    74    0.0000074    0.000407    5.4910971 56    58    0.0000058    0.0003248    5.4914219 57    45    0.0000045    0.0002565    5.4916784 58    40    0.000004    0.000232    5.4919104 59    34    0.0000034    0.0002006    5.492111 60    22    0.0000022    0.000132    5.492243 61    22    0.0000022    0.0001342    5.4923772 62    20    0.000002    0.000124    5.4925012 63    21    0.0000021    0.0001323    5.4926335 64    16    0.0000016    0.0001024    5.4927359 65    5    0.0000005    0.0000325    5.4927684 66    9    0.0000009    0.0000594    5.4928278 67    5    0.0000005    0.0000335    5.4928613 68    6    0.0000006    0.0000408    5.4929021 69    4    0.0000004    0.0000276    5.4929297 70    3    0.0000003    0.000021    5.4929507 71    5    0.0000005    0.0000355    5.4929862 72    2    0.0000002    0.0000144    5.4930006 73    2    0.0000002    0.0000146    5.4930152 74    5    0.0000005    0.000037    5.4930522 76    2    0.0000002    0.0000152    5.4930674 77    2    0.0000002    0.0000154    5.4930828 78    2    0.0000002    0.0000156    5.4930984 81    1    0.0000001    0.0000081    5.4931065

17bii | 32s | 32sii | 41c | 41cv | 41cx | 42s | 48g | 48g+ | 48gx | 50g | 30b

05-02-2023, 08:03 PM
Post: #5
 David Hayden Senior Member Posts: 461 Joined: Dec 2013
RE: little math/programming problems April 2023
(04-29-2023 06:50 PM)Allen Wrote:
• 5.490012 expected rolls on average
I'm surprised that we came up with such different answers. With 10 million trials, I got an average of 5.49603. Here's my C++ code. Does anyone spot a problem? Maybe this is a case of rand() showing it's limitations?
Code:
#include<iostream>                               //Required for cin, cout #include<cmath> using namespace std; int main() {     int sides[6] = {0,0,0,0,2,3};     int N, D;     int totalRolls = 0;     cout << "Enter starting number of dice and number of simulations: ";     cin >> D >> N;     for (int i=0; i<N; ++i) {         int sum = D, count=0;         int dice;         for (count=0; (dice = sum); ++count) {  // note that (dice=sum) is assignment, not comparison             for (int i=sum=0; i<dice; ++i) {                 int side = sides[rand() % 6];                 sum += side;             }         }         totalRolls += count;     }     cout << "Average of " << (double)totalRolls/N << " rolls starting with " << D << " dice\n"; }

\$ ./foo
Enter starting number of dice and number of simulations: 5 10000000
Average of 5.49603 rolls starting with 5 dice
05-03-2023, 07:10 AM (This post was last modified: 05-03-2023 07:11 AM by ThomasF.)
Post: #6
 ThomasF Member Posts: 229 Joined: Sep 2016
RE: little math/programming problems April 2023
(05-02-2023 08:03 PM)David Hayden Wrote:  Here's my C++ code. Does anyone spot a problem? Maybe this is a case of rand() showing it's limitations?

Hi David!

I assume this is due to the limitations of rand().
With different seeds I get the following result (ie your code but added seed initialization).
Code:
  srand (time(NULL));

Code:
Average of 5.48456 rolls starting with 5 dice Average of 5.48837 rolls starting with 5 dice Average of 5.48826 rolls starting with 5 dice

Using another random generator (PCG), I get the following result:
Code:
Average of 5.49419 rolls starting with 5 dice Average of 5.49255 rolls starting with 5 dice Average of 5.49275 rolls starting with 5 dice

Edit: All runs used 5 dice and 10M simulations.

Cheers,
Thomas

[35/45/55/65/67/97/80 21/25/29C 31E/32E/33E|C/34C/38E 41C|CV|CX 71B 10C/11C/12C/15C|CE/16C 32S|SII/42S 28C|S 48GX/49G/50G 35S 41X]
05-03-2023, 12:18 PM (This post was last modified: 05-04-2023 12:13 PM by Albert Chan.)
Post: #7
 Albert Chan Senior Member Posts: 2,679 Joined: Jul 2018
RE: little math/programming problems April 2023
(05-02-2023 08:03 PM)David Hayden Wrote:  Does anyone spot a problem? Maybe this is a case of rand() showing it's limitations?

rand(), and many other random number generators, have weaker random low bits.
We should use random high bits instead.

< int side = sides[rand() % 6];
> int side = sides[rand() * 6ull / (RANDMAX+1u)];

Bonus: (RAND_MAX+1u) is likely pow-of-2, division should compiled as fast bits shift.
Note: we can speed up a bit more, replacing 6ull by 6u, if it does not cause overflow.

I have ignored slight bias mapping rand() to 0 .. 5
see http://www.azillionmonkeys.com/qed/random.html
see https://lemire.me/blog/2019/06/06/nearly...s-systems/

However, because of long tails, it is hard to get good estimate, even with good random numbers.
For 5 dices, correct expected value = 5.491 819 ...
05-12-2023, 06:14 PM
Post: #8
 Albert Chan Senior Member Posts: 2,679 Joined: Jul 2018
RE: little math/programming problems April 2023
(05-02-2023 08:03 PM)David Hayden Wrote:  int sides[6] = {0,0,0,0,2,3};

For a rough idea of expected values, we may try a simple "random" generator.
5 dices, with cycled pattern, different initial state:

Code:
000023000023000023000023000023000023 ... Pattern 1111122333  111112222233444   11111222223333344555    1111122222333334444455666     111112222233333444445555566777      11111222333334444455666

Note that 6 throws add to sum of only 5, sequence will always terminate.
Expected values = (3 + 4 + 5 + 6 + 7 + 6) / 6 = 31/6 ≈ 5.2

This "random" generator is super lumpy, with maximal 4 0's in a row.
We expected true expected values slightly bigger.
05-14-2023, 12:16 PM
Post: #9
 ttw Member Posts: 282 Joined: Jun 2014
RE: little math/programming problems April 2023
The following pRNG may be helpful if sampling the dice throws is a problem. It's a specially constructed linear congruential generator designed for die rolling.

Letting X stand for the internal state, one has

X = X * 93847949523997033 Mod ( 216172782113784439 )

The die roll is X ( Mod 6 ) or X ( Mod 6 ) + 1

The only constraint is that the first X must be non-zero 0 < Xo < 93847949523997033.

If this version has distributional problems, I got millions more available.
05-15-2023, 09:26 PM (This post was last modified: 05-15-2023 11:09 PM by Gilles.)
Post: #10
 Gilles Member Posts: 226 Joined: Oct 2014
RE: little math/programming problems April 2023
My (New)RPL version :
Code:
« → Loop    « 0      1 Loop START       5 'NbDe' LSTO         0 'TotalJet' LSTO        WHILE NbDe REPEAT         1 'TotalJet' STO+  0 'NewDe' LSTO          1 NbDe START           RAND 6 *            CASE             DUP 5 ≥ THEN 3 'NewDe' STO+ END             DUP 4 ≥ THEN 2 'NewDe' STO+ END           END           DROP          NEXT         NewDe 'NbDe' STO        END       TotalJet +      NEXT     Loop /    » »
Size : 296 Bytes

Execution on my laptop and HP50g, NewRPL
Code:
Loops    Time(s)    Result  Time           Laptop             HP50g ---------------------------------- 1000     0.015      5.438    4.1 10000    0.285      5.4475   39 100000   2.912      5.48428  389 (ratio vs PC ~ 133 ) 1E6      30.99      5.49546  ++ 10E6     311.79     5.49432  ++
05-16-2023, 07:47 PM
Post: #11
 pier4r Senior Member Posts: 2,248 Joined: Nov 2014
RE: little math/programming problems April 2023
Finally a calculator version as well, thanks Gilles! (and all the others as well of course)

Admittedly I could produce a calculator version too, but I didn't prioritize the task.

Wikis are great, Contribute :)
05-16-2023, 09:27 PM (This post was last modified: 05-16-2023 09:30 PM by Ajaja.)
Post: #12
 Ajaja Junior Member Posts: 33 Joined: Feb 2021
RE: little math/programming problems April 2023
My HP-42s/Free42 version:
Code:
00 { 61-Byte Prgm } 01▸LBL "DICE" 02 CLRG 03 STO 11 04 STO 09 05 2 06 STO 04 07 3 08 STO 05 09▸LBL 04 10 0 11 STO 06 12 5 13 STO 07 14▸LBL 03 15 0 16 STO 08 17▸LBL 01 18 RAN 19 6 20 × 21 RCL IND ST X 22 STO+ 08 23 DSE 07 24 GTO 01 25 1 26 STO+ 06 27 RCL 08 28 STO 07 29 X>0? 30 GTO 03 31 RCL 06 32 STO+ 10 33 DSE 09 34 GTO 04 35 RCL 10 36 RCL 11 37 ÷ 38 END
10000: 5,4555
100000: 5,49717
1000000: 5,492177
10000000: 5,4922995
05-17-2023, 09:48 PM
Post: #13
 pier4r Senior Member Posts: 2,248 Joined: Nov 2014
RE: little math/programming problems April 2023
Interesting the difference between the newRPL and the Free42.

I thought that the math library (at least for basic stuff like RAND) of the HP42 and RPL was similar but likely it is not.

Wikis are great, Contribute :)
05-18-2023, 07:54 AM (This post was last modified: 05-18-2023 07:58 AM by Ajaja.)
Post: #14
 Ajaja Junior Member Posts: 33 Joined: Feb 2021
RE: little math/programming problems April 2023
I think it just has a large deviation.
I've changed the program a bit to calculate SDEV:
Code:
00 { 57-Byte Prgm } 01▸LBL "DICE" 02 CLRG 03 STO 09 04 2 05 STO 04 06 3 07 STO 05 08▸LBL 04 09 CLX 10 STO 06 11 5 12 STO 07 13▸LBL 03 14 CLX 15 STO 08 16▸LBL 01 17 RAN 18 6 19 × 20 RCL IND ST X 21 STO+ 08 22 DSE 07 23 GTO 01 24 1 25 STO+ 06 26 RCL 08 27 STO 07 28 X>0? 29 GTO 03 30 RCL 06 31 PRX 32 Σ+ 33 DSE 09 34 GTO 04 35 SDEV 36 END
It's ~5
So, standard error of the mean here (SDEV/SQRT(N)):
10000: ~0.05
100000: ~0.016
1000000: ~0.005
10000000: ~0.0016
 « Next Oldest | Next Newest »

User(s) browsing this thread: 1 Guest(s)