Here's one of my favorite ways to do it. The probability that two randomly selected integers are coprime is 6/π^2.
Code:
#include <stdio.h>
#include <math.h>
#ifndef __unix__
#warning /dev/urandom may not be available on this system!
#endif
typedef unsigned long ulong;
ulong sample()
{
static FILE* stream = NULL;
if(!stream)
stream = fopen("/dev/urandom", "rb");
ulong result;
fread(&result, 1, sizeof(result), stream);
return result;
}
ulong gcd(ulong n1,ulong n2)
{
if(!n2)
return n1;
return gcd(n2,n1 % n2);
}
int main(void)
{
ulong coprime = 0, total = 0;
for(;;)
{
if(gcd(sample(),sample()) == 1)
++coprime;
double estimate = 1 / sqrt(coprime/(6 * (double)++total));
printf(" %.10g\n", estimate);
}
}