dc

Write a program in dc!

"monte carlo estimator for pi":
[? 2^ ? 2^ + 1>i]su
[lx 1+ sx]si
[lu x lm 1+ d sm ln>z]sz

5k
?sn
lzx
lx ln / 4* p
q

It needs a stream of floats:
/* frand.c */
#include <u.h>
#include <libc.h>

void
main(void)
{
srand(time(0));

for(;;)
print("%f\n", frand());
}

Driver:
#!/bin/rc
# runpi <number of samples>

{ echo $1; frand } | dc pi.dc

Example:
% runpi 10000
3.14840

>>3

Now I got it working. Thanks for the idea!

It took a while to figure out because setting the scale factor k fucked up the integer operations for RANDU.

Also, I'm quite pleased with myself for coming up with the pseudocomment string at the top.

[
Registers:
u - routine : execute i if sum of squares less than 1
i - routine : increment register x
z - routine : iterator - execute u while n > m++
r - routine : RANDU PRNG
m - variable: number of samples
x - variable: number of samples inside circle
s - variable: seed for r
k - variable: scale for division
n - variable: number of iterations (user input)
]c
[lrx 2^ lrx 2^ + 1>i]su
[lx 1+ sx]si
[lu x lm 1+ d sm ln>z]sz
[0k ls 65539 * 2 31^ % d ss lkk 2 31 ^ /]sr
? sn
5dksk
1 ss
lzx
lx lm / 4*
p