monte carlo 1d radiative transfer code written in haskell using random-mersenne package and infinite list of random numbers.

mc1drad_mersenne_list.hs

import System.Random.Mersenne

np  = 1000000
tau = 1.0 -- cloud optical thickness
sal = 0.5 -- surface albedo
omg = 1.0 -- single scattering albedo
zu  = tau -- atmosphere upper boundary
zb  = 0.0 -- atmosphere lower boundary

main :: IO ()
main = do
    seed0 <- getStdGen
    (sumr, sumt) <- randoms seed0 >>= world 
    print $ [sumr / np, sumt / np]

world :: [Double] -> IO (Double, Double)
world rn = loop 0 0.0 0.0 rn
    where
        loop :: Double -> Double -> Double -> [Double] -> IO (Double, Double)
        loop t sumr sumt rn'
            | t == np = return (sumr, sumt)
            | otherwise = do
                ((wr, wt), rn'') <- solver 1.0 zu (-1.0) 0.0 0.0 rn'
                loop (t+1) (sumr + wr) (sumt + wt) rn''

solver :: Double -> Double -> Double
          -> Double -> Double
          -> [Double] -> IO ((Double, Double), [Double])
solver w z dir sumr sumt rn = do
    (n:rn') <- return rn
    path <- return $ -log(n)
    nextphase (z + dir * path) rn'
    where
        nextphase zn rn''
            | zn <= zb = solver (w * sal) zb 1.0 sumr (sumt + w) rn''
            | zn > zu = return ((sumr + w, sumt), rn'')
            | otherwise  = do
                (n2:rn''') <- return rn''
                dnew <- return $ if n2 <= 0.5
                                   then 1.0 * dir
                                   else -1.0 * dir
                solver (w * omg) zn dnew sumr sumt rn'''

benchmark

for 1000000 photons

[0.5996505903320313,0.8006988193359375]
./mc1drad_mwc  2.37s user 0.27s system 97% cpu 2.693 total

[0.6001902515869141,0.7996194968261718]
./mc1drad_io  6.98s user 0.22s system 95% cpu 7.543 total

[0.599968560546875,0.80006287890625]
./mc1drad_mersenne  1.45s user 0.23s system 94% cpu 1.777 total

[0.5998440072021485,0.8003119855957032]
./mc1drad_mersenne_list  1.50s user 0.22s system 92% cpu 1.849 total

[0.6000919125976563,0.7998161748046875]
./mc1drad_mwc  2.36s user 0.27s system 96% cpu 2.721 total