Sunday, June 5, 2011

Monte Carlo Simulation to Estimate pi

Imagine a person randomly throwing darts at a square board with a unit circle inscribed on it. The proportion of darts that land in the circle multiplied by 4 gives an estimate of the value of pi. The point (x, y) hit by the dart is in the circle if x^2 + y^2 <= 1.

What follows is my implementation of this idea in Teradata sql, with the value of pi estimated to two decimal places. The resource at the end of this post has a nice discourse on this topic, with Python implementations and a video on Monte Carlo simulation.
/* Create a table with a random point (x,y) on the board */
CREATE VOLATILE TABLE tblTemp AS
(SELECT RANDOM(-100,100)/100.00 AS x
,RANDOM(-100,100)/100.00 AS y
,CASE WHEN x**2+y**2 <= 1 THEN 1
ELSE 0
END AS inCirc
)WITH DATA PRIMARY INDEX (x)
ON COMMIT PRESERVE ROWS;

/* Generate a sufficiently large (but not too large!) number of random points (x,y) on the board */
CREATE VOLATILE TABLE tblRndPts AS
(WITH RECURSIVE tblRnd (lvl, x, y, inCirc) AS
(SELECT CAST(1 AS DECIMAL(18,0)) AS lvl, x, y, inCirc
FROM tblTemp

UNION ALL

SELECT lvl+1
,RANDOM(-100,100)/100.00 AS x1
,RANDOM(-100,100)/100.00 AS y1
,CASE WHEN x1**2+y1**2 <= 1 THEN 1
ELSE 0
END AS inCirc
FROM tblRnd AS t1
WHERE t1.lvl < 10000 /* total number of points (x, y) */
)
SELECT *
FROM tblRnd
) WITH DATA PRIMARY INDEX (lvl)
ON COMMIT PRESERVE ROWS;

/* Calculate pi and the diff from the actual value */
SELECT 4.00 * SUM(inCirc)/COUNT(*) AS pi, CAST(pi() - pi AS DECIMAL(4,2)) AS diff
FROM tblRndPts;

Click here to read more about this topic.