math-grads 0.1.5.1 → 0.1.6.2

raw patch · 11 files changed

+737/−286 lines, 11 filesdep +bimapdep +ilistdep +lensdep ~aesondep ~arraydep ~basePVP ok

version bump matches the API change (PVP)

Dependencies added: bimap, ilist, lens

Dependency ranges changed: aeson, array, base, containers, hspec, linear, matrix, mtl, random, vector

API changes (from Hackage documentation)

Files

• README.md +2/−2
• math-grads.cabal +40/−20
• src/Math/Grads/Algo/Cycles.hs +6/−3
• src/Math/Grads/Algo/Isomorphism.hs +10/−247
• src/Math/Grads/Algo/Isomorphism/RI.hs +127/−0
• src/Math/Grads/Algo/Isomorphism/Types.hs +32/−0
• src/Math/Grads/Algo/Isomorphism/Ullman.hs +244/−0
• src/Math/Grads/Algo/SSSR.hs +125/−0
• src/Math/Grads/GenericGraph.hs +14/−5
• test/Isomorphism.hs +45/−9
• test/SSSR.hs +92/−0

@@ -30,9 +30,9 @@  ## Algorithms -### Ullman's subgraph isomorphism algorithm+### Ullmann's subgraph isomorphism algorithm -Math.Grads contains implementation of Ullman's subgraph isomorphism +Math.Grads contains implementation of Ullmann's subgraph isomorphism  [algorithm](https://www.cs.bgu.ac.il/~dinitz/Course/SS-12/Ullman_Algorithm.pdf). There are several functions that one can find helpful in order to check two graphs for isomorphism or subgraph isomorphism: 

@@ -1,6 +1,6 @@ cabal-version: >=1.10 name: math-grads-version: 0.1.5.1+version: 0.1.6.2 license: BSD3 license-file: LICENSE copyright: 2017 Alexandr Sadovnikov@@ -36,7 +36,11 @@         Math.Grads.Algo.Cycles         Math.Grads.Algo.Interaction         Math.Grads.Algo.Isomorphism+        Math.Grads.Algo.Isomorphism.RI+        Math.Grads.Algo.Isomorphism.Types+        Math.Grads.Algo.Isomorphism.Ullman         Math.Grads.Algo.Paths+        Math.Grads.Algo.SSSR         Math.Grads.Algo.Traversals         Math.Grads.Drawing.Coords         Math.Grads.Graph@@ -56,14 +60,17 @@     default-language: Haskell2010     build-depends:         base >=4.7 && <5,-        aeson <1.5,-        array <0.6,-        containers <0.7,-        linear <1.21,-        matrix <0.4,-        mtl <2.3,-        random <1.2,-        vector <0.13+        aeson >=1.4.2.0 && <1.5,+        array >=0.5.3.0 && <0.6,+        bimap >=0.3.3 && <0.4,+        containers >=0.6.0.1 && <0.7,+        lens ==4.17.*,+        linear >=1.20.8 && <1.21,+        matrix >=0.3.6.1 && <0.4,+        mtl >=2.2.2 && <2.3,+        ilist >=0.3.1.0 && <0.4,+        random ==1.1.*,+        vector >=0.12.0.2 && <0.13  test-suite Coords-test     type: exitcode-stdio-1.0@@ -72,11 +79,11 @@     default-language: Haskell2010     ghc-options: -threaded -rtsopts -with-rtsopts=-N     build-depends:-        base <4.13,-        containers <0.7,-        hspec <2.7,+        base >=4.12.0.0 && <4.13,+        containers >=0.6.0.1 && <0.7,+        hspec >=2.6.1 && <2.7,         math-grads -any,-        random <1.2+        random ==1.1.*  test-suite Graph-test     type: exitcode-stdio-1.0@@ -85,9 +92,9 @@     default-language: Haskell2010     ghc-options: -threaded -rtsopts -with-rtsopts=-N     build-depends:-        base <4.13,-        containers <0.7,-        hspec <2.7,+        base >=4.12.0.0 && <4.13,+        containers >=0.6.0.1 && <0.7,+        hspec >=2.6.1 && <2.7,         math-grads -any  test-suite Isomorphism-test@@ -97,8 +104,21 @@     default-language: Haskell2010     ghc-options: -threaded -rtsopts -with-rtsopts=-N     build-depends:-        base <4.13,-        array <0.6,-        containers <0.7,-        hspec <2.7,+        base >=4.12.0.0 && <4.13,+        array >=0.5.3.0 && <0.6,+        containers >=0.6.0.1 && <0.7,+        hspec >=2.6.1 && <2.7,+        math-grads -any++test-suite SSSR-test+    type: exitcode-stdio-1.0+    main-is: SSSR.hs+    hs-source-dirs: test+    default-language: Haskell2010+    ghc-options: -threaded -rtsopts -with-rtsopts=-N+    build-depends:+        base >=4.12.0.0 && <4.13,+        array >=0.5.3.0 && <0.6,+        containers >=0.6.0.1 && <0.7,+        hspec >=2.6.1 && <2.7,         math-grads -any

