Examples#

Monte Carlo Pi#

Estimate Pi by counting random points inside a unit circle:

li x := :randu(-1.0, 1.0)
li y := :randu(-1.0, 1.0)
li inside := | x^2 + y^2 <= 1.0 -> inside + 1
             | _                -> inside

obs_r estimate := 4.0 * inside / :i

run() {1_000_000, 10}
:print("Pi estimate:", :mean(estimate))

Geometric Brownian Motion#

Simulate a stock price with drift and volatility:

S0 := 100.0
mu := 0.05
sigma := 0.2
dt := 1.0 / 252.0

li price :=
    [0] = S0
    | _ -> price * :exp((mu - 0.5*sigma^2)*dt + sigma*:sqrt(dt)*:randn(0,1))

obs_r final_price := price

run() {252, 10000}
:print("Expected price:", :mean(final_price))
:print("Std deviation:", :stddev(final_price))

Branching Process Extinction#

Model population dynamics with random reproduction:

fn reproduce(n: int): int
    result := 0
    for i = 1, n do
        result = result + :randpoisson(1.0)
    end
    return result
end

li pop :=
    [0] = 10
    | pop == 0 -> 0
    | _        -> reproduce(pop)

obs_r extinct := pop == 0
stop_r(pop == 0, "extinct")

run() {100, 1000}
:print("Extinction probability:", :probability(extinct))