Monte Carlo

open System

type Point =
    { X : float
      Y : float }

let isInsideUnitCircle point = point.X ** 2.0 + point.Y ** 2.0 <= 1.0

let randomPoints =
    let random = Random()
    Seq.initInfinite (fun _ ->
        { X = random.NextDouble()
          Y = random.NextDouble() })

let pi batch =
    let a =
        (Seq.unfold (fun state ->
             let (total, pointsInside, _) = state

             let newPoinstsInside =
                 pointsInside + (Seq.take batch randomPoints
                                 |> Seq.filter isInsideUnitCircle
                                 |> Seq.length)

             let ratio = (float) pointsInside / (float) total
             let newState = (total + batch, newPoinstsInside, ratio * 4.0)
             Some(state, newState))) (0, 0, 0.0)
    seq {
        for v in a do
            let (_, _, pi) = v
            yield pi
    }

[<EntryPoint>]
let main argv =
    let nDefault = 100

    let n =
        let argvLength = Array.length argv
        match argvLength with
        | _ when argvLength >= 1 ->
            let first = Array.head argv
            match Int32.TryParse first with
            | true, n -> n
            | false, _ -> nDefault
        | _ -> nDefault
    for estimate in Seq.take 1000 (pi n) do
        printfn "%f" estimate
    0 // return an integer exit code