@@ -2,6 +2,7 @@ -- module Math.Grads.Algo.Cycles   ( findCycles+  , findSimpleCycles   , findLocalCycles   , getCyclic   , isEdgeInCycle@@ -87,11 +88,13 @@                         else helperFilter (tail res) [head res]   where     -- TODO: We need to remove this filter.-    cycles = filter (\x -> length x < 21) (findLocalCycles' bonds)+    cycles = filter (\x -> length x < 21) (findSimpleCycles bonds)     res = filter (`filterBigCycles` cycles) cycles -findLocalCycles' :: Eq e => EdgeList e -> [EdgeList e]-findLocalCycles' bonds = concatMap (\(a, b, _) -> dfsAllPaths bonds a b) cycleBonds+-- | Finds all simple cycles in fused cycles system.+--+findSimpleCycles :: Eq e => EdgeList e -> [EdgeList e]+findSimpleCycles bonds = concatMap (\(a, b, _) -> dfsAllPaths bonds a b) cycleBonds   where     stBonds    = dfsSt bonds     cycleBonds = bonds \\ stBonds

@@ -6,7 +6,7 @@ -- module Math.Grads.Algo.Isomorphism   ( EComparator, VComparator-  , GComparable (..)+  , GComparable(..)   , VertexIndex   , getIso   , getMultiIso@@ -14,47 +14,17 @@   , isIsoSub   ) where -import           Control.Arrow           (second, (&&&), (***))-import qualified Data.Array              as A-import           Data.List               (delete, sortOn)-import           Data.Map                (Map)-import qualified Data.Map.Strict         as M-import           Data.Matrix             (Matrix (..), getElem, getRow, mapRow,-                                          matrix, multStd, ncols, nrows,-                                          setElem, transpose)-import           Data.Maybe              (isJust, listToMaybe)-import           Data.Tuple              (swap)-import qualified Data.Vector             as V-import           Math.Grads.GenericGraph (GenericGraph (..))-import           Math.Grads.Graph        (Graph, GraphEdge, changeIndsEdge,-                                          fromList, incidentIdx, toList, vCount,-                                          (!.))-import           Math.Grads.Utils        (nub)+import           Data.Map                          (Map)+import           Data.Maybe                        (isJust, listToMaybe) -type GenericGraphIso v e = GenericGraph Int e+import qualified Math.Grads.Algo.Isomorphism.RI    as RI+import           Math.Grads.Algo.Isomorphism.Types (EComparator,+                                                    GComparable (..),+                                                    VComparator, VertexIndex)+import           Math.Grads.GenericGraph           (GenericGraph (..))+import           Math.Grads.Graph                  (toList) --- | Type alias for 'Int'.----type VertexIndex = Int --- | Function that checks whether two vertices are identical.--- Due to properties related to index of vertex,--- like number of neighbors, we consider vertex indices instead of vertices.----type VComparator v1 v2 = VertexIndex -> VertexIndex -> Bool---- | Function that checks whether two edges are identical.--- Due to properties related to index of vertex,--- like belonging to a cycle, we consider GraphEdge (Int, Int, e) instead of e.----type EComparator e1 e2 = GraphEdge e1 -> GraphEdge e2 -> Bool---- | Type class for graphs that could be checked for isomorphism.----class (Graph g1, Graph g2) => GComparable g1 v1 e1 g2 v2 e2 where-  vComparator :: g1 v1 e1 -> g2 v2 e2 -> VComparator v1 v2-  eComparator :: g1 v1 e1 -> g2 v2 e2 -> EComparator e1 e2- -- | Checks whether two graphs are isomorphic. -- isIso :: (Ord v1, Ord v2, GComparable GenericGraph v1 e1 GenericGraph v2 e2, Eq e1, Eq e2)@@ -91,211 +61,4 @@             => GenericGraph v1 e1 -- ^ queryGraph             -> GenericGraph v2 e2 -- ^ targetGraph             -> [Map Int Int]-getMultiIso queryGraph' targetGraph' = matches-  where-    ((queryGraph, queryGraphWI), fromIsoToOldQ) = second inverseMap $ graphToGraphIso queryGraph'-    ((targetGraph, targetGraphWI), fromIsoToOldT) = second inverseMap $ graphToGraphIso targetGraph'--    vComp = vComparator queryGraphWI targetGraphWI-    eComp = eComparator queryGraphWI targetGraphWI--    isos = isoGraph vComp eComp queryGraph targetGraph-    matches = fmap (\x -> getMatchMap x fromIsoToOldQ fromIsoToOldT) isos--inverseMap :: Map Int Int -> Map Int Int-inverseMap = M.fromList . (swap <$>) . M.toList--getMatchMap :: Matrix Int -> Map Int Int -> Map Int Int -> Map Int Int-getMatchMap isoMatrix fromIsoToOldQ fromIsoToOldT = res-  where-    forMap = fmap (getMatchRow isoMatrix) [0 .. nrows isoMatrix - 1]-    res = M.fromList (fmap ((fromIsoToOldQ M.!) *** (fromIsoToOldT M.!)) forMap)--getMatchRow :: Matrix Int -> Int -> (Int, Int)-getMatchRow isoMatrix ind = (ind, helper 0)-  where-    row = getRow (ind + 1) isoMatrix--    helper :: Int -> Int-    helper counter = if row V.! counter == 1 then counter-                     else helper (counter + 1)--isoGraph :: (Eq e1, Eq e2) => VComparator v1 v2-                           -> EComparator e1 e2-                           -> GenericGraphIso v1 e1-                           -> GenericGraphIso v2 e2-                           -> [Matrix Int]-isoGraph vComp eComp queryGraph targetGraph = res-  where-    queryGraphEdges = (fst <$>) <$> gAdjacency queryGraph-    sizeOfQueryGraph = vCount queryGraph-    pMatrix = matrix sizeOfQueryGraph sizeOfQueryGraph (\(i, j) -> if i - 1 `elem` queryGraphEdges A.! (j - 1) then 1 else 0)--    targetGraphEdges = (fst <$>) <$> gAdjacency targetGraph-    sizeOfTargetGraph = vCount targetGraph-    gMatrix = matrix sizeOfTargetGraph sizeOfTargetGraph (\(i, j) -> if i - 1 `elem` targetGraphEdges A.! (j - 1) then 1 else 0)--    mMatrix = matrix sizeOfQueryGraph sizeOfTargetGraph (\(i, j) -> if fits vComp eComp queryGraph targetGraph i j then 1 else 0)--    currentRow = 0-    unusedColumns = [1 .. ncols mMatrix]--    res = recurse eComp queryGraph targetGraph unusedColumns currentRow gMatrix pMatrix mMatrix--fits :: (Eq e1, Eq e2) => VComparator v1 v2-                       -> EComparator e1 e2-                       -> GenericGraphIso v1 e1-                       -> GenericGraphIso v2 e2-                       -> Int-                       -> Int-                       -> Bool-fits vComp eComp queryGraph targetGraph i j = res-  where-    (vertex, edges) = (gIndex queryGraph A.! (i - 1), incidentIdx queryGraph $ i - 1)-    (vertex', edges') = (gIndex targetGraph A.! (j - 1), incidentIdx targetGraph $ j - 1)-    res = length edges <= length edges' && canBeSubset eComp edges edges' && vertex `vComp` vertex'--canBeSubset :: forall e1 e2. EComparator e1 e2 -> [GraphEdge e1] -> [GraphEdge e2] -> Bool-canBeSubset eComp query target = uniqueSeq maps-  where-    bondsInd = zip [0 ..] target-    maps = findMatches <$> query--    findMatches :: GraphEdge e1 -> [Int]-    findMatches thisEdge = fst <$> filter (\(_, otherEdge) -> eComp thisEdge otherEdge) bondsInd--uniqueSeq :: [[Int]] -> Bool-uniqueSeq maps = res-  where-    seqs = sequence maps--    res = any (\x -> length x == length (nub x)) seqs---- | Converts input graph into graph in which vertices with most amount of edges have lowest indices.----graphToGraphIso :: (Ord v) => GenericGraph v e-                           -> ((GenericGraphIso v e, GenericGraph v e), M.Map Int Int)-graphToGraphIso graph = res-  where-    (vertices, edges) = toList graph-    vArr = gIndex graph--    indsWithNCount = fmap (id &&& (length . (graph !.))) [0 .. length vertices - 1]-    sortedInds = fst <$> sortOn (\x -> - (snd x)) indsWithNCount-    changesMap = M.fromList (zip sortedInds [0 ..])--    sortedV = fmap (vArr A.!) sortedInds-    changedEdges = fmap (changeIndsEdge (changesMap M.!)) edges--    forGraphWI = (sortedV, changedEdges)-    forGraph = ([0 .. length sortedV - 1], changedEdges)--    res = ((fromList forGraph, fromList forGraphWI), changesMap)---- | Ullman's subgraph isomorphism algorithm itself.----recurse :: (Eq e1, Eq e2) => EComparator e1 e2-                          -> GenericGraphIso v1 e1-                          -> GenericGraphIso v2 e2-                          -> [Int]-                          -> Int-                          -> Matrix Int-                          -> Matrix Int-                          -> Matrix Int-                          -> [Matrix Int]-recurse eComp queryGraph targetGraph unusedColumns currentRow gMatrix pMatrix mMatrix = res-  where-    prunedM = prune eComp queryGraph targetGraph mMatrix currentRow--    recs = concatMap pruneNext unusedColumns--    res | hasEmptyRow mMatrix = []-        | currentRow == nrows mMatrix && isIsomorphism gMatrix pMatrix mMatrix = [mMatrix]-        | not (hasEmptyRow prunedM) = recs-        | otherwise = []--    pruneNext :: Int -> [Matrix Int]-    pruneNext x = recurse eComp queryGraph targetGraph newColumns newRow gMatrix pMatrix changedMatrix-      where-        newColumns = delete x unusedColumns-        newRow = currentRow + 1-        changedMatrix = changeRow prunedM newRow x--prune :: (Eq e1, Eq e2) => EComparator e1 e2-                        -> GenericGraphIso v1 e1-                        -> GenericGraphIso v2 e2-                        -> Matrix Int-                        -> Int-                        -> Matrix Int-prune eComp queryGraph targetGraph mMatrix currentRow | null indicesToChange = mMatrix-                                                      | hasEmptyRow mMatrix = mMatrix-                                                      | otherwise = res-  where-    numberOfMRows = nrows mMatrix-    numberOfMColumns = ncols mMatrix-    pairsOfindices = [(i, j) | i <- [1.. numberOfMRows], j <- [1.. numberOfMColumns], getElem i j mMatrix == 1]--    suitPair :: Int -> Int -> Bool-    suitPair = hasSuitableNeighbors eComp queryGraph targetGraph mMatrix--    indicesToChange = filter (not . uncurry suitPair) pairsOfindices-    changedMMatrix = foldl (flip (setElem 0)) mMatrix indicesToChange--    res = prune eComp queryGraph targetGraph changedMMatrix currentRow---- | Returns True if we can map all neighbors of query vertex to neighbors of target vertex in mMatrix.----hasSuitableNeighbors :: forall v1 v2 e1 e2. (Eq e1, Eq e2) => EComparator e1 e2-                                                           -> GenericGraphIso v1 e1-                                                           -> GenericGraphIso v2 e2-                                                           -> Matrix Int-                                                           -> Int-                                                           -> Int-                                                           -> Bool-hasSuitableNeighbors eComp queryGraph targetGraph mMatrix query target = doesSatisfy-  where-    iQ = query - 1-    iT = target - 1--    neighborsOfQ = (\(i, e) -> (iQ, i, e)) <$> queryGraph !. iQ-    neighborsOfT = (\(i, e) -> (iT, i, e)) <$> targetGraph !. iT--    hasProperNeighbor :: GraphEdge e1 -> Bool-    hasProperNeighbor edge = any (\edge' -> getProperElem edge edge' == 1 && eComp edge edge') neighborsOfT--    getProperElem :: GraphEdge e1 -> GraphEdge e2 -> Int-    getProperElem (_, b, _) (_, b', _) = getElem (b + 1) (b' + 1) mMatrix--    doesSatisfy = all hasProperNeighbor neighborsOfQ---- | Checks whether mMatrix encodes an isomorphism between pMatrix and gMatrix.----isIsomorphism :: Matrix Int -- ^ gMatrix-              -> Matrix Int -- ^ pMatrix-              -> Matrix Int -- ^ mMatrix-              -> Bool-isIsomorphism gMatrix pMatrix mMatrix = leqMatrices pMatrix check-  where-    check = multStd mMatrix (transpose (multStd mMatrix gMatrix))---- | Componentwise "less or equal" operation for matrices.----leqMatrices :: Matrix Int -> Matrix Int -> Bool-leqMatrices matrixA matrixB = nrows matrixA * ncols matrixA <= nrows matrixB * ncols matrixB && helper elems-  where-    numOfRows = nrows matrixA-    numOfColumns = ncols matrixB-    elems = [(i, j) | i <- [1..numOfRows], j <- [1..numOfColumns]]-    helper = foldr (\x -> (&&) (uncurry getElem x matrixA <= uncurry getElem x matrixB)) True---- | Replace all elements in row with 0 apart from chosen one.----changeRow :: Matrix Int -> Int -> Int -> Matrix Int-changeRow mMatrix row column = mapRow helper row mMatrix-  where helper column' a = if column' /= column then 0 else a--hasEmptyRow :: Matrix Int -> Bool-hasEmptyRow prunedMatrix = cond-  where-    numberOfRows = nrows prunedMatrix-    cond = any (\x -> all (== 0) (getRow x prunedMatrix)) [1 .. numberOfRows]+getMultiIso = RI.getMultiIso

@@ -0,0 +1,127 @@+{-# LANGUAGE FlexibleContexts    #-}+{-# LANGUAGE RankNTypes          #-}+{-# LANGUAGE ScopedTypeVariables #-}+{-# LANGUAGE ViewPatterns        #-}++module Math.Grads.Algo.Isomorphism.RI+  ( getMultiIso+  ) where++import           Data.Array                        (Array)+import qualified Data.Array                        as A+import           Data.Bimap                        (Bimap)+import qualified Data.Bimap                        as BM+import           Data.List                         (delete, find, maximumBy)+import           Data.Map                          (Map)+import           Data.Ord                          (comparing)++import           Math.Grads.Algo.Isomorphism.Types (GComparable (..))+import           Math.Grads.GenericGraph           (GenericGraph, gAdjacency,+                                                    getEdge)+++-- | RI isomorphism algorithm.+-- <https://bmcbioinformatics.biomedcentral.com/articles/10.1186/1471-2105-14-S7-S13>++type AdjArr  = Array Int [Int]+type VComp   = (Int -> Int -> Bool)+type EComp   = ((Int, Int) -> (Int, Int) -> Bool)+++getMultiIso :: Ord v1 => Eq e1 => Ord v2 => Eq e2+            => GComparable GenericGraph v1 e1 GenericGraph v2 e2+            => GenericGraph v1 e1 -> GenericGraph v2 e2+            -> [Map Int Int]+getMultiIso queryGraph@(toAdjArr -> queryAdjArr) targetGraph@(toAdjArr -> targetAdjArr) = isos+  where+    (order, parents) = buildMatchingOrder queryAdjArr+    isos             = findIsos queryAdjArr targetAdjArr vComp eComp order parents++    vComp = vComparator queryGraph targetGraph+    eComp = eIndComp+      where+        comp                   = eComparator queryGraph targetGraph+        eIndComp (a, b) (x, y) = comp (a, b, getEdge queryGraph a b) (x, y, getEdge targetGraph x y)++buildMatchingOrder :: AdjArr -> ([Int], [Maybe Int])+buildMatchingOrder graph = buildOrder [] [] $ A.indices graph+  where+    buildOrder :: [Int] -> [Maybe Int] -> [Int] -> ([Int], [Maybe Int])+    buildOrder visited parents unvisited+      | null unvisited = (visited, parents)+      | otherwise      = buildOrder nextVisited nextParents nextUnvisited+      where+        maxVertex       = maximumBy (comparing vertexRank) unvisited+        maxVertexParent = find ((maxVertex `elem`) . (graph !.)) visited++        nextVisited   = visited ++ [maxVertex]+        nextParents   = parents ++ [maxVertexParent]+        nextUnvisited = delete maxVertex unvisited++        vertexRank :: Int -> (Int, Int, Int)+        vertexRank ind+          | null visited = (degree,  0,       0)+          | otherwise    = (visRank, neiRank, unvisRank)+          where+            neis      = graph !. ind+            neisUnvis = filter (`notElem` visited) neis++            degree    = length neis+            visRank   = count neisWithInd visited+            neiRank   = count neisWithInd1Unvis visited+            unvisRank = count notNeiWithVis neisUnvis++            neisWithInd       = (ind `elem`) . (graph !.)+            neisWithInd1Unvis = any (`elem` neisUnvis) . (graph !.)+            notNeiWithVis     = all (`notElem` visited) . (graph !.)++findIsos :: AdjArr -> AdjArr -> VComp -> EComp+         -> [Int] -> [Maybe Int]+         -> [Map Int Int]+findIsos queryGraph targetGraph vComp eComp order parents =+  BM.toMap <$> goRI order parents BM.empty+  where+    goRI :: [Int] -> [Maybe Int] -> Bimap Int Int -> [Bimap Int Int]+    goRI []            _             match = [match]+    goRI _             []            match = [match]+    goRI (queryV : vs) (parent : ps) match = concatMap (goRI vs ps) matches+      where+        matches :: [Bimap Int Int]+        matches = (\targetV -> BM.insert queryV targetV match) <$> targetVs++        targetVs :: [Int]+        targetVs = filter isValidMatch $ maybe targetVsUnmatched findCandidates parent++        targetVsUnmatched :: [Int]+        targetVsUnmatched = filter (`BM.notMemberR` match) $ A.indices targetGraph++        findCandidates :: Int -> [Int]+        findCandidates = filter (`BM.notMemberR` match) . (targetGraph !.) . (match BM.!)++        isValidMatch :: Int -> Bool+        isValidMatch targetV =+          vComp queryV targetV &&+          length queryNeis <= length targetNeis &&+          all (`elem` targetNeisMatched) queryNeisMatchedProjection &&+          all (uncurry eComp) matchedEdges+          where+            queryNeis  = queryGraph  !. queryV+            targetNeis = targetGraph !. targetV++            queryNeisMatched           = filter (`BM.member` match) queryNeis+            queryNeisMatchedProjection = (match BM.!) <$> queryNeisMatched+            targetNeisMatched          = filter (`BM.memberR` match) targetNeis++            matchedEdges :: [((Int, Int), (Int, Int))]+            matchedEdges = matchedEdgePair <$> queryNeisMatched+              where+                matchedEdgePair queryNei = ((queryV, queryNei), (targetV, match BM.! queryNei))++toAdjArr :: GenericGraph v b -> Array Int [Int]+toAdjArr = fmap (fmap fst) . gAdjacency++(!.) :: AdjArr -> Int -> [Int]+(!.) = (A.!)++count :: (a -> Bool) -> [a] -> Int+count p = length . filter p

@@ -0,0 +1,32 @@+{-# LANGUAGE MultiParamTypeClasses #-}+module Math.Grads.Algo.Isomorphism.Types+  ( VertexIndex+  , VComparator+  , EComparator+  , GComparable(..)+  ) where++import           Math.Grads.Graph (Graph, GraphEdge)+++-- | Type alias for 'Int'.+--+type VertexIndex = Int++-- | Function that checks whether two vertices are identical.+-- Due to properties related to index of vertex,+-- like number of neighbors, we consider vertex indices instead of vertices.+--+type VComparator v1 v2 = VertexIndex -> VertexIndex -> Bool++-- | Function that checks whether two edges are identical.+-- Due to properties related to index of vertex,+-- like belonging to a cycle, we consider GraphEdge (Int, Int, e) instead of e.+--+type EComparator e1 e2 = GraphEdge e1 -> GraphEdge e2 -> Bool++-- | Type class for graphs that could be checked for isomorphism.+--+class (Graph g1, Graph g2) => GComparable g1 v1 e1 g2 v2 e2 where+  vComparator :: g1 v1 e1 -> g2 v2 e2 -> VComparator v1 v2+  eComparator :: g1 v1 e1 -> g2 v2 e2 -> EComparator e1 e2

@@ -0,0 +1,244 @@+{-# LANGUAGE FlexibleContexts      #-}+{-# LANGUAGE MultiParamTypeClasses #-}+{-# LANGUAGE ScopedTypeVariables   #-}++module Math.Grads.Algo.Isomorphism.Ullman+  ( getMultiIso+  ) where++import           Control.Arrow                     (second, (&&&), (***))+import qualified Data.Array                        as A+import           Data.List                         (delete, sortOn)+import           Data.Map                          (Map)+import qualified Data.Map.Strict                   as M+import           Data.Matrix                       (Matrix (..), getElem,+                                                    getRow, mapRow, matrix,+                                                    multStd, ncols, nrows,+                                                    setElem, transpose)+import           Data.Tuple                        (swap)+import qualified Data.Vector                       as V++import           Math.Grads.Algo.Isomorphism.Types (EComparator,+                                                    GComparable (..),+                                                    VComparator)+import           Math.Grads.GenericGraph           (GenericGraph (..))+import           Math.Grads.Graph                  (GraphEdge, changeIndsEdge,+                                                    fromList, incidentIdx,+                                                    toList, vCount, (!.))+import           Math.Grads.Utils                  (nub)+++type GenericGraphIso v e = GenericGraph Int e++getMultiIso :: (Ord v1, Ord v2, GComparable GenericGraph v1 e1 GenericGraph v2 e2, Eq e1, Eq e2)+            => GenericGraph v1 e1 -- ^ queryGraph+            -> GenericGraph v2 e2 -- ^ targetGraph+            -> [Map Int Int]+getMultiIso queryGraph' targetGraph' = matches+  where+    ((queryGraph, queryGraphWI), fromIsoToOldQ) = second inverseMap $ graphToGraphIso queryGraph'+    ((targetGraph, targetGraphWI), fromIsoToOldT) = second inverseMap $ graphToGraphIso targetGraph'++    vComp = vComparator queryGraphWI targetGraphWI+    eComp = eComparator queryGraphWI targetGraphWI++    isos = isoGraph vComp eComp queryGraph targetGraph+    matches = fmap (\x -> getMatchMap x fromIsoToOldQ fromIsoToOldT) isos++inverseMap :: Map Int Int -> Map Int Int+inverseMap = M.fromList . (swap <$>) . M.toList++getMatchMap :: Matrix Int -> Map Int Int -> Map Int Int -> Map Int Int+getMatchMap isoMatrix fromIsoToOldQ fromIsoToOldT = res+  where+    forMap = fmap (getMatchRow isoMatrix) [0 .. nrows isoMatrix - 1]+    res = M.fromList (fmap ((fromIsoToOldQ M.!) *** (fromIsoToOldT M.!)) forMap)++getMatchRow :: Matrix Int -> Int -> (Int, Int)+getMatchRow isoMatrix ind = (ind, helper 0)+  where+    row = getRow (ind + 1) isoMatrix++    helper :: Int -> Int+    helper counter = if row V.! counter == 1 then counter+                     else helper (counter + 1)++isoGraph :: (Eq e1, Eq e2) => VComparator v1 v2+                           -> EComparator e1 e2+                           -> GenericGraphIso v1 e1+                           -> GenericGraphIso v2 e2+                           -> [Matrix Int]+isoGraph vComp eComp queryGraph targetGraph = res+  where+    queryGraphEdges = (fst <$>) <$> gAdjacency queryGraph+    sizeOfQueryGraph = vCount queryGraph+    pMatrix = matrix sizeOfQueryGraph sizeOfQueryGraph (\(i, j) -> if i - 1 `elem` queryGraphEdges A.! (j - 1) then 1 else 0)++    targetGraphEdges = (fst <$>) <$> gAdjacency targetGraph+    sizeOfTargetGraph = vCount targetGraph+    gMatrix = matrix sizeOfTargetGraph sizeOfTargetGraph (\(i, j) -> if i - 1 `elem` targetGraphEdges A.! (j - 1) then 1 else 0)++    mMatrix = matrix sizeOfQueryGraph sizeOfTargetGraph (\(i, j) -> if fits vComp eComp queryGraph targetGraph i j then 1 else 0)++    currentRow = 0+    unusedColumns = [1 .. ncols mMatrix]++    res = recurse eComp queryGraph targetGraph unusedColumns currentRow gMatrix pMatrix mMatrix++fits :: (Eq e1, Eq e2) => VComparator v1 v2+                       -> EComparator e1 e2+                       -> GenericGraphIso v1 e1+                       -> GenericGraphIso v2 e2+                       -> Int+                       -> Int+                       -> Bool+fits vComp eComp queryGraph targetGraph i j = res+  where+    (vertex, edges) = (gIndex queryGraph A.! (i - 1), incidentIdx queryGraph $ i - 1)+    (vertex', edges') = (gIndex targetGraph A.! (j - 1), incidentIdx targetGraph $ j - 1)+    res = length edges <= length edges' && canBeSubset eComp edges edges' && vertex `vComp` vertex'++canBeSubset :: forall e1 e2. EComparator e1 e2 -> [GraphEdge e1] -> [GraphEdge e2] -> Bool+canBeSubset eComp query target = uniqueSeq maps+  where+    bondsInd = zip [0 ..] target+    maps = findMatches <$> query++    findMatches :: GraphEdge e1 -> [Int]+    findMatches thisEdge = fst <$> filter (\(_, otherEdge) -> eComp thisEdge otherEdge) bondsInd++uniqueSeq :: [[Int]] -> Bool+uniqueSeq maps = res+  where+    seqs = sequence maps++    res = any (\x -> length x == length (nub x)) seqs++-- | Converts input graph into graph in which vertices with most amount of edges have lowest indices.+--+graphToGraphIso :: (Ord v) => GenericGraph v e+                           -> ((GenericGraphIso v e, GenericGraph v e), M.Map Int Int)+graphToGraphIso graph = res+  where+    (vertices, edges) = toList graph+    vArr = gIndex graph++    indsWithNCount = fmap (id &&& (length . (graph !.))) [0 .. length vertices - 1]+    sortedInds = fst <$> sortOn (\x -> - (snd x)) indsWithNCount+    changesMap = M.fromList (zip sortedInds [0 ..])++    sortedV = fmap (vArr A.!) sortedInds+    changedEdges = fmap (changeIndsEdge (changesMap M.!)) edges++    forGraphWI = (sortedV, changedEdges)+    forGraph = ([0 .. length sortedV - 1], changedEdges)++    res = ((fromList forGraph, fromList forGraphWI), changesMap)++-- | Ullmann's subgraph isomorphism algorithm itself.+--+recurse :: (Eq e1, Eq e2) => EComparator e1 e2+                          -> GenericGraphIso v1 e1+                          -> GenericGraphIso v2 e2+                          -> [Int]+                          -> Int+                          -> Matrix Int+                          -> Matrix Int+                          -> Matrix Int+                          -> [Matrix Int]+recurse eComp queryGraph targetGraph unusedColumns currentRow gMatrix pMatrix mMatrix = res+  where+    prunedM = prune eComp queryGraph targetGraph mMatrix currentRow++    recs = concatMap pruneNext unusedColumns++    res | hasEmptyRow mMatrix = []+        | currentRow == nrows mMatrix && isIsomorphism gMatrix pMatrix mMatrix = [mMatrix]+        | not (hasEmptyRow prunedM) = recs+        | otherwise = []++    pruneNext :: Int -> [Matrix Int]+    pruneNext x = recurse eComp queryGraph targetGraph newColumns newRow gMatrix pMatrix changedMatrix+      where+        newColumns = delete x unusedColumns+        newRow = currentRow + 1+        changedMatrix = changeRow prunedM newRow x++prune :: (Eq e1, Eq e2) => EComparator e1 e2+                        -> GenericGraphIso v1 e1+                        -> GenericGraphIso v2 e2+                        -> Matrix Int+                        -> Int+                        -> Matrix Int+prune eComp queryGraph targetGraph mMatrix currentRow | null indicesToChange = mMatrix+                                                      | hasEmptyRow mMatrix = mMatrix+                                                      | otherwise = res+  where+    numberOfMRows = nrows mMatrix+    numberOfMColumns = ncols mMatrix+    pairsOfindices = [(i, j) | i <- [1.. numberOfMRows], j <- [1.. numberOfMColumns], getElem i j mMatrix == 1]++    suitPair :: Int -> Int -> Bool+    suitPair = hasSuitableNeighbors eComp queryGraph targetGraph mMatrix++    indicesToChange = filter (not . uncurry suitPair) pairsOfindices+    changedMMatrix = foldl (flip (setElem 0)) mMatrix indicesToChange++    res = prune eComp queryGraph targetGraph changedMMatrix currentRow++-- | Returns True if we can map all neighbors of query vertex to neighbors of target vertex in mMatrix.+--+hasSuitableNeighbors :: forall v1 v2 e1 e2. (Eq e1, Eq e2) => EComparator e1 e2+                                                           -> GenericGraphIso v1 e1+                                                           -> GenericGraphIso v2 e2+                                                           -> Matrix Int+                                                           -> Int+                                                           -> Int+                                                           -> Bool+hasSuitableNeighbors eComp queryGraph targetGraph mMatrix query target = doesSatisfy+  where+    iQ = query - 1+    iT = target - 1++    neighborsOfQ = (\(i, e) -> (iQ, i, e)) <$> queryGraph !. iQ+    neighborsOfT = (\(i, e) -> (iT, i, e)) <$> targetGraph !. iT++    hasProperNeighbor :: GraphEdge e1 -> Bool+    hasProperNeighbor edge = any (\edge' -> getProperElem edge edge' == 1 && eComp edge edge') neighborsOfT++    getProperElem :: GraphEdge e1 -> GraphEdge e2 -> Int+    getProperElem (_, b, _) (_, b', _) = getElem (b + 1) (b' + 1) mMatrix++    doesSatisfy = all hasProperNeighbor neighborsOfQ++-- | Checks whether mMatrix encodes an isomorphism between pMatrix and gMatrix.+--+isIsomorphism :: Matrix Int -- ^ gMatrix+              -> Matrix Int -- ^ pMatrix+              -> Matrix Int -- ^ mMatrix+              -> Bool+isIsomorphism gMatrix pMatrix mMatrix = leqMatrices pMatrix check+  where+    check = multStd mMatrix (transpose (multStd mMatrix gMatrix))++-- | Componentwise "less or equal" operation for matrices.+--+leqMatrices :: Matrix Int -> Matrix Int -> Bool+leqMatrices matrixA matrixB = nrows matrixA * ncols matrixA <= nrows matrixB * ncols matrixB && helper elems+  where+    numOfRows = nrows matrixA+    numOfColumns = ncols matrixB+    elems = [(i, j) | i <- [1..numOfRows], j <- [1..numOfColumns]]+    helper = foldr (\x -> (&&) (uncurry getElem x matrixA <= uncurry getElem x matrixB)) True++-- | Replace all elements in row with 0 apart from chosen one.+--+changeRow :: Matrix Int -> Int -> Int -> Matrix Int+changeRow mMatrix row column = mapRow helper row mMatrix+  where helper column' a = if column' /= column then 0 else a++hasEmptyRow :: Matrix Int -> Bool+hasEmptyRow prunedMatrix = cond+  where+    numberOfRows = nrows prunedMatrix+    cond = any (\x -> all (== 0) (getRow x prunedMatrix)) [1 .. numberOfRows]

@@ -0,0 +1,125 @@+{-# LANGUAGE MultiWayIf #-}+{-# LANGUAGE RankNTypes #-}++module Math.Grads.Algo.SSSR+  ( findSSSR+  ) where++import           Prelude                 hiding (map)++import           Control.Arrow           ((***))+import           Control.Lens            (over, _1, _2)+import           Data.Bimap              ((!>))+import           Data.List               (intersect, nub, sort)+import           Data.List.Index         (ifoldl)+import           Data.Map.Strict         (Map)+import qualified Data.Map.Strict         as M (empty, insert, member, (!))+import           Data.Matrix             (Matrix, matrix, unsafeGet, unsafeSet)+import qualified Data.Set                as S++import           Math.Grads.Algo.Cycles  (getCyclic)+import           Math.Grads.GenericGraph (GenericGraph, subgraphWithReindex)+import           Math.Grads.Graph        (EdgeList, toList)+++-- | RP-Path algorithm for searching the smallest set of smallest rings.+-- <https://www.ncbi.nlm.nih.gov/pubmed/19805142>++findSSSR :: Ord v => Ord e => GenericGraph v e -> [EdgeList e]+findSSSR graph = sssr+  where+    (reindex, cyclicGraph) = subgraphWithReindex graph . S.toList $ getCyclic graph+    g@(_, edges)           = toList cyclicGraph+    (n, m)                 = (length *** length) g+    maxSSSRs               = m - n + 1++    edgeIndex :: Map (Int, Int) Int+    edgeIndex = ifoldl insertEdge M.empty edges+      where+        insertEdge map ix (x, y, _) = M.insert (y, x) ix $ M.insert (x, y) ix map++    (pid, pid') = calculatePidMatrices n m edgeIndex+    sssrEdges   = takeSSSR n maxSSSRs pid pid'+    sssr        = fmap (backToOriginIndex . (edges !!)) <$> sssrEdges+      where+        backToOriginIndex = over _1 (reindex !>) . over _2 (reindex !>)++takeSSSR :: Int -> Int -> Matrix (Int, [[Int]]) -> Matrix [[Int]] -> [[Int]]+takeSSSR n maxSSSRs pid pid' = takeSSSR' [] [] 3 (1, 1)+  where+    takeSSSR' cycles edges len (i, j)+      | (i, j) > (n, n) || length cycles >= maxSSSRs =+          cycles+      | curLen /= len || length ij < 1 || length ij == 1 && null ij' =+          takeSSSR' cycles edges nextLen nextInd+      | otherwise =+          takeSSSR' nextCycles nextEdges nextLen nextInd+      where+        (ijDist, ij) = unsafeGet i j pid+        ij'          = unsafeGet i j pid'++        curLen  = ijDist * 2 + len `mod` 2++        newCycles = take (maxSSSRs - length cycles) $+          filter (`notElem` cycles) $+          filter (\c -> all (c `notContains`) cycles) $+          nub $+          cartesianProduct ij (if even len then ij else ij') notIntersects concatSort+          where+            notIntersects a b = null $ intersect a b+            concatSort a b    = sort $ a ++ b++            notContains a b       = length (a `intersect` b) < length b - 1++        nextCycles = cycles ++ newCycles+        nextEdges  = nub $ edges ++ concat newCycles++        (nextLen, nextInd)+          | (i, j) == (n, n) = (len + 1, (1,     1))+          | j == n           = (len,     (i + 1, 1))+          | otherwise        = (len,     (i,     j + 1))++calculatePidMatrices :: Int -> Int -> Map (Int, Int) Int+                     -> (Matrix (Int, [[Int]]), Matrix [[Int]])+calculatePidMatrices n m edgeIndex = calcPids initPid initPid' (1, 1, 1)+  where+    initPid :: Matrix(Int, [[Int]])+    initPid = matrix n n $ \(x, y) -> if+      | (x - 1, y - 1) `M.member` edgeIndex -> (1, [[edgeIndex M.! (x - 1, y - 1)]])+      | otherwise                           -> (m + 1, [])++    initPid' :: Matrix [[Int]]+    initPid' = matrix n n $ const []++    calcPids :: Matrix(Int, [[Int]]) -> Matrix [[Int]]+             -> (Int, Int, Int)+             -> (Matrix(Int, [[Int]]), Matrix [[Int]])+    calcPids pid pid' ind@(k, i, j)+      | ind > (n, n, n) = (pid, pid')+      | otherwise       = calcPids nextPid nextPid' nextInd+      where+        (ijDist, ij) = unsafeGet i j pid+        (ikDist, ik) = unsafeGet i k pid+        (kjDist, kj) = unsafeGet k j pid+        ikjDist      = ikDist + kjDist+        ikj          = cartesianProduct ik kj (/=) (++)+        ij'          = unsafeGet i j pid'++        nextPid+          | ijDist >  ikjDist = unsafeSet (ikjDist, ikj)       (i, j) pid+          | ijDist == ikjDist = unsafeSet (ijDist,  ij ++ ikj) (i, j) pid+          | otherwise         = pid++        nextPid'+          | ijDist == ikjDist + 1 = unsafeSet ij           (i, j) pid'+          | ijDist == ikjDist - 1 = unsafeSet (ij' ++ ikj) (i, j) pid'+          | ijDist >  ikjDist + 1 = unsafeSet []           (i, j) pid'+          | otherwise             = pid'++        nextInd+          | (i, j) == (n, n) = (k + 1, 1,     1)+          | j == n           = (k,     i + 1, i + 1)+          | otherwise        = (k,          i,j + 1)++cartesianProduct :: Eq a => [a] -> [a] -> (a -> a -> Bool) -> (a -> a -> b) -> [b]+cartesianProduct xs ys f g = [ g x y | x <- xs, y <- ys, f x y ]

@@ -18,7 +18,7 @@   , removeVertices   , safeAt   , safeIdx-  , subgraph+  , subgraph, subgraphWithReindex   , sumGraphs   , typeOfEdge   ) where@@ -28,6 +28,8 @@                                    genericParseJSON, genericToJSON) import           Data.Array       (Array) import qualified Data.Array       as A+import           Data.Bimap       (Bimap)+import qualified Data.Bimap       as BM import           Data.List        (find, groupBy, sortBy) import           Data.Map.Strict  (Map, mapKeys, member, (!)) import qualified Data.Map.Strict  as M@@ -36,6 +38,7 @@ import           GHC.Generics     (Generic) import           Math.Grads.Graph (Graph (..)) + -- | Generic undirected graph which stores elements of type v in its vertices (e.g. labels, atoms, states etc) -- and elements of type e in its edges (e.g. weights, bond types, functions over states etc). -- Note that loops and multiple edges between two vertices are allowed.@@ -131,7 +134,13 @@ -- Be careful with !. and ?. operators. -- subgraph :: Ord v => GenericGraph v e -> [Int] -> GenericGraph v e-subgraph graph toKeep = fromList (newVertices, newEdges)+subgraph graph = snd . subgraphWithReindex graph++-- | Get subgraph on given vertices and mapping from old `toKeep` indices to+-- new indices of resulting subgraph.+--+subgraphWithReindex :: Ord v => GenericGraph v e -> [Int] -> (Bimap Int Int, GenericGraph v e)+subgraphWithReindex graph toKeep = (vMap, fromList (newVertices, newEdges))   where     vSet :: S.Set Int     vSet = S.fromList toKeep@@ -142,10 +151,10 @@     (oldVertices, edges)  = filter eRemain <$> toList graph     (newVertices, oldIdx) = unzip . filter (\(_, ix) -> ix `S.member` vSet) $ zip oldVertices [0..] -    vMap :: Map Int Int-    vMap = M.fromList $ zip oldIdx [0 ..]+    vMap :: Bimap Int Int+    vMap = BM.fromList $ zip oldIdx [0 ..] -    newEdges = map (\(at, other, bond) -> (vMap ! at, vMap ! other, bond)) edges+    newEdges = map (\(at, other, bond) -> (vMap BM.! at, vMap BM.! other, bond)) edges  -- | Add given vertices to graph. --

@@ -1,16 +1,23 @@-{-# OPTIONS_GHC -fno-warn-orphans #-}+{-# OPTIONS_GHC -fno-warn-orphans  #-}+{-# LANGUAGE FlexibleContexts      #-} {-# LANGUAGE MultiParamTypeClasses #-}  module Main where -import qualified Data.Array                  as A-import           Data.Map.Strict             (Map)-import qualified Data.Map.Strict             as M-import           Math.Grads.Algo.Isomorphism (GComparable (..), isIsoSub)-import           Math.Grads.GenericGraph     (GenericGraph, gIndex)-import           Math.Grads.Graph            (fromList)+import           Control.Arrow                      ((***))+import qualified Data.Array                         as A+import           Data.List                          (sort)+import           Data.Map.Strict                    (Map)+import qualified Data.Map.Strict                    as M+import           Math.Grads.Graph                   (fromList, toList) import           Test.Hspec +import           Math.Grads.Algo.Isomorphism        (GComparable (..), isIsoSub)+import qualified Math.Grads.Algo.Isomorphism.RI     as RI+import qualified Math.Grads.Algo.Isomorphism.Ullman as UI+import           Math.Grads.GenericGraph            (GenericGraph, gIndex,+                                                     getEdge)+ instance GComparable GenericGraph Int Int GenericGraph Int Int where   vComparator g1 g2 ind1 ind2 = gIndex g1 A.! ind1 == gIndex g2 A.! ind2   eComparator _ _ (_, _, t) (_, _, t') = t == t'@@ -69,28 +76,57 @@     it "Path" $ do         graph <- fmap (M.! "only_path") testMap         graph `shouldSatisfy` isIsoSub pathGraph+        uiIsoEqRiIso pathGraph graph     it "Conjugated cycles" $ do         graph <- fmap (M.! "only_cycles") testMap         graph `shouldSatisfy` isIsoSub conjugatedCycles+        uiIsoEqRiIso conjugatedCycles graph     it "Connected cycles" $ do         graph <- fmap (M.! "simple_drawing") testMap         graph `shouldSatisfy` isIsoSub connectedCycles+        uiIsoEqRiIso connectedCycles graph     it "Conjugated cycles again" $ do         graph <- fmap (M.! "hard_drawing") testMap         graph `shouldSatisfy` isIsoSub conjugatedCycles+        uiIsoEqRiIso conjugatedCycles graph     it "Cycle and triangle" $ do         graph <- fmap (M.! "paths_through_conjugated_cycles") testMap         graph `shouldSatisfy` isIsoSub cycleAndTriangle+        uiIsoEqRiIso cycleAndTriangle graph     it "Big graph" $ do         graph <- fmap (M.! "takes_long_if_done_wrong") testMap         graph `shouldSatisfy` isIsoSub bigSubGraph+        uiIsoEqRiIso bigSubGraph graph     it "Triangle and triangle. No match" $ do         graph <- fmap (M.! "paths_through_conjugated_cycles") testMap         graph `shouldNotSatisfy` isIsoSub triangleAndTriangle+        uiIsoEqRiIso triangleAndTriangle graph     it "Cycle and triangle. No match" $ do         graph <- fmap (M.! "simple_drawing") testMap         graph `shouldNotSatisfy` isIsoSub cycleAndTriangle+        uiIsoEqRiIso cycleAndTriangle graph +uiIsoEqRiIso :: GenericGraph Int Int -> GenericGraph Int Int -> Expectation+uiIsoEqRiIso query target = do+    mapM_ (`shouldSatisfy` isValidIso query target) uiIsos+    mapM_ (`shouldSatisfy` isValidIso query target) riIsos+    length uiIsos `shouldBe` length riIsos+    toIsoList uiIsos `shouldBe` toIsoList riIsos+  where+    uiIsos = UI.getMultiIso query target+    riIsos = RI.getMultiIso query target++    toIsoList = sort . fmap (sort . M.toList)++isValidIso :: GenericGraph Int Int -> GenericGraph Int Int -> Map Int Int -> Bool+isValidIso query target iso = vsEq && esEq+  where+    (queryVs, queryEs) = toList query+    targetVs           = fst $ toList target++    vsEq = all (uncurry (==) . ((queryVs !!) *** (targetVs !!))) $ M.toList iso+    esEq = all (\(v1, v2, t) -> t == getEdge target (iso M.! v1) (iso M.! v2)) queryEs++ main :: IO ()-main = hspec $ do-  testIsIsoSub+main = hspec testIsIsoSub

@@ -0,0 +1,92 @@+{-# LANGUAGE ViewPatterns #-}+module Main where++import           Data.List               (group, nub, sort)+import           Data.Map                (Map)+import qualified Data.Map                as M++import           Test.Hspec              (Expectation, Spec, describe, hspec,+                                          it, shouldBe, shouldMatchList,+                                          shouldSatisfy)++import           Math.Grads.Algo.SSSR    (findSSSR)+import           Math.Grads.GenericGraph (GenericGraph)+import           Math.Grads.Graph        (GraphEdge, fromList)+++main :: IO ()+main = hspec $ describe "SSSR" $ do+    graphSpec "graphA" graphA resultA+    graphSpec "graphB" graphB resultB+    graphSpec "graphC" graphC resultC+    graphSpec "graphD" graphD resultD+++graphSpec :: String -> GenericGraph Int Int -> Map Int Int -> Spec+graphSpec name (sort . fmap sort . findSSSR -> sssr) cycleMap = it name $ do+    sssr `shouldMatchList` nub sssr+    mapM_ checkSimpleCycle sssr+    sssr `shouldSatisfy` (==) (sum . M.elems $ cycleMap) . length+    mapM_ (\(len, count) -> numCyclesOfLen len `shouldBe` count) $ M.toList cycleMap+  where+    checkSimpleCycle :: [GraphEdge Int] -> Expectation+    checkSimpleCycle = mapM_ (`shouldSatisfy` (==) 2 . length) . group . sort . concatMap (\(x, y, _) -> [x, y])++    numCyclesOfLen :: Int -> Int+    numCyclesOfLen n = length . filter ((==) n . length) $ sssr+++-- | This 4 graphs comes from figure 5 page 5 of the original paper.+-- <https://www.ncbi.nlm.nih.gov/pubmed/19805142>+--++-- | Note that this graph image wrong in the paper - it has 3 cycles of length 4.+-- We fix it with 3 additional edges to match with results from table 4 page 5 of the paper.+--+graphA :: GenericGraph Int Int+graphA = fromList (+  [0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20],++  [(0,1,1),(0,5,1),(0,10,1),(0,13,1),(0,20,1),(1,2,1),(1,3,1),(1,4,1),(1,5,1),(1,20,1),(2,3,1),+  (2,8,1),(2,10,1),(2,20,1),(3,4,1),(3,8,1),(3,14,1),(3,15,1),(4,5,1),(4,15,1),(4,16,1),(4,19,1),+  (5,6,1),(5,13,1),(5,19,1),(6,7,1),(6,12,1),(6,13,1),(6,18,1),(6,19,1),(7,8,1),(7,9,1),(7,11,1),+  (7,12,1),(7,14,1),(7,17,1),(7,18,1),(8,9,1),(8,10,1),(8,14,1),(9,10,1),(9,11,1),(10,11,1),+  (10,12,1),(10,13,1),(11,12,1),(12,13,1),(14,15,1),(14,17,1),(15,16,1),(15,17,1),(16,17,1),+  (16,18,1),(16,19,1),(17,18,1),(18,19,1)])++resultA :: Map Int Int+resultA = M.fromList [(3, 36)]+++graphB :: GenericGraph Int Int+graphB = fromList (+  [0,1,2,3,4,5,6,7,8,9,10,11,12],++  [(0,1,1),(0,10,1),(1,2,1),(1,8,1),(1,12,1),(2,3,1),(2,7,1),(3,4,1),(3,12,1),(4,5,1),(4,11,1),+  (5,6,1),(5,10,1),(6,7,1),(6,9,1),(7,8,1),(8,9,1),(9,10,1),(10,11,1),(11,12,1)])++resultB :: Map Int Int+resultB = M.fromList [(4, 6), (5, 2)]+++graphC :: GenericGraph Int Int+graphC = fromList (+  [0,1,2,3,4,5,6,7,8,9,10,11,12],++  [(0,1,1),(0,12,1),(1,2,1),(1,8,1),(1,9,1),(1,10,1),(2,3,1),(2,4,1),(2,9,1),(3,4,1),(3,5,1),+  (3,9,1),(4,5,1),(4,12,1),(5,6,1),(5,12,1),(6,7,1),(6,11,1),(6,12,1),(7,8,1),(7,10,1),(7,11,1),+  (8,9,1),(8,10,1),(10,11,1),(11,12,1)])++resultC :: Map Int Int+resultC = M.fromList [(3, 12), (5, 2)]+++graphD :: GenericGraph Int Int+graphD = fromList (+  [0,1,2,3,4,5,6,7,8,9],++  [(0,1,1),(0,3,1),(1,2,1),(1,9,1),(2,3,1),(2,9,1),(3,4,1),(3,5,1),(4,5,1),(4,6,1),(5,6,1),(6,7,1),+  (6,8,1),(7,8,1),(7,9,1),(8,9,1)])++resultD :: Map Int Int+resultD = M.fromList [(3, 5), (4, 1), (6, 1)]