[ QUOTE ]
Great post. Could you share the programm?
[/ QUOTE ]
Sure. The program is written in C, and uses math.h, so remember to add the -lm flag when compiling. This version uses random() and gettimeofday() for a seed. If your C compiler doesnt have these functions just use rand(), and seed it with srand(time(NULL)). rand() is not a great random number generator, as the low order bits tend to cycle, though that shouldn't matter much here.
The box_muller function implements the Box-Muller transform, a method of generating random numbers along a normally distributed function. The main program is pretty simple, for each mean (bb/100) we want to examine, we set bankroll to +0 bets and just start simulating lots of hands. Each time our bankroll passes the previous high point, the difference between it and the lowest point reached since then is counted as a downswing. So, if we start at 0, lose 300 bbs, gain 100 bbs, lose 300 bbs, when we eventually get back to 0, that counts as 1 500 bb swing, not 2 300 bb swings.
<font class="small">Code:</font><hr /><pre>
//************************************************** ******************************
// random_walk.c
// Ronald Nussbaum
//************************************************** ******************************
#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <sys/time.h>
#define ITERATIONS 2000000000 // number of hands
#define STD 15 // std/100
// http://www.taygeta.com/random/boxmuller.html
double box_muller(double m, double s) // normal random variate generator
{ // mean m, standard deviation s
double x1, x2, y1, w;
static double y2;
static int use_last = 0;
if (use_last) // use value from previous call
{
y1 = y2;
use_last = 0;
}
else
{
do
{
x1 = 2.0 * ((double) random() / (double) RAND_MAX) - 1.0;
x2 = 2.0 * ((double) random() / (double) RAND_MAX) - 1.0;
w = x1 * x1 + x2 * x2;
}
while (w >= 1.0);
w = sqrt((-2.0 * log(w)) / w);
y1 = x1 * w;
y2 = x2 * w;
use_last = 1;
}
return (m + y1 * s);
}
int main(int argc, char **argv)
{
long int i, t, downswings[10];
long double mean, bankroll, high, low;
struct timeval tv;
gettimeofday(&tv, NULL);
srandom(tv.tv_usec);
printf("For std/100 = %d:\n", STD);
printf(" 100-199 200-299 300-399 400-499 500-599 600-699 700-799 800+\n");
for (mean = .5; mean <= 4; mean += .5)
{
high = low = bankroll = 0;
for (i = 0; i < 10; ++i)
downswings[i] = 0;
for (i = 0; i < ITERATIONS; i += 100)
{
bankroll += box_muller(mean, STD);
if (bankroll > high)
{
t = (int) (high - low) / 100;
if (t > 8)
t = 8;
++downswings[t];
high = low = bankroll;
}
else if (bankroll < low)
low = bankroll;
}
printf("%.1f ", (double) mean);
for (i = 1; i <= 8; ++i)
if (downswings[i] > 0)
{
if ((ITERATIONS / downswings[i]) < 9999999)
printf("%7ld ", ITERATIONS / downswings[i]);
else printf(" %6.0E ", (double) (ITERATIONS / downswings[i]));
}
else printf(" ");
printf("\n");
}
return EXIT_SUCCESS;
}
</pre><hr />
On a 32 bit machine, ITERATIONS cannot exceed the maximum size of a long int, which is just over 2 billion. To get more accurate approximations, you either need to use a large number package, or a 64 bit machine, which has a limit in the quintillions. Here are some results from a 64 bit machine simulating 100 billion hands per trial:
<font class="small">Code:</font><hr /><pre>
For std/100 = 15:
100-199 200-299 300-399 400-499 500-599 600-699 700-799 800+
0.5 54087 152788 323449 603343 1046725 1733763 2810962 1651936
1.0 30426 104171 289054 745901 1868530 4599816 1E+07 2E+07
1.5 24429 112099 452988 1784630 7061149 3E+07 1E+08 3E+08
2.0 23419 153301 962343 6127826 4E+07 2E+08 1E+09 5E+09
2.5 25085 245873 2487809 3E+07 3E+08 3E+09 2E+10 1E+11
3.0 29107 438392 7304601 1E+08 3E+09 2E+10 1E+11
3.5 35756 838475 2E+07 7E+08 2E+10 1E+11
4.0 45889 1731871 8E+07 4E+09 1E+11
For std/100 = 20:
100-199 200-299 300-399 400-499 500-599 600-699 700-799 800+
0.5 57068 149643 290361 490186 760121 1115436 1598210 544988
1.0 29780 83583 177262 338975 606152 1053363 1798011 1200984
1.5 21257 65961 162519 370933 821631 1774119 3829510 4298302
2.0 17481 62315 184495 526803 1465502 4100209 1E+07 2E+07
2.5 15666 65967 243236 880871 3221338 1E+07 4E+07 1E+08
3.0 14944 76072 354111 1681491 8045699 4E+07 2E+08 7E+08
3.5 14927 93370 559665 3482864 2E+07 1E+08 9E+08 3E+09
4.0 15461 120526 938086 7792410 6E+07 6E+08 5E+09 3E+10
</pre><hr />
I hope this is useful. Please let me know if I made any serious mathematical goofs in my code (hopefully not).
Senni