Wykonanie Floyd-Warshall w Haskell - Naprawianie wycieku przestrzeni

Question

Chciałem napisać wydajną implementację algorytmu najkrótszej ścieżki Floyd-Warshall dla wszystkich par w Haskell za pomocą Vector s, aby uzyskać dobrą wydajność.

Implementacja jest dość prosta, ale zamiast trójwymiarowego | V | × | V | × | V | matrix, używany jest wektor dwuwymiarowy, ponieważ tylko czytamy poprzednią wartość k.

Tak więc algorytm jest tak naprawdę serią kroków, w których wektor 2D jest przekazywany i generowany jest nowy wektor 2D. Ostateczny wektor 2D zawiera najkrótsze ścieżki między wszystkimi węzłami (i, j).

Moja intuicja podpowiadała mi, że ważne byłoby, aby upewnić się, że poprzedni 2D vector oceniano przed każdym kroku, więc użyłem BangPatterns na prev argument funkcji fw i ścisłym foldl':

{-# Language BangPatterns #-} 

import   Control.DeepSeq 
import   Control.Monad  (forM_) 
import   Data.List   (foldl') 
import qualified Data.Map.Strict  as M 
import   Data.Vector   (Vector, (!), (//)) 
import qualified Data.Vector   as V 
import qualified Data.Vector.Mutable as V hiding (length, replicate, take) 

type Graph = Vector (M.Map Int Double) 
type TwoDVector = Vector (Vector Double) 

infinity :: Double 
infinity = 1/0 

-- calculate shortest path between all pairs in the given graph, if there are 
-- negative cycles, return Nothing 
allPairsShortestPaths :: Graph -> Int -> Maybe TwoDVector 
allPairsShortestPaths g v = 
    let initial = fw g v V.empty 0 
     results = foldl' (fw g v) initial [1..v] 
    in if negCycle results 
     then Nothing 
     else Just results 
    where -- check for negative elements along the diagonal 
     negCycle a = any not $ map (\i -> a ! i ! i >= 0) [0..(V.length a-1)] 

-- one step of the Floyd-Warshall algorithm 
fw :: Graph -> Int -> TwoDVector -> Int -> TwoDVector 
fw g v !prev k = V.create $ do           -- ← bang 
    curr <- V.new v 
    forM_ [0..(v-1)] $ \i -> 
    V.write curr i $ V.create $ do 
     ivec <- V.new v 
     forM_ [0..(v-1)] $ \j -> do 
     let d = distance g prev i j k 
     V.write ivec j d 
     return ivec 
    return curr 

distance :: Graph -> TwoDVector -> Int -> Int -> Int -> Double 
distance g _ i j 0 -- base case; 0 if same vertex, edge weight if neighbours 
    | i == j = 0.0 
    | otherwise = M.findWithDefault infinity j (g ! i) 
distance _ a i j k = let c1 = a ! i ! j 
         c2 = (a ! i ! (k-1))+(a ! (k-1) ! j) 
         in min c1 c2 

Jednak podczas uruchamiania tego programu z grafem o 1000 węzłach z 47978 krawędziami, rzeczy nie wyglądają dobrze. Wykorzystanie pamięci jest bardzo wysokie, a uruchomienie programu trwa zbyt długo. Program został skompilowany z ghc -O2.

I przebudowany program do profilowania i ogranicza liczbę powtórzeń do 50:

results = foldl' (fw g v) initial [1..50] 

Potem prowadził program z +RTS -p -hc i +RTS -p -hd:

To jest ... interesujące, ale myślę, że to pokazuje, że jest zgodne umiłowanie ton uderzeń. Niedobrze.

Ok, więc po kilku ujęć w ciemności, dodałem deepseq w fw aby upewnić prevnaprawdę jest evaluted:

let d = prev `deepseq` distance g prev i j k 

Teraz wszystko wygląda lepiej, a ja mogę faktycznie uruchomić program do zakończenia ze stałym użyciem pamięci. To oczywiste, że huk na argumencie prev nie wystarczył.

Dla porównania z poprzednich wykresów, tutaj jest wykorzystanie pamięci przez 50 powtórzeń po dodaniu deepseq:

Ok, więc rzeczy są lepsze, ale nadal mam pewne pytania:

1. Czy to właściwe rozwiązanie dla tego wycieku przestrzeni kosmicznej? Nie mam racji, czując, że wstawianie deepseq jest trochę brzydkie?
2. Czy korzystanie z Vector s tutaj idiomatyczne/poprawne? Buduję zupełnie nowy wektor dla każdej iteracji i mam nadzieję, że śmieciarz usunie stare s. Vector.
3. Czy są inne rzeczy, które mogę zrobić, aby przyspieszyć ten proces?

Dla odniesienia, tutaj jest graph.txt: http://sebsauvage.net/paste/?45147f7caf8c5f29#7tiCiPovPHWRm1XNvrSb/zNl3ujF3xB3yehrxhEdVWw=

Oto main:

main = do 
    ls <- fmap lines $ readFile "graph.txt" 
    let numVerts = head . map read . words . head $ ls 
    let edges = map (map read . words) (tail ls) 
    let g = V.create $ do 
     g' <- V.new numVerts 
     forM_ [0..(numVerts-1)] (\idx -> V.write g' idx M.empty) 
     forM_ edges $ \[f,t,w] -> do 
      -- subtract one from vertex IDs so we can index directly 
      curr <- V.read g' (f-1) 
      V.write g' (f-1) $ M.insert (t-1) (fromIntegral w) curr 
     return g' 
    let a = allPairsShortestPaths g numVerts 
    case a of 
    Nothing -> putStrLn "Negative cycle detected." 
    Just a' -> do 
     putStrLn $ "The shortest, shortest path has length " 
       ++ show ((V.minimum . V.map V.minimum) a') 

Answer 1

Po pierwsze, pewne ogólne porządki Kod:

W funkcji fw, jawnie przeznaczyć i wypełnić zmienne wektory. Istnieje jednak gotowa funkcja do tego celu, a mianowicie generate. fw można zatem zapisać jako

V.generate v (\i -> V.generate v (\j -> distance g prev i j k)) 

Podobnie kod generacji wykres można zastąpić replicate i accum:

let parsedEdges = map (\[f,t,w] -> (f - 1, (t - 1, fromIntegral w))) edges 
let g = V.accum (flip (uncurry M.insert)) (V.replicate numVerts M.empty) parsedEdges 

Zauważ, że to całkowicie usuwa wszystkie potrzeby mutacji, bez utraty wydajności.

teraz do rzeczywistych pytań:

1. Z mojego doświadczenia wynika, deepseq jest bardzo przydatna, ale tylko jak szybko naprawić nieszczelności kosmicznych, takich jak ten. Podstawowym problemem nie jest to, że musisz wymusić wyniki po ich wyprodukowaniu. Zamiast tego użycie deepseq sugeruje, że powinieneś budować strukturę w bardziej ścisły sposób. W rzeczywistości, jeśli dodać wzór huk w kodzie tworzenia wektora tak:

   let !d = distance g prev i j k 
   

   Wtedy problem zostanie rozwiązany bez deepseq. Zauważ, że to nie działa z kodem generate, ponieważ z jakiegoś powodu (mogę utworzyć żądanie funkcji), vector nie zapewnia ścisłych funkcji dla wektorów pudełkowych. Jednak gdy dojdę do rozpakowanych wektorów w odpowiedzi na pytanie 3, które są surowe, oba podejścia działają bez adnotacji ścisłości.

2. O ile mi wiadomo, wzorzec wielokrotnego generowania nowych wektorów jest idiomatyczny. Jedyną rzeczą, która nie jest idiomatyczna, jest użycie zmienności - z wyjątkiem przypadków, gdy są one bezwzględnie konieczne, zmienne wektory są generalnie odradzane.

3. Istnieje kilka rzeczy do zrobienia:

   • Najprościej można zastąpić Map Int z IntMap. Ponieważ to nie jest powolne działanie funkcji, nie ma to większego znaczenia, ale może być znacznie szybsze przy dużych obciążeniach.

   • Możesz przełączyć się do korzystania z niezakodowanych wektorów. Chociaż wektor zewnętrzny musi pozostać w pudełku, ponieważ wektory wektorów nie mogą być rozpakowane, może to być wektor wewnętrzny. To rozwiązuje również twój problem związany ze ścisłością - ponieważ unboxowane wektory są surowe w swoich elementach, nie dostajesz przecieku przestrzeni. Zauważ, że na moim komputerze poprawia to wydajność z 4,1 sekundy do 1,3 sekundy, więc rozpakowywanie jest bardzo pomocne.

   • Możesz spłaszczyć wektor w jeden i użyć mnożenia i dzielenia, aby przełączać się między dwuwymiarowymi wskaźnikami i jednym wymiarowym wskaźnikiem. Nie polecam tego, ponieważ jest trochę zaangażowany, dość brzydki, a dzięki podziałowi faktycznie spowalnia kod na moim komputerze.

   • Można użyć repa. Ma to ogromną zaletę polegającą na automatycznym równolegywaniu kodu. Zauważ, że od repa spłaszczam jego tablice i najwyraźniej nie pozbywam się poprawnie podziałów potrzebnych do ładnego wypełnienia (można to zrobić z zagnieżdżonymi pętlami, ale myślę, że używa on jednej pętli i podziału), ma taką samą wydajność kara, jak wspomniałem powyżej, przynosząc czas wykonywania z 1,3 sekundy do 1,8. Jeśli jednak włączysz równoległość i użyjesz maszyny wielordzeniowej, zaczniesz dostrzegać pewne korzyści. Niestety, obecny przypadek testowy jest zbyt mały, by dostrzec wiele korzyści, więc na mojej 6 rdzeniowej maszynie widzę, jak spada z powrotem do 1,2 sekundy. Jeśli podniosę rozmiar z powrotem do [1..v] zamiast [1..50], paralelizm przenosi go z 32 do 13. Przypuszczalnie, jeśli nadasz temu programowi większe wejście, możesz zobaczyć więcej korzyści.

     Jeśli jesteś zainteresowany, zamieściłem moją repa -ified version here.

   • EDYCJA: Użyj -fllvm. Testowanie na moim komputerze, używając repa, uzyskuję 14,7 sekund bez równoległości, co jest prawie tak dobre, jak bez -fllvm i równolegle. Ogólnie rzecz biorąc, LLVM może po prostu obsługiwać kod oparty na tablicy, jak to bardzo dobrze.