Permalink
Browse files

fragmentation is now possible!

  • Loading branch information...
1 parent c8d51f7 commit cf51d5380d4cbf8dfec3f65d6eb535e108f26424 @rhz committed Aug 7, 2012
Showing with 444 additions and 120 deletions.
  1. +48 −0 Frag.hs
  2. +104 −0 Fragmentation.hs
  3. +12 −16 KappaParser.hs
  4. +66 −75 Matching.hs
  5. +58 −27 Mixture.hs
  6. +7 −0 README.md
  7. +141 −0 Rule.hs
  8. +8 −2 Utils.hs
View
@@ -0,0 +1,48 @@
+module Main where
+
+import qualified KappaParser as KP
+import qualified Env as E
+import qualified Mixture as M
+import qualified Rule as R
+import qualified Fragmentation as F
+import qualified PlainKappa as K
+import Utils
+
+import qualified Data.Map as Map
+import qualified Data.Set as Set
+import System.Environment (getArgs)
+import Data.List (intercalate, nubBy)
+
+main :: IO ()
+main = do inputFilename : nStr : _ <- getArgs
+ m <- KP.parseFromFile inputFilename
+
+ let env = E.createEnv m
+ rules = [ R.evalRule env rule | (_, rule) <- KP.rules m ]
+ obss = [ (M.evalKExpr env True obs, name) | KP.KExprWithName name obs <- KP.obss m ]
+
+ odes = take (read nStr) $ F.fragment rules (map fst obss)
+ fragments = map fst odes ++ (snd =<< snd =<< odes)
+
+ names = reverse $ addFrags 1 (reverse obss) fragments -- reverse forth and back
+ addFrags _ obss [] = obss
+ addFrags n obss (obs:todo)
+ | (obs', name):_ <- filter (F.ccIso obs . fst) obss = addFrags n ((obs, name):obss) todo -- reuse name
+ | otherwise = addFrags (n+1) ((obs, "F" ++ show n):obss) todo -- assign a name
+
+ mapM_ (printFrag env) (nubBy ((==) `on` snd) names)
+ putStrLn ""
+ mapM_ (printODE (Map.fromList names)) odes
+ where
+ printFrag :: E.Env -> (F.Obs, String) -> IO ()
+ printFrag env (obs, name) = putStrLn $ name ++ " := " ++ M.toKappa env obs
+
+ printODE :: Map.Map F.Obs String -> (F.Obs, F.ODE) -> IO ()
+ printODE names (obs, ode) = putStrLn $ names Map.! obs ++ " = " ++ odeRhs
+ where
+ odeRhs | null ode = "0"
+ | otherwise = intercalate " + " (map showODET ode)
+
+ showODET :: F.ODET -> String
+ showODET (rate, obss) = K.showAExpr rate ++ " " ++ intercalate " " (map (names Map.!) obss)
+
View
@@ -0,0 +1,104 @@
+{-# LANGUAGE TupleSections #-}
+
+module Fragmentation where
+
+import qualified Data.Vector as Vec
+import qualified Data.Map as Map
+import qualified Data.Set as Set
+
+import qualified KappaParser as KP
+import qualified Mixture as M
+import qualified Rule as R
+import Matching
+import Utils
+
+type Obs = M.Mixture
+type ODET = (R.Rate, [Obs]) -- ODE term
+type ODE = [ODET]
+type ODES = [(Obs, ODE)] -- ODE system, that is [( Obs, [(R.Rate, [Obs])] )]
+
+-- TODO check if the first two tests are any good for performance
+-- the first probably is, but the second one I'm not that sure
+ccIso :: M.Mixture -> M.Mixture -> Bool
+ccIso m1 m2 = agentCount m1 == agentCount m2 -- same amount of agents
+ && typeCount m1 == typeCount m2 -- same 'typing map'
+ && (agentCount m1 == 0 || any (match Set.empty Set.empty) anchors)
+ where
+ agentCount = Vec.length . M.agents
+ typeCount = frequencies . map M.agentName . agentList
+ agentList = Vec.toList . M.agents
+
+ anchors = return . (0, ) <$> M.agentIds m2
+
+ -- This function assumes the given m1 and m2 are connected components, it doesn't check if that's true.
+ -- That's why we return True when we don't have anything else in the todo list,
+ -- because that means we have visited all agents in both m1 and m2
+ match :: Set.Set M.AgentId -> Set.Set M.AgentId -> [(M.AgentId, M.AgentId)] -> Bool
+ match visited1 visited2 [] = True
+ match visited1 visited2 ((id1, id2) : todo)
+ | Set.member id1 visited1 && Set.member id2 visited2 = match visited1 visited2 todo -- TODO should I check that all nbs are visited?
+ | Set.member id1 visited1 = error "Fragmentation.ccIso.match: id1 visited but not id2"
+ | Set.member id2 visited2 = error "Fragmentation.ccIso.match: id2 visited but not id1"
+ | a1 == a2 = match (Set.insert id1 visited1) (Set.insert id2 visited2) (todo ++ nbs)
+ | otherwise = False
+ where
+ nbs = do (sId, (M.Site{ M.bindingState = M.Bound}, M.Site{ M.bindingState = M.Bound })) <- indexedList $ Vec.zip (M.interface a1) (M.interface a2)
+ let (nb1, _) = M.follow m1 (id1, sId) ? "Fragmentation.ccIso.match: disconnected graph (1)"
+ (nb2, _) = M.follow m2 (id2, sId) ? "Fragmentation.ccIso.match: disconnected graph (2)"
+ return (nb1, nb2)
+ -- nbsVisited = all (visited1 `Set.member`) (map fst nbs) && all (visited2 `Set.member`) (map snd nbs)
+
+ a1 = M.agents m1 Vec.!? id1 ? "Fragmentation.ccIso.match: agent id not found"
+ a2 = M.agents m2 Vec.!? id2 ? "Fragmentation.ccIso.match: agent id not found"
+
+
+fragment :: [R.Rule] -> [Obs] -> ODES
+fragment rules obss = fragment' obss []
+ where
+ fragment' :: [Obs] -> [Obs] -> ODES
+ fragment' [] _ = []
+ fragment' (obs:todo) visited
+ | any (ccIso obs) visited = fragment' todo visited
+ | otherwise = (obs, ode) : fragment' (todo ++ newFragments) (obs : visited)
+ where
+ ode = getODE obs =<< rules
+ newFragments = snd =<< ode
+
+ getODE :: Obs -> R.Rule -> ODE
+ getODE obs rule = lhsTerms ++ rhsTerms
+ where
+ (msitesLhs, msitesRhs) = R.modifiedSites rule
+ minglueingsLhs = minimalGlueings obs (R.lhs rule)
+ minglueingsRhs = minimalGlueings obs (R.rhs rule)
+
+ -- you have to filter out the minimal glueings that don't intersect in the pullback with the set of modified sites
+ -- when you glue on the left, the new observable is just the minimal glueing
+ -- when you glue on the right, the new observable is the inverse of the rule applied to the minimal glueing
+
+ lhsRelevantMG = filter (isRelevant msitesLhs (R.lhs rule)) minglueingsLhs
+ rhsRelevantMG = filter (isRelevant msitesRhs (R.rhs rule)) minglueingsRhs
+
+ isRelevant :: [M.Endpoint] -> M.Mixture -> (M.Mixture, Injection, Injection, M.Mixture) -> Bool
+ isRelevant msites m2 (m0, _, m2Inj, _) = any inPullback msites
+ where
+ inPullback (aId, sId) = agentInPullback && s0 `siteMatch` s2 -- TODO why should I check if s2 matches s0?
+ where id0 = m2Inj Map.! aId
+ s0 = M.interface (M.agents m0 Vec.! id0) Vec.! sId
+ s2 = M.interface (M.agents m2 Vec.! aId) Vec.! sId
+ agentInPullback = id0 < Vec.length (M.agents m0)
+
+ lhsTerms = (neg $ R.rate rule, ) <$> M.split <$> codomain <$> lhsRelevantMG
+ rhsTerms = ( R.rate rule, ) <$> M.split <$> invert rule <$> rhsRelevantMG
+
+ codomain (_, _, _, m3) = m3
+
+ invert :: R.Rule -> (M.Mixture, Injection, Injection, M.Mixture) -> M.Mixture
+ invert rule (_, _, rhsInj, m3) = R.apply invertedScript rhsInj m3
+ where invertedScript = R.actionScript (R.rhs rule) (R.lhs rule)
+
+
+neg :: KP.AExpr -> KP.AExpr
+neg (KP.Integer n) = KP.Integer (negate n)
+neg (KP.Float n) = KP.Float (negate n)
+neg x = KP.Duo KP.Mult (KP.Integer (-1)) x
+
View
@@ -108,7 +108,8 @@ createChain first@(Agent fname fintf) second@(Agent sname sintf) last@(Agent lna
error $ "KappaParser.createChain: all agents in a chain must have the same sites in their interface"
| firstLink /= firstLink' =
error $ "KappaParser.createChain: first and second agents in chain must be bound by sites '" ++ rightSite ++ "' and '" ++ leftSite ++ "', respectively"
- | otherwise = first : take n agentsInChain ++ [last]
+ | otherwise =
+ first : take n agentsInChain ++ [last]
where
agentsInChain = iterate nextAgentInChain second
n = (lastLink - firstLink) `quot` step
@@ -135,17 +136,15 @@ createChain first@(Agent fname fintf) second@(Agent sname sintf) last@(Agent lna
hasSameSites :: Interface -> Interface -> Bool
hasSameSites i1 i2 = map siteName i1 == map siteName i2
-kexpr :: Parser KExpr
-kexpr = reverse . unpackChains [] [] <$> commaSep1 (liftM Right agent <|> liftM Left ellipsis) <?> "kappa expression"
- where
- ellipsis = reserved "..."
+unpackChains :: KExpr -> [Either () Agent] -> KExpr
+unpackChains acc [] = reverse acc
+unpackChains acc ((Right a1):(Right a2):(Left ()):(Right a3):xs) = unpackChains (reverse (createChain a1 a2 a3) ++ acc) xs
+unpackChains acc ((Right a):xs) = unpackChains (a:acc) xs
+unpackChains acc xs = error $ "malformed chain expression"
- unpackChains :: KExpr -> KExpr -> [Either () Agent] -> KExpr
- unpackChains acc [b2,b1] [] = b1:b2:acc
- unpackChains acc [b2,b1] ((Right a):xs) = unpackChains (b1:acc) [a,b2] xs
- unpackChains acc buf ((Right a):xs) = unpackChains acc (a:buf) xs
- unpackChains acc [b2,b1] ((Left _):(Right a):xs) = unpackChains (reverse (createChain b1 b2 a) ++ acc) [] xs
- unpackChains _ _ _ = error "malformed chain expression"
+kexpr :: Parser KExpr
+kexpr = unpackChains [] <$> commaSep1 (liftM Right agent <|> liftM Left ellipsis) <?> "kappa expression"
+ where ellipsis = reserved "..."
rule :: Parser Rule
rule = do lhs <- kexpr
@@ -310,11 +309,8 @@ obsP = do reserved "obs:"
varP :: Parser Var
varP = do name <- identifier
reservedOp "="
- ke <- kexpr
- if null ke
- then do ae <- aexpr
- return (name, Right ae)
- else return (name, Left ke)
+ expr <- liftM Left kexpr <|> liftM Right aexpr
+ return (name, expr)
energyShape :: Parser Shape
energyShape = do expr <- kexpr
Oops, something went wrong.

0 comments on commit cf51d53

Please sign in to comment.