I wrote my own greatest common divisor (GCD) function, but just used the one in Emacs.
(require 'calc) (require 'calc-ext) ;; greatest common divisor (defun my-gcd (x y) (if (eq y 0) x (my-gcd y (% x y)))) ;; x and y are coprime if gcd(x, y) = 1 (defun coprimep (x y) (eq 1 (math-gcd x y))) (defun coprime-probability (n maxint) "Probability of coprime for N random numbers of MAXINT." (random t) (/ (apply '+ (mapcar (function (lambda (p) (if (coprimep (car p) (cdr p)) 1.0 0.0))) (mapcar* 'cons (mapcar 'random (make-vector n maxint)) (mapcar 'random (make-vector n maxint))))) n))
Here's an example of running the coprime generation using a million pairs and using as much precision Emacs can handle for integers, and then using the result to calculate pi.
(let* ((n #xfffff) ;; 1048575 (maxint #xfffffff) ;; 268435455 (p (coprime-probability n maxint)) (pi-approx (if (eq 0 p) 0.0 (sqrt (/ 6 p))))) (princ (format "Probability that %d pairs of integers less than %d are prime: %f\n" n maxint p)) (princ (format "pi (%s) is approximately " (math-format-number (math-pi)))) pi-approx)
Here's a result of running this once:
Probability that 1048575 pairs of integers less than 268435455 are prime: 0.607957
pi (3.14159265359) is approximately 3.141516663694778
It only takes a few seconds to run. Though, if you increase the number of pairs, Emacs quickly started choking my computer.
What was most intriguing to me about this factoid, was not the relationship to pi. I was surprised that over 60 percent of all points in two-dimensional space are coprime. I didn't expect that to be that high.
(require 'calc) (require 'calc-ext) (math-format-number (math-div 6 (math-pow (math-pi) 2))) >> "0.607927101854"
The probability of selecting a generally prime number is zero. "Apples to oranges", I know.