## Posts Tagged ‘pi’

### Programming Praxis – Buffon’s Needle

March 15, 2013

In today’s Programming Praxis exercise, our goal is to approximate pi by using Georges-Louis Leclerc’s method of dropping needles on a board. Let’s get started, shall we?

```import Control.Applicative
import System.Random```

We’re going to be needing two sets of random numbers; one for the position and one for the angle of the needles. In order to save some code and to stop Haskell complaining about ambiguous types we make a function to generate an infinite amount of numbers in the [0,1) range.

```rnds :: IO [Double]
rnds = fmap randoms newStdGen```

We approximate pi by dividing the total amount of needles dropped by the number of needles that hit a line.

```buffon :: Int -> IO Double
buffon n = (fromIntegral n /) . sum . take n <\$>
(zipWith (\y t -> if y < sin (t*pi/2) / 2 then 1 else 0) <\$> rnds <*> rnds)```

Running the simulation reveals that this isn’t a very practical way of approximating pi: after one million needles it generally only has two correct digits.

```main :: IO ()
main = print =<< buffon 1000000```

### Programming Praxis – Calculating Pi

October 9, 2009

In today’s Programming Praxis problem we have to implement two ways of calculating pi. Let’s get started, shall we?

Some imports:

```import Control.Monad
import System.Random
```

The first algorithm is the Monte Carlo method. Sadly, the need for monads and the int-float conversion make this code less concise than it could be.

```monteCarloPi :: Int -> IO ()
monteCarloPi n = do
hits <- fmap sum \$ liftM2 (zipWith checkHit) rs rs
print \$ fromIntegral hits / fromIntegral n
where rs = replicateM n \$ randomRIO (0,1 :: Double)
checkHit x y = if x*x + y*y < 1 then 4 else 0
```

Next up is Archimedes’ algorithm.

```boundPi :: Int -> (Double, Double)
boundPi n = iterate f (3/2 * sqrt 3, 3 * sqrt 3) !! (n - 1)
where f (b, a) = let x = 2 / (1 / a + 1 / b)
in (sqrt \$ b * x, x)
```

A quick test to see if everything is working correctly:

```main :: IO ()
main = do monteCarloPi 10000
print \$ boundPi 6
```

Everything seems to be in order. Of course, we could just say

`main = print pi`

🙂