msakai / toysolver / 496

Committed 10 Nov 2024 11:05AM UTC coverage: 69.994% (-1.1%) from 71.113%

Build # 496

Build Type

push

github

Committed by

web-flow

Commit Message

Merge pull request #117 from msakai/update-coveralls-and-haddock

GitHub Actions: Update coveralls and haddock configuration

Run Details

9872 of 14104 relevant lines covered (69.99%)

0.7 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

57.93

/src/ToySolver/Graph/ShortestPath.hs

{-# OPTIONS_GHC -Wall #-}
{-# OPTIONS_HADDOCK show-extensions #-}
{-# LANGUAGE ExistentialQuantification #-}
{-# LANGUAGE ScopedTypeVariables #-}
--------------------------------------------------------------------------
-- |
-- Module      :  ToySolver.Graph.ShortestPath
-- Copyright   :  (c) Masahiro Sakai 2016-2017
-- License     :  BSD-style
--
-- Maintainer  :  masahiro.sakai@gmail.com
-- Stability   :  provisional
-- Portability :  non-portable
--
-- This module provides functions for shotest path computation.
--
-- Reference:
--
-- * Friedrich Eisenbrand. “Linear and Discrete Optimization”.
--   <https://www.coursera.org/course/linearopt>
--
--------------------------------------------------------------------------
module ToySolver.Graph.ShortestPath
  (
  -- * Graph data types
    Graph
  , Edge
  , OutEdge
  , InEdge

  -- * Fold data type
  , Fold (..)
  , monoid'
  , monoid
  , unit
  , pair
  , path
  , firstOutEdge
  , lastInEdge
  , cost

  -- * Path data types
  , Path (..)
  , pathFrom
  , pathTo
  , pathCost
  , pathEmpty
  , pathAppend
  , pathEdges
  , pathEdgesBackward
  , pathEdgesSeq
  , pathVertexes
  , pathVertexesBackward
  , pathVertexesSeq
  , pathFold
  , pathMin

  -- * Shortest-path algorithms
  , bellmanFord
  , dijkstra
  , floydWarshall

  -- * Utility functions
  , bellmanFordDetectNegativeCycle
  ) where

import Control.Monad
import Control.Monad.ST
import Control.Monad.Trans
import Control.Monad.Trans.Except
import Data.Hashable
import qualified Data.HashTable.Class as H
import qualified Data.HashTable.ST.Cuckoo as C
import Data.IntMap.Strict (IntMap)
import qualified Data.IntMap.Strict as IntMap
import qualified Data.IntSet as IntSet
import qualified Data.Heap as Heap -- http://hackage.haskell.org/package/heaps
import Data.List (foldl')
import Data.Maybe (fromJust)
import Data.Monoid
import Data.Ord
import Data.Sequence (Seq)
import qualified Data.Sequence as Seq
import Data.STRef

-- ------------------------------------------------------------------------

-- | Graph represented as a map from vertexes to their outgoing edges
type Graph cost label = IntMap [OutEdge cost label]

-- | Vertex data type
type Vertex = Int

-- | Edge data type
type Edge cost label = (Vertex, Vertex, cost, label)

-- | Outgoing edge data type (source vertex is implicit)
type OutEdge cost label = (Vertex, cost, label)

-- | Incoming edge data type (target vertex is implicit)
type InEdge cost label = (Vertex, cost, label)

-- | path data type.
data Path cost label
  = Empty Vertex
    -- ^ empty path
  | Singleton (Edge cost label)
    -- ^ path with single edge
  | Append (Path cost label) (Path cost label) !cost
    -- ^ concatenation of two paths
  deriving (Eq, Show)

-- | Source vertex
pathFrom :: Path cost label -> Vertex
pathFrom (Empty v) = v
pathFrom (Singleton (from,_,_,_)) = from
pathFrom (Append p1 _ _) = pathFrom p1

-- | Target vertex
pathTo :: Path cost label -> Vertex
pathTo (Empty v) = v
pathTo (Singleton (_,to,_,_)) = to
pathTo (Append _ p2 _) = pathTo p2

-- | Cost of a path
pathCost :: Num cost => Path cost label -> cost
pathCost (Empty _) = 0
pathCost (Singleton (_,_,c,_)) = c
pathCost (Append _ _ c) = c

-- | Empty path
pathEmpty :: Vertex -> Path cost label
pathEmpty = Empty

-- | Concatenation of two paths
pathAppend :: (Num cost) => Path cost label -> Path cost label -> Path cost label
pathAppend p1 p2
  | pathTo p1 /= pathFrom p2 = error "ToySolver.Graph.ShortestPath.pathAppend: pathTo/pathFrom mismatch"
  | otherwise =
      case (p1, p2) of
        (Empty _, _) -> p2
        (_, Empty _) -> p1
        _ -> Append p1 p2 (pathCost p1 + pathCost p2)

-- | Edges of a path
pathEdges :: Path cost label -> [Edge cost label]
pathEdges p = f p []
  where
    f (Empty _) xs = xs
    f (Singleton e) xs = e : xs
    f (Append p1 p2 _) xs = f p1 (f p2 xs)

-- | Edges of a path, but in the reverse order
pathEdgesBackward :: Path cost label -> [Edge cost label]
pathEdgesBackward p = f p []
  where
    f (Empty _) xs = xs
    f (Singleton e) xs = e : xs
    f (Append p1 p2 _) xs = f p2 (f p1 xs)

-- | Edges of a path, but as `Seq`
pathEdgesSeq :: Path cost label -> Seq (Edge cost label)
pathEdgesSeq (Empty _) = Seq.empty
pathEdgesSeq (Singleton e) = Seq.singleton e
pathEdgesSeq (Append p1 p2 _) = pathEdgesSeq p1 <> pathEdgesSeq p2

-- | Vertexes of a path
pathVertexes :: Path cost label -> [Vertex]
pathVertexes p = pathFrom p : f p []
  where
    f (Empty _) xs = xs
    f (Singleton (_,v2,_,_)) xs = v2 : xs
    f (Append p1 p2 _) xs = f p1 (f p2 xs)

-- | Vertexes of a path, but in the reverse order
pathVertexesBackward :: Path cost label -> [Vertex]
pathVertexesBackward p = pathTo p : f p []
  where
    f (Empty _) xs = xs
    f (Singleton (v1,_,_,_)) xs = v1 : xs
    f (Append p1 p2 _) xs = f p2 (f p1 xs)

-- | Vertex of a path, but as `Seq`
pathVertexesSeq :: Path cost label -> Seq Vertex
pathVertexesSeq p = f True p
  where
    f True  (Empty v) = Seq.singleton v
    f False (Empty _) = mempty
    f True  (Singleton (v1,v2,_,_)) = Seq.fromList [v1, v2]
    f False (Singleton (v1,_,_,_))  = Seq.singleton v1
    f b (Append p1 p2 _) = f False p1 <> f b p2

pathMin :: (Num cost, Ord cost) => Path cost label -> Path cost label -> Path cost label
pathMin p1 p2
  | pathCost p1 <= pathCost p2 = p1
  | otherwise = p2

-- | Fold a path using a given 'Fold` operation.
pathFold :: Fold cost label a -> Path cost label -> a
pathFold (Fold fV fE fC fD) p = fD (h p)
  where
    h (Empty v) = fV v
    h (Singleton e) = fE e
    h (Append p1 p2 _) = fC (h p1) (h p2)

-- ------------------------------------------------------------------------

-- | Strict pair type
data Pair a b = Pair !a !b

-- ------------------------------------------------------------------------

-- | Operations for folding edge information along with a path into a @r@ value.
--
-- @Fold cost label r@ consists of three operations
--
-- * @Vertex -> a@ corresponds to an empty path,
--
-- * @Edge cost label -> a@ corresponds to a single edge,
--
-- * @a -> a -> a@ corresponds to path concatenation,
--
-- and @a -> r@ to finish the computation.
data Fold cost label r
  = forall a. Fold (Vertex -> a) (Edge cost label -> a) (a -> a -> a) (a -> r)

instance Functor (Fold cost label) where
  {-# INLINE fmap #-}
  fmap f (Fold fV fE fC fD) = Fold fV fE fC (f . fD)

instance Applicative (Fold cost label) where
  {-# INLINE pure #-}
  pure a = Fold (\_ -> ()) (\_ -> ()) (\_ _ -> ()) (const a)

  {-# INLINE (<*>) #-}
  Fold fV1 fE1 fC1 fD1 <*> Fold fV2 fE2 fC2 fD2 =
    Fold (\v -> Pair (fV1 v) (fV2 v))
         (\e -> Pair (fE1 e) (fE2 e))
         (\(Pair a1 b1) (Pair a2 b2) -> Pair (fC1 a1 a2) (fC2 b1 b2))
         (\(Pair a b) -> fD1 a (fD2 b))

-- | Project `Edge` into a monoid value and fold using monoidal operation.
monoid' :: Monoid m => (Edge cost label -> m) -> Fold cost label m
monoid' f = Fold (\_ -> mempty) f mappend id

-- | Similar to 'monoid'' but a /label/ is directly used as a monoid value.
monoid :: Monoid m => Fold cost m m
monoid = monoid' (\(_,_,_,m) -> m)

-- | Ignore contents and return a unit value.
unit :: Fold cost label ()
unit = monoid' (\_ -> ())

-- | Pairs two `Fold` into one that produce a tuple.
pair :: Fold cost label a -> Fold cost label b -> Fold cost label (a,b)
pair (Fold fV1 fE1 fC1 fD1) (Fold fV2 fE2 fC2 fD2) =
  Fold (\v -> Pair (fV1 v) (fV2 v))
       (\e -> Pair (fE1 e) (fE2 e))
       (\(Pair a1 b1) (Pair a2 b2) -> Pair (fC1 a1 a2) (fC2 b1 b2))
       (\(Pair a b) -> (fD1 a, fD2 b))

-- | Construct a `Path` value.
path :: (Num cost) => Fold cost label (Path cost label)
path = Fold pathEmpty Singleton pathAppend id

-- | Compute cost of a path.
cost :: Num cost => Fold cost label cost
cost = Fold (\_ -> 0) (\(_,_,c,_) -> c) (+) id

-- | Get the first `OutEdge` of a path.
firstOutEdge :: Fold cost label (First (OutEdge cost label))
firstOutEdge = monoid' (\(_,v,c,l) -> First (Just (v,c,l)))

-- | Get the last `InEdge` of a path.
-- This is useful for constructing a /parent/ map of a spanning tree.
lastInEdge :: Fold cost label (Last (InEdge cost label))
lastInEdge = monoid' (\(v,_,c,l) -> Last (Just (v,c,l)))

-- ------------------------------------------------------------------------

-- | Bellman-Ford algorithm for finding shortest paths from source vertexes
-- to all of the other vertices in a weighted graph with negative weight
-- edges allowed.
--
-- It compute shortest-paths from given source vertexes, and folds edge
-- information along shortest paths using a given 'Fold' operation.
bellmanFord
  :: Real cost
  => Fold cost label a
     -- ^ Operation used to fold shotest paths
  -> Graph cost label
  -> [Vertex]
     -- ^ List of source vertexes @vs@
  -> IntMap (cost, a)
bellmanFord (Fold fV fE fC fD) g ss = runST $ do
  let n = IntMap.size g
  d <- C.newSized n
  forM_ ss $ \s -> H.insert d s (Pair 0 (fV s))

  updatedRef <- newSTRef (IntSet.fromList ss)
  _ <- runExceptT $ replicateM_ n $ do
    us <- lift $ readSTRef updatedRef
    when (IntSet.null us) $ throwE ()
    lift $ do
      writeSTRef updatedRef IntSet.empty
      forM_ (IntSet.toList us) $ \u -> do
        -- modifySTRef' updatedRef (IntSet.delete u) -- possible optimization
        Pair du a <- liftM fromJust $ H.lookup d u
        forM_ (IntMap.findWithDefault [] u g) $ \(v, c, l) -> do
          m <- H.lookup d v
          case m of
            Just (Pair c0 _) | c0 <= du + c -> return ()
            _ -> do
              H.insert d v (Pair (du + c) (a `fC` fE (u,v,c,l)))
              modifySTRef' updatedRef (IntSet.insert v)

  liftM (fmap (\(Pair c x) -> (c, fD x))) $ freezeHashTable d

freezeHashTable :: H.HashTable h => h s Int v -> ST s (IntMap v)
freezeHashTable h = H.foldM (\m (k,v) -> return $! IntMap.insert k v m) IntMap.empty h

-- | Utility function for detecting a negative cycle.
bellmanFordDetectNegativeCycle
  :: forall cost label a. Real cost
  => Fold cost label a
     -- ^ Operation used to fold a cycle
  -> Graph cost label
  -> IntMap (cost, Last (InEdge cost label))
     -- ^ Result of @'bellmanFord' 'lastInEdge' vs@
  -> Maybe a
bellmanFordDetectNegativeCycle (Fold fV fE fC fD) g d = either (Just . fD) (const Nothing) $ do
  forM_ (IntMap.toList d) $ \(u,(du,_)) ->
    forM_ (IntMap.findWithDefault [] u g) $ \(v, c, l) -> do
      let (dv,_) = d IntMap.! v
      when (du + c < dv) $ do
        -- a negative cycle is detected
        let d' = IntMap.insert v (du + c, Last (Just (u, c, l))) d
            parent u = do
              case IntMap.lookup u d' of
                Just (_, Last (Just (v,_,_))) -> v
                _ -> undefined
            u0 = go (parent (parent v)) (parent v)
              where
                go hare tortoise
                  | hare == tortoise = hare
                  | otherwise = go (parent (parent hare)) (parent tortoise)
        let go u result = do
              let Just (_, Last (Just (v,c,l))) = IntMap.lookup u d'
              if v == u0 then
                fC (fE (v,u,c,l)) result
              else
                go v (fC (fE (v,u,c,l)) result)
        Left $ go u0 (fV u0)

-- ------------------------------------------------------------------------

-- | Dijkstra's algorithm for finding shortest paths from source vertexes
-- to all of the other vertices in a weighted graph with non-negative edge
-- weight.
--
-- It compute shortest-paths from given source vertexes, and folds edge
-- information along shortest paths using a given 'Fold' operation.
dijkstra
  :: forall cost label a. Real cost
  => Fold cost label a
     -- ^ Operation used to fold shotest paths
  -> Graph cost label
  -> [Vertex]
     -- ^ List of source vertexes
  -> IntMap (cost, a)
dijkstra (Fold fV fE fC (fD :: x -> a)) g ss =
  fmap (\(Pair c x) -> (c, fD x)) $
    loop (Heap.fromList [Heap.Entry 0 (Pair s (fV s)) | s <- ss]) IntMap.empty
  where
    loop
      :: Heap.Heap (Heap.Entry cost (Pair Vertex x))
      -> IntMap (Pair cost x)
      -> IntMap (Pair cost x)
    loop q visited =
      case Heap.viewMin q of
        Nothing -> visited
        Just (Heap.Entry c (Pair v a), q')
          | v `IntMap.member` visited -> loop q' visited
          | otherwise ->
              let q2 = Heap.fromList
                       [ Heap.Entry (c+c') (Pair ch (a `fC` fE (v,ch,c',l)))
                       | (ch,c',l) <- IntMap.findWithDefault [] v g
                       , not (ch `IntMap.member` visited)
                       ]
              in loop (Heap.union q' q2) (IntMap.insert v (Pair c a) visited)

-- ------------------------------------------------------------------------

-- | Floyd-Warshall algorithm for finding shortest paths in a weighted graph
-- with positive or negative edge weights (but with no negative cycles).
--
-- It compute shortest-paths between each pair of vertexes, and folds edge
-- information along shortest paths using a given 'Fold' operation.
floydWarshall
  :: forall cost label a. Real cost
  => Fold cost label a
     -- ^ Operation used to fold shotest paths
  -> Graph cost label
  -> IntMap (IntMap (cost, a))
floydWarshall (Fold fV fE fC (fD :: x -> a)) g =
  fmap (fmap (\(Pair c x) -> (c, fD x))) $
    IntMap.unionWith (IntMap.unionWith minP) (foldl' f tbl0 vs) paths0
  where
    vs = IntMap.keys g

    paths0 :: IntMap (IntMap (Pair cost x))
    paths0 = IntMap.mapWithKey (\v _ -> IntMap.singleton v (Pair 0 (fV v))) g

    tbl0 :: IntMap (IntMap (Pair cost x))
    tbl0 = IntMap.mapWithKey (\v es -> IntMap.fromListWith minP [(u, (Pair c (fE (v,u,c,l)))) | (u,c,l) <- es]) g

    minP :: Pair cost x -> Pair cost x -> Pair cost x
    minP = minBy (comparing (\(Pair c _) -> c))

    f :: IntMap (IntMap (Pair cost x))
      -> Vertex
      -> IntMap (IntMap (Pair cost x))
    f tbl vk =
      case IntMap.lookup vk tbl of
        Nothing -> tbl
        Just hk -> IntMap.map h tbl
          where
            h :: IntMap (Pair cost x) -> IntMap (Pair cost x)
            h m =
              case IntMap.lookup vk m of
                Nothing -> m
                Just (Pair c1 x1) -> IntMap.unionWith minP m (IntMap.map (\(Pair c2 x2) -> (Pair (c1+c2) (fC x1 x2))) hk)

minBy :: (a -> a -> Ordering) -> a -> a -> a
minBy f a b =
  case f a b of
    LT -> a
    EQ -> a
    GT -> b

-- ------------------------------------------------------------------------

1	{-# OPTIONS_GHC -Wall #-}
2	{-# OPTIONS_HADDOCK show-extensions #-}
3	{-# LANGUAGE ExistentialQuantification #-}
4	{-# LANGUAGE ScopedTypeVariables #-}
5	--------------------------------------------------------------------------
6	-- \|
7	-- Module : ToySolver.Graph.ShortestPath
8	-- Copyright : (c) Masahiro Sakai 2016-2017
9	-- License : BSD-style
10	--
11	-- Maintainer : masahiro.sakai@gmail.com
12	-- Stability : provisional
13	-- Portability : non-portable
14	--
15	-- This module provides functions for shotest path computation.
16	--
17	-- Reference:
18	--
19	-- * Friedrich Eisenbrand. “Linear and Discrete Optimization”.
20	-- <https://www.coursera.org/course/linearopt>
21	--
22	--------------------------------------------------------------------------
23	module ToySolver.Graph.ShortestPath
24	(
25	-- * Graph data types
26	Graph
27	, Edge
28	, OutEdge
29	, InEdge
30
31	-- * Fold data type
32	, Fold (..)
33	, monoid'
34	, monoid
35	, unit
36	, pair
37	, path
38	, firstOutEdge
39	, lastInEdge
40	, cost
41
42	-- * Path data types
43	, Path (..)
44	, pathFrom
45	, pathTo
46	, pathCost
47	, pathEmpty
48	, pathAppend
49	, pathEdges
50	, pathEdgesBackward
51	, pathEdgesSeq
52	, pathVertexes
53	, pathVertexesBackward
54	, pathVertexesSeq
55	, pathFold
56	, pathMin
57
58	-- * Shortest-path algorithms
59	, bellmanFord
60	, dijkstra
61	, floydWarshall
62
63	-- * Utility functions
64	, bellmanFordDetectNegativeCycle
65	) where
66
67	import Control.Monad
68	import Control.Monad.ST
69	import Control.Monad.Trans
70	import Control.Monad.Trans.Except
71	import Data.Hashable
72	import qualified Data.HashTable.Class as H
73	import qualified Data.HashTable.ST.Cuckoo as C
74	import Data.IntMap.Strict (IntMap)
75	import qualified Data.IntMap.Strict as IntMap
76	import qualified Data.IntSet as IntSet
77	import qualified Data.Heap as Heap -- http://hackage.haskell.org/package/heaps
78	import Data.List (foldl')
79	import Data.Maybe (fromJust)
80	import Data.Monoid
81	import Data.Ord
82	import Data.Sequence (Seq)
83	import qualified Data.Sequence as Seq
84	import Data.STRef
85
86	-- ------------------------------------------------------------------------
87
88	-- \| Graph represented as a map from vertexes to their outgoing edges
89	type Graph cost label = IntMap [OutEdge cost label]
90
91	-- \| Vertex data type
92	type Vertex = Int
93
94	-- \| Edge data type
95	type Edge cost label = (Vertex, Vertex, cost, label)
96
97	-- \| Outgoing edge data type (source vertex is implicit)
98	type OutEdge cost label = (Vertex, cost, label)
99
100	-- \| Incoming edge data type (target vertex is implicit)
101	type InEdge cost label = (Vertex, cost, label)
102
103	-- \| path data type.
104	data Path cost label
105	= Empty Vertex
106	-- ^ empty path
107	\| Singleton (Edge cost label)
108	-- ^ path with single edge
109	\| Append (Path cost label) (Path cost label) !cost
110	-- ^ concatenation of two paths
111	deriving (Eq, Show)	×
112
113	-- \| Source vertex
114	pathFrom :: Path cost label -> Vertex
115	pathFrom (Empty v) = v	1✔
116	pathFrom (Singleton (from,_,_,_)) = from	1✔
117	pathFrom (Append p1 _ _) = pathFrom p1	1✔
118
119	-- \| Target vertex
120	pathTo :: Path cost label -> Vertex
121	pathTo (Empty v) = v	1✔
122	pathTo (Singleton (_,to,_,_)) = to	1✔
123	pathTo (Append _ p2 _) = pathTo p2	1✔
124
125	-- \| Cost of a path
126	pathCost :: Num cost => Path cost label -> cost
127	pathCost (Empty _) = 0	1✔
128	pathCost (Singleton (_,_,c,_)) = c	1✔
129	pathCost (Append _ _ c) = c	1✔
130
131	-- \| Empty path
132	pathEmpty :: Vertex -> Path cost label
133	pathEmpty = Empty	1✔
134
135	-- \| Concatenation of two paths
136	pathAppend :: (Num cost) => Path cost label -> Path cost label -> Path cost label
137	pathAppend p1 p2	1✔
138	\| pathTo p1 /= pathFrom p2 = error "ToySolver.Graph.ShortestPath.pathAppend: pathTo/pathFrom mismatch"	×
139	\| otherwise =	×
140	case (p1, p2) of	1✔
141	(Empty _, _) -> p2	1✔
142	(_, Empty _) -> p1	1✔
143	_ -> Append p1 p2 (pathCost p1 + pathCost p2)	1✔
144
145	-- \| Edges of a path
146	pathEdges :: Path cost label -> [Edge cost label]
147	pathEdges p = f p []	×
148	where
149	f (Empty _) xs = xs	×
150	f (Singleton e) xs = e : xs	×
151	f (Append p1 p2 _) xs = f p1 (f p2 xs)	×
152
153	-- \| Edges of a path, but in the reverse order
154	pathEdgesBackward :: Path cost label -> [Edge cost label]
155	pathEdgesBackward p = f p []	×
156	where
157	f (Empty _) xs = xs	×
158	f (Singleton e) xs = e : xs	×
159	f (Append p1 p2 _) xs = f p2 (f p1 xs)	×
160
161	-- \| Edges of a path, but as `Seq`
162	pathEdgesSeq :: Path cost label -> Seq (Edge cost label)
163	pathEdgesSeq (Empty _) = Seq.empty	×
164	pathEdgesSeq (Singleton e) = Seq.singleton e	×
165	pathEdgesSeq (Append p1 p2 _) = pathEdgesSeq p1 <> pathEdgesSeq p2	×
166
167	-- \| Vertexes of a path
168	pathVertexes :: Path cost label -> [Vertex]
169	pathVertexes p = pathFrom p : f p []	1✔
170	where
171	f (Empty _) xs = xs	×
172	f (Singleton (_,v2,_,_)) xs = v2 : xs	1✔
173	f (Append p1 p2 _) xs = f p1 (f p2 xs)	1✔
174
175	-- \| Vertexes of a path, but in the reverse order
176	pathVertexesBackward :: Path cost label -> [Vertex]
177	pathVertexesBackward p = pathTo p : f p []	×
178	where
179	f (Empty _) xs = xs	×
180	f (Singleton (v1,_,_,_)) xs = v1 : xs	×
181	f (Append p1 p2 _) xs = f p2 (f p1 xs)	×
182
183	-- \| Vertex of a path, but as `Seq`
184	pathVertexesSeq :: Path cost label -> Seq Vertex
185	pathVertexesSeq p = f True p	×
186	where
187	f True (Empty v) = Seq.singleton v	×
188	f False (Empty _) = mempty	×
189	f True (Singleton (v1,v2,_,_)) = Seq.fromList [v1, v2]	×
190	f False (Singleton (v1,_,_,_)) = Seq.singleton v1	×
191	f b (Append p1 p2 _) = f False p1 <> f b p2	×
192
193	pathMin :: (Num cost, Ord cost) => Path cost label -> Path cost label -> Path cost label
194	pathMin p1 p2	×
195	\| pathCost p1 <= pathCost p2 = p1	×
196	\| otherwise = p2	×
197
198	-- \| Fold a path using a given 'Fold` operation.
199	pathFold :: Fold cost label a -> Path cost label -> a
200	pathFold (Fold fV fE fC fD) p = fD (h p)	×
201	where
202	h (Empty v) = fV v	×
203	h (Singleton e) = fE e	×
204	h (Append p1 p2 _) = fC (h p1) (h p2)	×
205
206	-- ------------------------------------------------------------------------
207
208	-- \| Strict pair type
209	data Pair a b = Pair !a !b
210
211	-- ------------------------------------------------------------------------
212
213	-- \| Operations for folding edge information along with a path into a @r@ value.
214	--
215	-- @Fold cost label r@ consists of three operations
216	--
217	-- * @Vertex -> a@ corresponds to an empty path,
218	--
219	-- * @Edge cost label -> a@ corresponds to a single edge,
220	--
221	-- * @a -> a -> a@ corresponds to path concatenation,
222	--
223	-- and @a -> r@ to finish the computation.
224	data Fold cost label r
225	= forall a. Fold (Vertex -> a) (Edge cost label -> a) (a -> a -> a) (a -> r)
226
227	instance Functor (Fold cost label) where	×
228	{-# INLINE fmap #-}
229	fmap f (Fold fV fE fC fD) = Fold fV fE fC (f . fD)	×
230
231	instance Applicative (Fold cost label) where	×
232	{-# INLINE pure #-}
233	pure a = Fold (\_ -> ()) (\_ -> ()) (\_ _ -> ()) (const a)	×
234
235	{-# INLINE (<*>) #-}
236	Fold fV1 fE1 fC1 fD1 <*> Fold fV2 fE2 fC2 fD2 =	×
237	Fold (\v -> Pair (fV1 v) (fV2 v))	×
238	(\e -> Pair (fE1 e) (fE2 e))	×
239	(\(Pair a1 b1) (Pair a2 b2) -> Pair (fC1 a1 a2) (fC2 b1 b2))	×
240	(\(Pair a b) -> fD1 a (fD2 b))	×
241
242	-- \| Project `Edge` into a monoid value and fold using monoidal operation.
243	monoid' :: Monoid m => (Edge cost label -> m) -> Fold cost label m
244	monoid' f = Fold (\_ -> mempty) f mappend id	1✔
245
246	-- \| Similar to 'monoid'' but a /label/ is directly used as a monoid value.
247	monoid :: Monoid m => Fold cost m m
248	monoid = monoid' (\(_,_,_,m) -> m)	×
249
250	-- \| Ignore contents and return a unit value.
251	unit :: Fold cost label ()
252	unit = monoid' (\_ -> ())	×
253
254	-- \| Pairs two `Fold` into one that produce a tuple.
255	pair :: Fold cost label a -> Fold cost label b -> Fold cost label (a,b)
256	pair (Fold fV1 fE1 fC1 fD1) (Fold fV2 fE2 fC2 fD2) =	×
257	Fold (\v -> Pair (fV1 v) (fV2 v))	×
258	(\e -> Pair (fE1 e) (fE2 e))	×
259	(\(Pair a1 b1) (Pair a2 b2) -> Pair (fC1 a1 a2) (fC2 b1 b2))	×
260	(\(Pair a b) -> (fD1 a, fD2 b))	×
261
262	-- \| Construct a `Path` value.
263	path :: (Num cost) => Fold cost label (Path cost label)
264	path = Fold pathEmpty Singleton pathAppend id	1✔
265
266	-- \| Compute cost of a path.
267	cost :: Num cost => Fold cost label cost
268	cost = Fold (\_ -> 0) (\(_,_,c,_) -> c) (+) id	×
269
270	-- \| Get the first `OutEdge` of a path.
271	firstOutEdge :: Fold cost label (First (OutEdge cost label))
272	firstOutEdge = monoid' (\(_,v,c,l) -> First (Just (v,c,l)))	×
273
274	-- \| Get the last `InEdge` of a path.
275	-- This is useful for constructing a /parent/ map of a spanning tree.
276	lastInEdge :: Fold cost label (Last (InEdge cost label))
277	lastInEdge = monoid' (\(v,_,c,l) -> Last (Just (v,c,l)))	1✔
278
279	-- ------------------------------------------------------------------------
280
281	-- \| Bellman-Ford algorithm for finding shortest paths from source vertexes
282	-- to all of the other vertices in a weighted graph with negative weight
283	-- edges allowed.
284	--
285	-- It compute shortest-paths from given source vertexes, and folds edge
286	-- information along shortest paths using a given 'Fold' operation.
287	bellmanFord
288	:: Real cost
289	=> Fold cost label a
290	-- ^ Operation used to fold shotest paths
291	-> Graph cost label
292	-> [Vertex]
293	-- ^ List of source vertexes @vs@
294	-> IntMap (cost, a)
295	bellmanFord (Fold fV fE fC fD) g ss = runST $ do	1✔
296	let n = IntMap.size g	1✔
297	d <- C.newSized n	1✔
298	forM_ ss $ \s -> H.insert d s (Pair 0 (fV s))	1✔
299
300	updatedRef <- newSTRef (IntSet.fromList ss)	1✔
301	_ <- runExceptT $ replicateM_ n $ do	1✔
302	us <- lift $ readSTRef updatedRef	1✔
303	when (IntSet.null us) $ throwE ()	×
304	lift $ do	1✔
305	writeSTRef updatedRef IntSet.empty	1✔
306	forM_ (IntSet.toList us) $ \u -> do	1✔
307	-- modifySTRef' updatedRef (IntSet.delete u) -- possible optimization
308	Pair du a <- liftM fromJust $ H.lookup d u	1✔
309	forM_ (IntMap.findWithDefault [] u g) $ \(v, c, l) -> do	×
310	m <- H.lookup d v	1✔
311	case m of	1✔
312	Just (Pair c0 _) \| c0 <= du + c -> return ()	×
313	_ -> do	1✔
314	H.insert d v (Pair (du + c) (a `fC` fE (u,v,c,l)))	1✔
315	modifySTRef' updatedRef (IntSet.insert v)	1✔
316
317	liftM (fmap (\(Pair c x) -> (c, fD x))) $ freezeHashTable d	1✔
318
319	freezeHashTable :: H.HashTable h => h s Int v -> ST s (IntMap v)
320	freezeHashTable h = H.foldM (\m (k,v) -> return $! IntMap.insert k v m) IntMap.empty h	1✔
321
322	-- \| Utility function for detecting a negative cycle.
323	bellmanFordDetectNegativeCycle
324	:: forall cost label a. Real cost
325	=> Fold cost label a
326	-- ^ Operation used to fold a cycle
327	-> Graph cost label
328	-> IntMap (cost, Last (InEdge cost label))
329	-- ^ Result of @'bellmanFord' 'lastInEdge' vs@
330	-> Maybe a
331	bellmanFordDetectNegativeCycle (Fold fV fE fC fD) g d = either (Just . fD) (const Nothing) $ do	1✔
332	forM_ (IntMap.toList d) $ \(u,(du,_)) ->	1✔
333	forM_ (IntMap.findWithDefault [] u g) $ \(v, c, l) -> do	×
334	let (dv,_) = d IntMap.! v	1✔
335	when (du + c < dv) $ do	1✔
336	-- a negative cycle is detected
337	let d' = IntMap.insert v (du + c, Last (Just (u, c, l))) d	×
338	parent u = do	1✔
339	case IntMap.lookup u d' of	1✔
340	Just (_, Last (Just (v,_,_))) -> v	1✔
341	_ -> undefined	×
342	u0 = go (parent (parent v)) (parent v)	1✔
343	where
344	go hare tortoise	1✔
345	\| hare == tortoise = hare	1✔
346	\| otherwise = go (parent (parent hare)) (parent tortoise)	×
347	let go u result = do	1✔
348	let Just (_, Last (Just (v,c,l))) = IntMap.lookup u d'	1✔
349	if v == u0 then	1✔
350	fC (fE (v,u,c,l)) result	1✔
351	else
352	go v (fC (fE (v,u,c,l)) result)	1✔
353	Left $ go u0 (fV u0)	1✔
354
355	-- ------------------------------------------------------------------------
356
357	-- \| Dijkstra's algorithm for finding shortest paths from source vertexes
358	-- to all of the other vertices in a weighted graph with non-negative edge
359	-- weight.
360	--
361	-- It compute shortest-paths from given source vertexes, and folds edge
362	-- information along shortest paths using a given 'Fold' operation.
363	dijkstra
364	:: forall cost label a. Real cost
365	=> Fold cost label a
366	-- ^ Operation used to fold shotest paths
367	-> Graph cost label
368	-> [Vertex]
369	-- ^ List of source vertexes
370	-> IntMap (cost, a)
371	dijkstra (Fold fV fE fC (fD :: x -> a)) g ss =	1✔
372	fmap (\(Pair c x) -> (c, fD x)) $	×
373	loop (Heap.fromList [Heap.Entry 0 (Pair s (fV s)) \| s <- ss]) IntMap.empty	1✔
374	where
375	loop
376	:: Heap.Heap (Heap.Entry cost (Pair Vertex x))
377	-> IntMap (Pair cost x)
378	-> IntMap (Pair cost x)
379	loop q visited =	1✔
380	case Heap.viewMin q of	1✔
381	Nothing -> visited	1✔
382	Just (Heap.Entry c (Pair v a), q')
383	\| v `IntMap.member` visited -> loop q' visited	1✔
384	\| otherwise ->	×
385	let q2 = Heap.fromList	1✔
386	[ Heap.Entry (c+c') (Pair ch (a `fC` fE (v,ch,c',l)))	1✔
387	\| (ch,c',l) <- IntMap.findWithDefault [] v g	×
388	, not (ch `IntMap.member` visited)	1✔
389	]
390	in loop (Heap.union q' q2) (IntMap.insert v (Pair c a) visited)	1✔
391
392	-- ------------------------------------------------------------------------
393
394	-- \| Floyd-Warshall algorithm for finding shortest paths in a weighted graph
395	-- with positive or negative edge weights (but with no negative cycles).
396	--
397	-- It compute shortest-paths between each pair of vertexes, and folds edge
398	-- information along shortest paths using a given 'Fold' operation.
399	floydWarshall
400	:: forall cost label a. Real cost
401	=> Fold cost label a
402	-- ^ Operation used to fold shotest paths
403	-> Graph cost label
404	-> IntMap (IntMap (cost, a))
405	floydWarshall (Fold fV fE fC (fD :: x -> a)) g =	1✔
406	fmap (fmap (\(Pair c x) -> (c, fD x))) $	1✔
407	IntMap.unionWith (IntMap.unionWith minP) (foldl' f tbl0 vs) paths0	1✔
408	where
409	vs = IntMap.keys g	1✔
410
411	paths0 :: IntMap (IntMap (Pair cost x))
412	paths0 = IntMap.mapWithKey (\v _ -> IntMap.singleton v (Pair 0 (fV v))) g	1✔
413
414	tbl0 :: IntMap (IntMap (Pair cost x))
415	tbl0 = IntMap.mapWithKey (\v es -> IntMap.fromListWith minP [(u, (Pair c (fE (v,u,c,l)))) \| (u,c,l) <- es]) g	1✔
416
417	minP :: Pair cost x -> Pair cost x -> Pair cost x
418	minP = minBy (comparing (\(Pair c _) -> c))	1✔
419
420	f :: IntMap (IntMap (Pair cost x))
421	-> Vertex
422	-> IntMap (IntMap (Pair cost x))
423	f tbl vk =	1✔
424	case IntMap.lookup vk tbl of	1✔
425	Nothing -> tbl	×
426	Just hk -> IntMap.map h tbl	1✔
427	where
428	h :: IntMap (Pair cost x) -> IntMap (Pair cost x)
429	h m =	1✔
430	case IntMap.lookup vk m of	1✔
431	Nothing -> m	1✔
432	Just (Pair c1 x1) -> IntMap.unionWith minP m (IntMap.map (\(Pair c2 x2) -> (Pair (c1+c2) (fC x1 x2))) hk)	1✔
433
434	minBy :: (a -> a -> Ordering) -> a -> a -> a
435	minBy f a b =	1✔
436	case f a b of	1✔
437	LT -> a	1✔
438	EQ -> a	1✔
439	GT -> b	1✔
440
441	-- ------------------------------------------------------------------------

msakai / toysolver / 496

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous