Skip to content

Instantly share code, notes, and snippets.

@rhz rhz/gillespie.hs
Created Apr 26, 2011

Embed
What would you like to do?
import qualified Data.Map as Map
import Random (randomRIO)
type Mol = String
type Mixture = Map.Map Mol Integer
type Rxn = (Mixture, (Mixture -> Mixture)) -- (Lhs, Action)
type State = (Double, Integer, Mixture) -- (Time, NumSteps, Mixture)
getActivities :: [Rxn] -> Mixture -> [Integer]
getActivities rxns mixture = []
selectRxn :: [Rxn] -> [Integer] -> Integer -> IO Rxn
selectRxn rxns activities totalActivity = return $ head rxns
-- do rnd <- getStdRandom $ randomR (0.0, totalActivity)
step :: [Rxn] -> State -> IO State
step rxns (t, ns, mixture) =
let activities = getActivities rxns mixture
totalActivity = sum activities
in if totalActivity == 0
then return (t, ns+1, mixture)
else do r <- selectRxn rxns activities totalActivity
rnd <- randomRIO (0.0, 1.0)
let dt = log (1.0/rnd) / totalActivity
return (t+dt, ns+1, (snd r) mixture)
simulate :: [Rxn] -> Mixture -> Integer -> IO [State]
simulate rxns mixture numSteps = take numSteps (iterate (step rxns) (0.0, 0, mixture))
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.