{-# LANGUAGE LambdaCase #-} ----------------------------------------------------------------------------- -- | -- Module : Algebra.Graph.AdjacencyMap.Algorithm -- Copyright : (c) Andrey Mokhov 2016-2025 -- License : MIT (see the file LICENSE) -- Maintainer : andrey.mokhov@gmail.com -- Stability : unstable -- -- __Alga__ is a library for algebraic construction and manipulation of graphs -- in Haskell. See for the -- motivation behind the library, the underlying theory, and implementation details. -- -- This module provides basic graph algorithms, such as /depth-first search/, -- implemented for the "Algebra.Graph.AdjacencyMap" data type. ----------------------------------------------------------------------------- module Algebra.Graph.AdjacencyMap.Algorithm ( -- * Algorithms bfsForest, bfs, dfsForest, dfsForestFrom, dfs, reachable, topSort, isAcyclic, scc, -- * Correctness properties isDfsForestOf, isTopSortOf, -- * Type synonyms Cycle ) where import Control.Monad import Control.Monad.Trans.Cont import Control.Monad.Trans.State.Strict import Data.Foldable (for_) import Data.Either import Data.List.NonEmpty (NonEmpty(..), (<|)) import Data.Maybe import Data.Tree (Forest, Tree (..), flatten, levels, unfoldForestM_BF) import Algebra.Graph.AdjacencyMap import qualified Algebra.Graph.NonEmpty.AdjacencyMap as NonEmpty import qualified Data.Array as Array import qualified Data.List as List import qualified Data.Map.Strict as Map import qualified Data.Set as Set -- | Compute the /breadth-first search/ forest of a graph, such that adjacent -- vertices are explored in increasing order according to their 'Ord' instance. -- The search is seeded by a list of vertices that will become the roots of the -- resulting forest. Duplicates in the list will have their first occurrence -- explored and subsequent ones ignored. The seed vertices that do not belong to -- the graph are also ignored. -- -- Complexity: /O((L + m) * log n)/ time and /O(n)/ space, where /L/ is the -- number of seed vertices. -- -- @ -- 'forest' $ bfsForest ('edge' 1 2) [0] == 'empty' -- 'forest' $ bfsForest ('edge' 1 2) [1] == 'edge' 1 2 -- 'forest' $ bfsForest ('edge' 1 2) [2] == 'vertex' 2 -- 'forest' $ bfsForest ('edge' 1 2) [0,1,2] == 'vertices' [1,2] -- 'forest' $ bfsForest ('edge' 1 2) [2,1,0] == 'vertices' [1,2] -- 'forest' $ bfsForest ('edge' 1 1) [1] == 'vertex' 1 -- 'isSubgraphOf' ('forest' $ bfsForest x vs) x == True -- bfsForest x ('vertexList' x) == 'map' (\\v -> Node v []) ('Data.List.nub' $ 'vertexList' x) -- bfsForest x [] == [] -- bfsForest 'empty' vs == [] -- bfsForest (3 * (1 + 4) * (1 + 5)) [1,4] == [ Node { rootLabel = 1 -- , subForest = [ Node { rootLabel = 5 -- , subForest = [] }]} -- , Node { rootLabel = 4 -- , subForest = [] }] -- 'forest' $ bfsForest ('circuit' [1..5] + 'circuit' [5,4..1]) [3] == 'path' [3,2,1] + 'path' [3,4,5] -- -- @ bfsForest :: Ord a => AdjacencyMap a -> [a] -> Forest a bfsForest x vs = evalState (explore [ v | v <- vs, hasVertex v x ]) Set.empty where explore = filterM discovered >=> unfoldForestM_BF walk walk v = (v,) <$> adjacentM v adjacentM v = filterM discovered $ Set.toList (postSet v x) discovered v = do new <- gets (not . Set.member v) when new $ modify' (Set.insert v) return new -- | A version of 'bfsForest' where the resulting forest is converted to a level -- structure. Adjacent vertices are explored in the increasing order according -- to their 'Ord' instance. Flattening the result via @'concat'@ @.@ @'bfs'@ @x@ -- gives an enumeration of reachable vertices in the breadth-first search order. -- -- Complexity: /O((L + m) * min(n,W))/ time and /O(n)/ space, where /L/ is the -- number of seed vertices. -- -- @ -- bfs ('edge' 1 2) [0] == [] -- bfs ('edge' 1 2) [1] == [[1], [2]] -- bfs ('edge' 1 2) [2] == [[2]] -- bfs ('edge' 1 2) [1,2] == [[1,2]] -- bfs ('edge' 1 2) [2,1] == [[2,1]] -- bfs ('edge' 1 1) [1] == [[1]] -- bfs 'empty' vs == [] -- bfs x [] == [] -- bfs (1 * 2 + 3 * 4 + 5 * 6) [1,2] == [[1,2]] -- bfs (1 * 2 + 3 * 4 + 5 * 6) [1,3] == [[1,3], [2,4]] -- bfs (3 * (1 + 4) * (1 + 5)) [3] == [[3], [1,4,5]] -- -- bfs ('circuit' [1..5] + 'circuit' [5,4..1]) [3] == [[2], [1,3], [5,4]] -- 'concat' $ bfs ('circuit' [1..5] + 'circuit' [5,4..1]) [3] == [3,2,4,1,5] -- 'map' 'concat' . 'List.transpose' . 'map' 'levels' . 'bfsForest' x == bfs x -- @ bfs :: Ord a => AdjacencyMap a -> [a] -> [[a]] bfs x = map concat . List.transpose . map levels . bfsForest x dfsForestFromImpl :: Ord a => AdjacencyMap a -> [a] -> Forest a dfsForestFromImpl g vs = evalState (explore vs) Set.empty where explore (v:vs) = discovered v >>= \case True -> (:) <$> walk v <*> explore vs False -> explore vs explore [] = return [] walk v = Node v <$> explore (adjacent v) adjacent v = Set.toList (postSet v g) discovered v = do new <- gets (not . Set.member v) when new $ modify' (Set.insert v) return new -- | Compute the /depth-first search/ forest of a graph, where adjacent vertices -- are explored in the increasing order according to their 'Ord' instance. -- -- Complexity: /O((n + m) * min(n,W))/ time and /O(n)/ space. -- -- @ -- 'forest' $ dfsForest 'empty' == 'empty' -- 'forest' $ dfsForest ('edge' 1 1) == 'vertex' 1 -- 'forest' $ dfsForest ('edge' 1 2) == 'edge' 1 2 -- 'forest' $ dfsForest ('edge' 2 1) == 'vertices' [1,2] -- 'isSubgraphOf' ('forest' $ dfsForest x) x == True -- 'isDfsForestOf' (dfsForest x) x == True -- dfsForest . 'forest' . dfsForest == dfsForest -- dfsForest ('vertices' vs) == 'map' (\\v -> Node v []) ('Data.List.nub' $ 'Data.List.sort' vs) -- dfsForest $ 3 * (1 + 4) * (1 + 5) == [ Node { rootLabel = 1 -- , subForest = [ Node { rootLabel = 5 -- , subForest = [] }]} -- , Node { rootLabel = 3 -- , subForest = [ Node { rootLabel = 4 -- , subForest = [] }]}] -- 'forest' (dfsForest $ 'circuit' [1..5] + 'circuit' [5,4..1]) == 'path' [1,2,3,4,5] -- @ dfsForest :: Ord a => AdjacencyMap a -> Forest a dfsForest g = dfsForestFromImpl g (vertexList g) -- | Compute the /depth-first search/ forest of a graph starting from the given -- seed vertices, where adjacent vertices are explored in the increasing order -- according to their 'Ord' instance. Note that the resulting forest does not -- necessarily span the whole graph, as some vertices may be unreachable. The -- seed vertices which do not belong to the graph are ignored. -- -- Complexity: /O((L + m) * log n)/ time and /O(n)/ space, where /L/ is the -- number of seed vertices. -- -- @ -- 'forest' $ dfsForestFrom 'empty' vs == 'empty' -- 'forest' $ dfsForestFrom ('edge' 1 1) [1] == 'vertex' 1 -- 'forest' $ dfsForestFrom ('edge' 1 2) [0] == 'empty' -- 'forest' $ dfsForestFrom ('edge' 1 2) [1] == 'edge' 1 2 -- 'forest' $ dfsForestFrom ('edge' 1 2) [2] == 'vertex' 2 -- 'forest' $ dfsForestFrom ('edge' 1 2) [1,2] == 'edge' 1 2 -- 'forest' $ dfsForestFrom ('edge' 1 2) [2,1] == 'vertices' [1,2] -- 'isSubgraphOf' ('forest' $ dfsForestFrom x vs) x == True -- 'isDfsForestOf' (dfsForestFrom x ('vertexList' x)) x == True -- dfsForestFrom x ('vertexList' x) == 'dfsForest' x -- dfsForestFrom x [] == [] -- dfsForestFrom (3 * (1 + 4) * (1 + 5)) [1,4] == [ Node { rootLabel = 1 -- , subForest = [ Node { rootLabel = 5 -- , subForest = [] } -- , Node { rootLabel = 4 -- , subForest = [] }] -- 'forest' $ dfsForestFrom ('circuit' [1..5] + 'circuit' [5,4..1]) [3] == 'path' [3,2,1,5,4] -- @ dfsForestFrom :: Ord a => AdjacencyMap a -> [a] -> Forest a dfsForestFrom g vs = dfsForestFromImpl g [ v | v <- vs, hasVertex v g ] -- | Return the list vertices visited by the /depth-first search/ in a graph, -- starting from the given seed vertices. Adjacent vertices are explored in the -- increasing order according to their 'Ord' instance. -- -- Complexity: /O((L + m) * log n)/ time and /O(n)/ space, where /L/ is the -- number of seed vertices. -- -- @ -- dfs 'empty' vs == [] -- dfs ('edge' 1 1) [1] == [1] -- dfs ('edge' 1 2) [0] == [] -- dfs ('edge' 1 2) [1] == [1,2] -- dfs ('edge' 1 2) [2] == [2] -- dfs ('edge' 1 2) [1,2] == [1,2] -- dfs ('edge' 1 2) [2,1] == [2,1] -- dfs x [] == [] -- -- 'Data.List.and' [ 'hasVertex' v x | v <- dfs x vs ] == True -- dfs (3 * (1 + 4) * (1 + 5)) [1,4] == [1,5,4] -- dfs ('circuit' [1..5] + 'circuit' [5,4..1]) [3] == [3,2,1,5,4] -- @ dfs :: Ord a => AdjacencyMap a -> [a] -> [a] dfs x = concatMap flatten . dfsForestFrom x -- | Return the list of vertices /reachable/ from a source vertex in a graph. -- The vertices in the resulting list appear in the /depth-first search order/. -- -- Complexity: /O(m * log n)/ time and /O(n)/ space. -- -- @ -- reachable 'empty' x == [] -- reachable ('vertex' 1) 1 == [1] -- reachable ('edge' 1 1) 1 == [1] -- reachable ('edge' 1 2) 0 == [] -- reachable ('edge' 1 2) 1 == [1,2] -- reachable ('edge' 1 2) 2 == [2] -- reachable ('path' [1..8] ) 4 == [4..8] -- reachable ('circuit' [1..8] ) 4 == [4..8] ++ [1..3] -- reachable ('clique' [8,7..1]) 8 == [8] ++ [1..7] -- -- 'Data.List.and' [ 'hasVertex' v x | v <- reachable x y ] == True -- @ reachable :: Ord a => AdjacencyMap a -> a -> [a] reachable x y = dfs x [y] type Cycle = NonEmpty type Result a = Either (Cycle a) [a] data NodeState = Entered | Exited data S a = S { parent :: Map.Map a a , entry :: Map.Map a NodeState , order :: [a] } topSortImpl :: Ord a => AdjacencyMap a -> StateT (S a) (Cont (Result a)) (Result a) topSortImpl g = liftCallCC' callCC $ \cyclic -> do let vertices = map fst $ Map.toDescList $ adjacencyMap g adjacent = Set.toDescList . flip postSet g dfsRoot x = nodeState x >>= \case Nothing -> enterRoot x >> dfs x >> exit x _ -> return () dfs x = forM_ (adjacent x) $ \y -> nodeState y >>= \case Nothing -> enter x y >> dfs y >> exit y Just Exited -> return () Just Entered -> cyclic . Left . retrace x y =<< gets parent forM_ vertices dfsRoot Right <$> gets order where nodeState v = gets (Map.lookup v . entry) enter u v = modify' (\(S m n vs) -> S (Map.insert v u m) (Map.insert v Entered n) vs) enterRoot v = modify' (\(S m n vs) -> S m (Map.insert v Entered n) vs) exit v = modify' (\(S m n vs) -> S m (Map.alter (fmap leave) v n) (v:vs)) where leave = \case Entered -> Exited Exited -> error "Internal error: dfs search order violated" retrace curr head parent = aux (curr :| []) where aux xs@(curr :| _) | head == curr = xs | otherwise = aux (parent Map.! curr <| xs) -- | Compute a topological sort of a graph or discover a cycle. -- -- Vertices are explored in the decreasing order according to their 'Ord' -- instance. This gives the lexicographically smallest topological ordering in -- the case of success. In the case of failure, the cycle is characterized by -- being the lexicographically smallest up to rotation with respect to -- @Ord@ @(Dual@ @Int)@ in the first connected component of the graph containing -- a cycle, where the connected components are ordered by their largest vertex -- with respect to @Ord a@. -- -- Complexity: /O((n + m) * min(n,W))/ time and /O(n)/ space. -- -- @ -- topSort (1 * 2 + 3 * 1) == Right [3,1,2] -- topSort ('path' [1..5]) == Right [1..5] -- topSort (3 * (1 * 4 + 2 * 5)) == Right [3,1,2,4,5] -- topSort (1 * 2 + 2 * 1) == Left (2 ':|' [1]) -- topSort ('path' [5,4..1] + 'edge' 2 4) == Left (4 ':|' [3,2]) -- topSort ('circuit' [1..3]) == Left (3 ':|' [1,2]) -- topSort ('circuit' [1..3] + 'circuit' [3,2,1]) == Left (3 ':|' [2]) -- topSort (1 * 2 + (5 + 2) * 1 + 3 * 4 * 3) == Left (1 ':|' [2]) -- fmap ('flip' 'isTopSortOf' x) (topSort x) /= Right False -- 'isRight' . topSort == 'isAcyclic' -- topSort . 'vertices' == Right . 'nub' . 'sort' -- @ topSort :: Ord a => AdjacencyMap a -> Either (Cycle a) [a] topSort g = runCont (evalStateT (topSortImpl g) initialState) id where initialState = S Map.empty Map.empty [] -- | Check if a given graph is /acyclic/. -- -- Complexity: /O((n+m)*log n)/ time and /O(n)/ space. -- -- @ -- isAcyclic (1 * 2 + 3 * 1) == True -- isAcyclic (1 * 2 + 2 * 1) == False -- isAcyclic . 'circuit' == 'null' -- isAcyclic == 'isRight' . 'topSort' -- @ isAcyclic :: Ord a => AdjacencyMap a -> Bool isAcyclic = isRight . topSort -- | Compute the /condensation/ of a graph, where each vertex corresponds to a -- /strongly-connected component/ of the original graph. Note that component -- graphs are non-empty, and are therefore of type -- "Algebra.Graph.NonEmpty.AdjacencyMap". -- -- Details about the implementation can be found at -- . -- -- Complexity: /O((n+m)*log n)/ time and /O(n+m)/ space. -- -- @ -- scc 'empty' == 'empty' -- scc ('vertex' x) == 'vertex' (NonEmpty.'NonEmpty.vertex' x) -- scc ('vertices' xs) == 'vertices' ('map' 'NonEmpty.vertex' xs) -- scc ('edge' 1 1) == 'vertex' (NonEmpty.'NonEmpty.edge' 1 1) -- scc ('edge' 1 2) == 'edge' (NonEmpty.'NonEmpty.vertex' 1) (NonEmpty.'NonEmpty.vertex' 2) -- scc ('circuit' (1:xs)) == 'vertex' (NonEmpty.'NonEmpty.circuit1' (1 'Data.List.NonEmpty.:|' xs)) -- scc (3 * 1 * 4 * 1 * 5) == 'edges' [ (NonEmpty.'NonEmpty.vertex' 3 , NonEmpty.'NonEmpty.vertex' 5 ) -- , (NonEmpty.'NonEmpty.vertex' 3 , NonEmpty.'NonEmpty.clique1' [1,4,1]) -- , (NonEmpty.'NonEmpty.clique1' [1,4,1], NonEmpty.'NonEmpty.vertex' 5 ) ] -- 'isAcyclic' . scc == 'const' True -- 'isAcyclic' x == (scc x == 'gmap' NonEmpty.'NonEmpty.vertex' x) -- @ scc :: Ord a => AdjacencyMap a -> AdjacencyMap (NonEmpty.AdjacencyMap a) scc g = condense g $ execState (gabowSCC g) initialState where initialState = SCC 0 0 [] [] Map.empty Map.empty [] [] [] data StateSCC a = SCC { _preorder :: {-# unpack #-} !Int , _component :: {-# unpack #-} !Int , boundaryStack :: [(Int,a)] , _pathStack :: [a] , preorders :: Map.Map a Int , components :: Map.Map a Int , _innerGraphs :: [AdjacencyMap a] , _innerEdges :: [(Int,(a,a))] , _outerEdges :: [(a,a)] } deriving (Show) gabowSCC :: Ord a => AdjacencyMap a -> State (StateSCC a) () gabowSCC g = do let dfs u = do p_u <- enter u for_ (postSet u g) $ \v -> do preorderId v >>= \case Nothing -> do updated <- dfs v if updated then outedge (u,v) else inedge (p_u,(u,v)) Just p_v -> do scc_v <- hasComponent v if scc_v then outedge (u,v) else popBoundary p_v >> inedge (p_u,(u,v)) exit u forM_ (vertexList g) $ \v -> do assigned <- hasPreorderId v unless assigned $ void $ dfs v where -- called when visiting vertex v. assigns preorder number to v, -- adds the (id, v) pair to the boundary stack b, and adds v to -- the path stack s. enter v = do SCC pre scc bnd pth pres sccs gs es_i es_o <- get let pre' = pre+1 bnd' = (pre,v):bnd pth' = v:pth pres' = Map.insert v pre pres put $! SCC pre' scc bnd' pth' pres' sccs gs es_i es_o return pre -- called on back edges. pops the boundary stack while the top -- vertex has a larger preorder number than p_v. popBoundary p_v = modify' (\(SCC pre scc bnd pth pres sccs gs es_i es_o) -> SCC pre scc (dropWhile ((>p_v).fst) bnd) pth pres sccs gs es_i es_o) -- called when exiting vertex v. if v is the bottom of a scc -- boundary, we add a new SCC, otherwise v is part of a larger scc -- being constructed and we continue. exit v = do boundaryStack <- gets boundaryStack case boundaryStack of (p_top, top) : newBoundaryStack | v == top -> do insertComponent p_top top newBoundaryStack return True _ -> return False insertComponent p_v v newBoundaryStack = modify' (\(SCC pre scc _oldBoundaryStack pth pres sccs gs es_i es_o) -> let (curr,v_pth') = span (/=v) pth pth' = drop 1 v_pth' -- Here we know that v_pth' starts with v (es,es_i') = span ((>=p_v).fst) es_i g_i | null es = vertex v | otherwise = edges (snd <$> es) scc' = scc + 1 sccs' = List.foldl' (\sccs x -> Map.insert x scc sccs) sccs (v:curr) gs' = g_i:gs in SCC pre scc' newBoundaryStack pth' pres sccs' gs' es_i' es_o) inedge uv = modify' (\(SCC pre scc bnd pth pres sccs gs es_i es_o) -> SCC pre scc bnd pth pres sccs gs (uv:es_i) es_o) outedge uv = modify' (\(SCC pre scc bnd pth pres sccs gs es_i es_o) -> SCC pre scc bnd pth pres sccs gs es_i (uv:es_o)) hasPreorderId v = gets (Map.member v . preorders) preorderId v = gets (Map.lookup v . preorders) hasComponent v = gets (Map.member v . components) condense :: Ord a => AdjacencyMap a -> StateSCC a -> AdjacencyMap (NonEmpty.AdjacencyMap a) condense g (SCC _ n _ _ _ assignment inner _ outer) | n == 1 = vertex $ convert g | otherwise = gmap (\c -> inner' Array.! (n-1-c)) outer' where inner' = Array.listArray (0,n-1) (convert <$> inner) outer' = es `overlay` vs vs = vertices [0..n-1] es = edges [ (sccid x, sccid y) | (x,y) <- outer ] sccid v = assignment Map.! v convert = fromJust . NonEmpty.toNonEmpty -- | Check if a given forest is a correct /depth-first search/ forest of a graph. -- The implementation is based on the paper "Depth-First Search and Strong -- Connectivity in Coq" by François Pottier. -- -- @ -- isDfsForestOf [] 'empty' == True -- isDfsForestOf [] ('vertex' 1) == False -- isDfsForestOf [Node 1 []] ('vertex' 1) == True -- isDfsForestOf [Node 1 []] ('vertex' 2) == False -- isDfsForestOf [Node 1 [], Node 1 []] ('vertex' 1) == False -- isDfsForestOf [Node 1 []] ('edge' 1 1) == True -- isDfsForestOf [Node 1 []] ('edge' 1 2) == False -- isDfsForestOf [Node 1 [], Node 2 []] ('edge' 1 2) == False -- isDfsForestOf [Node 2 [], Node 1 []] ('edge' 1 2) == True -- isDfsForestOf [Node 1 [Node 2 []]] ('edge' 1 2) == True -- isDfsForestOf [Node 1 [], Node 2 []] ('vertices' [1,2]) == True -- isDfsForestOf [Node 2 [], Node 1 []] ('vertices' [1,2]) == True -- isDfsForestOf [Node 1 [Node 2 []]] ('vertices' [1,2]) == False -- isDfsForestOf [Node 1 [Node 2 [Node 3 []]]] ('path' [1,2,3]) == True -- isDfsForestOf [Node 1 [Node 3 [Node 2 []]]] ('path' [1,2,3]) == False -- isDfsForestOf [Node 3 [], Node 1 [Node 2 []]] ('path' [1,2,3]) == True -- isDfsForestOf [Node 2 [Node 3 []], Node 1 []] ('path' [1,2,3]) == True -- isDfsForestOf [Node 1 [], Node 2 [Node 3 []]] ('path' [1,2,3]) == False -- @ isDfsForestOf :: Ord a => Forest a -> AdjacencyMap a -> Bool isDfsForestOf f am = case go Set.empty f of Just seen -> seen == vertexSet am Nothing -> False where go seen [] = Just seen go seen (t:ts) = do let root = rootLabel t guard $ root `Set.notMember` seen guard $ and [ hasEdge root (rootLabel subTree) am | subTree <- subForest t ] newSeen <- go (Set.insert root seen) (subForest t) guard $ postSet root am `Set.isSubsetOf` newSeen go newSeen ts -- | Check if a given list of vertices is a correct /topological sort/ of a graph. -- -- @ -- isTopSortOf [3,1,2] (1 * 2 + 3 * 1) == True -- isTopSortOf [1,2,3] (1 * 2 + 3 * 1) == False -- isTopSortOf [] (1 * 2 + 3 * 1) == False -- isTopSortOf [] 'empty' == True -- isTopSortOf [x] ('vertex' x) == True -- isTopSortOf [x] ('edge' x x) == False -- @ isTopSortOf :: Ord a => [a] -> AdjacencyMap a -> Bool isTopSortOf xs m = go Set.empty xs where go seen [] = seen == Map.keysSet (adjacencyMap m) go seen (v:vs) = postSet v m `Set.intersection` newSeen == Set.empty && go newSeen vs where newSeen = Set.insert v seen