Zapamiętywanie w Haskell?

136

Wszelkie wskazówki, jak skutecznie rozwiązać następującą funkcję w Haskellu dla dużych liczb (n > 108)

f(n) = max(n, f(n/2) + f(n/3) + f(n/4))

Widziałem przykłady zapamiętywania w Haskellu w celu rozwiązania liczb Fibonacciego, które obejmowały (leniwie) obliczanie wszystkich liczb Fibonacciego do wymaganego n. Ale w tym przypadku dla danego n wystarczy obliczyć bardzo niewiele wyników pośrednich.

Dzięki

Angel de Vicente
źródło
110
Tylko w tym sensie, że to jakaś praca, którą wykonuję w domu :-)
Angel de Vicente

Odpowiedzi:

256

Możemy to zrobić bardzo efektywnie, tworząc strukturę, którą możemy indeksować w czasie nieliniowym.

Ale najpierw,

{-# LANGUAGE BangPatterns #-}

import Data.Function (fix)

Zdefiniujmy f, ale niech używa „otwartej rekursji” zamiast wywoływać się bezpośrednio.

f :: (Int -> Int) -> Int -> Int
f mf 0 = 0
f mf n = max n $ mf (n `div` 2) +
                 mf (n `div` 3) +
                 mf (n `div` 4)

Możesz uzyskać niezapamiętać f, używającfix f

Pozwoli ci to sprawdzić, fczy to, co masz na myśli, dla małych wartości f, dzwoniąc, na przykład:fix f 123 = 144

Moglibyśmy to zapamiętać, definiując:

f_list :: [Int]
f_list = map (f faster_f) [0..]

faster_f :: Int -> Int
faster_f n = f_list !! n

To działa dość dobrze i zastępuje to, co miało zająć O (n ^ 3) czas, czymś, co zapamiętuje wyniki pośrednie.

Ale wciąż potrzeba czasu liniowego, aby indeksować, aby znaleźć zapamiętaną odpowiedź mf. Oznacza to, że wyniki takie jak:

*Main Data.List> faster_f 123801
248604

są do zaakceptowania, ale wynik nie skaluje się dużo lepiej. Możemy zrobić lepiej!

Najpierw zdefiniujmy nieskończone drzewo:

data Tree a = Tree (Tree a) a (Tree a)
instance Functor Tree where
    fmap f (Tree l m r) = Tree (fmap f l) (f m) (fmap f r)

Następnie zdefiniujemy sposób indeksowania do niego, abyśmy mogli zamiast tego znaleźć węzeł z indeksem nw czasie O (log n) :

index :: Tree a -> Int -> a
index (Tree _ m _) 0 = m
index (Tree l _ r) n = case (n - 1) `divMod` 2 of
    (q,0) -> index l q
    (q,1) -> index r q

... i możemy znaleźć drzewo pełne liczb naturalnych dla wygody, więc nie będziemy musieli bawić się tymi indeksami:

nats :: Tree Int
nats = go 0 1
    where
        go !n !s = Tree (go l s') n (go r s')
            where
                l = n + s
                r = l + s
                s' = s * 2

Ponieważ możemy indeksować, możesz po prostu przekonwertować drzewo na listę:

toList :: Tree a -> [a]
toList as = map (index as) [0..]

Możesz sprawdzić dotychczasową pracę, weryfikując, która toList natsdaje[0..]

Teraz,

f_tree :: Tree Int
f_tree = fmap (f fastest_f) nats

fastest_f :: Int -> Int
fastest_f = index f_tree

działa tak samo, jak w przypadku powyższej listy, ale zamiast potrzebować czasu liniowego na znalezienie każdego węzła, można go ścigać w czasie logarytmicznym.

Wynik jest znacznie szybszy:

*Main> fastest_f 12380192300
67652175206

*Main> fastest_f 12793129379123
120695231674999

W rzeczywistości jest to o wiele szybsze, że możesz przejść przez powyższe i zastąpić Intje Integerpowyższymi i niemal natychmiast uzyskać absurdalnie duże odpowiedzi

*Main> fastest_f' 1230891823091823018203123
93721573993600178112200489

*Main> fastest_f' 12308918230918230182031231231293810923
11097012733777002208302545289166620866358
Edward KMETT
źródło
3
Wypróbowałem ten kod i, co ciekawe, f_faster wydawał się wolniejszy niż f. Myślę, że te odniesienia do listy naprawdę spowolniły sprawę. Definicja nats i indeksu wydawała mi się dość tajemnicza, więc dodałem własną odpowiedź, która może wyjaśnić sprawę.
Pitarou
5
Nieskończona lista musi zajmować się połączoną listą o długości 111111111 elementów. Przypadek drzewa dotyczy logu n * liczby osiągniętych węzłów.
Edward KMETT,
2
tj. wersja listy musi tworzyć thunks dla wszystkich węzłów na liście, podczas gdy wersja drzewka unika tworzenia ich wielu.
Tom Ellis,
7
Wiem, że to dość stary post, ale nie powinien f_treebyć definiowany w whereklauzuli, aby uniknąć zapisywania niepotrzebnych ścieżek w drzewie między wywołaniami?
dfeuer
17
Powodem umieszczania go w CAF było to, że można było zapamiętywać połączenia. Gdybym miał kosztowną rozmowę, którą zapisywałem, prawdopodobnie zostawiłbym go w CAF, stąd technika pokazana tutaj. W rzeczywistej aplikacji istnieje oczywiście kompromis między korzyściami a kosztami trwałego zapamiętywania. Chociaż, biorąc pod uwagę pytanie, jak osiągnąć zapamiętywanie, myślę, że mylące byłoby odpowiadanie techniką, która celowo unika zapamiętywania w wywołaniach, a jeśli nic innego, ten komentarz tutaj wskaże ludziom fakt, że istnieją subtelności. ;)
Edward KMETT
17

Odpowiedź Edwarda jest tak wspaniałą perełką, że ją zduplikowałem i dostarczyłem implementacje memoListi memoTreekombinatory, które zapamiętują funkcję w formie otwarto-rekurencyjnej.

{-# LANGUAGE BangPatterns #-}

import Data.Function (fix)

f :: (Integer -> Integer) -> Integer -> Integer
f mf 0 = 0
f mf n = max n $ mf (div n 2) +
                 mf (div n 3) +
                 mf (div n 4)


-- Memoizing using a list

-- The memoizing functionality depends on this being in eta reduced form!
memoList :: ((Integer -> Integer) -> Integer -> Integer) -> Integer -> Integer
memoList f = memoList_f
  where memoList_f = (memo !!) . fromInteger
        memo = map (f memoList_f) [0..]

faster_f :: Integer -> Integer
faster_f = memoList f


-- Memoizing using a tree

data Tree a = Tree (Tree a) a (Tree a)
instance Functor Tree where
    fmap f (Tree l m r) = Tree (fmap f l) (f m) (fmap f r)

index :: Tree a -> Integer -> a
index (Tree _ m _) 0 = m
index (Tree l _ r) n = case (n - 1) `divMod` 2 of
    (q,0) -> index l q
    (q,1) -> index r q

nats :: Tree Integer
nats = go 0 1
    where
        go !n !s = Tree (go l s') n (go r s')
            where
                l = n + s
                r = l + s
                s' = s * 2

toList :: Tree a -> [a]
toList as = map (index as) [0..]

-- The memoizing functionality depends on this being in eta reduced form!
memoTree :: ((Integer -> Integer) -> Integer -> Integer) -> Integer -> Integer
memoTree f = memoTree_f
  where memoTree_f = index memo
        memo = fmap (f memoTree_f) nats

fastest_f :: Integer -> Integer
fastest_f = memoTree f
Tom Ellis
źródło
12

Nie jest to najbardziej efektywny sposób, ale zapamiętuje:

f = 0 : [ g n | n <- [1..] ]
    where g n = max n $ f!!(n `div` 2) + f!!(n `div` 3) + f!!(n `div` 4)

przy żądaniu f !! 144sprawdza się, czy f !! 143istnieje, ale jego dokładna wartość nie jest obliczana. Nadal jest to nieznany wynik obliczeń. Jedynymi dokładnymi obliczonymi wartościami są te potrzebne.

Tak więc początkowo, o ile zostało obliczone, program nic nie wie.

f = .... 

Kiedy wysyłamy żądanie f !! 12, zaczyna ono dopasowywać wzorce:

f = 0 : g 1 : g 2 : g 3 : g 4 : g 5 : g 6 : g 7 : g 8 : g 9 : g 10 : g 11 : g 12 : ...

Teraz zaczyna obliczać

f !! 12 = g 12 = max 12 $ f!!6 + f!!4 + f!!3

To rekurencyjnie tworzy kolejne żądanie na f, więc obliczamy

f !! 6 = g 6 = max 6 $ f !! 3 + f !! 2 + f !! 1
f !! 3 = g 3 = max 3 $ f !! 1 + f !! 1 + f !! 0
f !! 1 = g 1 = max 1 $ f !! 0 + f !! 0 + f !! 0
f !! 0 = 0

Teraz możemy trochę cofnąć

f !! 1 = g 1 = max 1 $ 0 + 0 + 0 = 1

Co oznacza, że ​​program teraz wie:

f = 0 : 1 : g 2 : g 3 : g 4 : g 5 : g 6 : g 7 : g 8 : g 9 : g 10 : g 11 : g 12 : ...

Nadal sączy się:

f !! 3 = g 3 = max 3 $ 1 + 1 + 0 = 3

Co oznacza, że ​​program teraz wie:

f = 0 : 1 : g 2 : 3 : g 4 : g 5 : g 6 : g 7 : g 8 : g 9 : g 10 : g 11 : g 12 : ...

Teraz kontynuujemy obliczenia f!!6:

f !! 6 = g 6 = max 6 $ 3 + f !! 2 + 1
f !! 2 = g 2 = max 2 $ f !! 1 + f !! 0 + f !! 0 = max 2 $ 1 + 0 + 0 = 2
f !! 6 = g 6 = max 6 $ 3 + 2 + 1 = 6

Co oznacza, że ​​program teraz wie:

f = 0 : 1 : 2 : 3 : g 4 : g 5 : 6 : g 7 : g 8 : g 9 : g 10 : g 11 : g 12 : ...

Teraz kontynuujemy obliczenia f!!12:

f !! 12 = g 12 = max 12 $ 6 + f!!4 + 3
f !! 4 = g 4 = max 4 $ f !! 2 + f !! 1 + f !! 1 = max 4 $ 2 + 1 + 1 = 4
f !! 12 = g 12 = max 12 $ 6 + 4 + 3 = 13

Co oznacza, że ​​program teraz wie:

f = 0 : 1 : 2 : 3 : 4 : g 5 : 6 : g 7 : g 8 : g 9 : g 10 : g 11 : 13 : ...

Więc obliczenia są wykonywane dość leniwie. Program wie, że f !! 8istnieje jakaś wartość for , że jest równa g 8, ale nie ma pojęcia, co nią g 8jest.

rampion
źródło
Dziękuję za to. Jak stworzyłbyś i wykorzystał dwuwymiarową przestrzeń rozwiązań? Czy to byłaby lista list? ig n m = (something with) f!!a!!b
vikingsteve
1
Jasne, że możesz. Dla prawdziwego rozwiązania, choć i pewnie korzystać z biblioteki memoization, jak memocombinators
zerwa
Niestety, to O (n ^ 2).
Qumeric
8

To jest dodatek do doskonałej odpowiedzi Edwarda Kmetta.

Kiedy wypróbowałem jego kod, definicje natsi indexwydawały się dość tajemnicze, więc piszę alternatywną wersję, która była dla mnie łatwiejsza do zrozumienia.

Definiuję indexi natsw kategoriach index'i nats'.

index' t njest zdefiniowany w zakresie [1..]. (Przypomnij sobie, że index tjest to zdefiniowane w zakresie [0..].) To działa przeszukuje drzewo traktując je njako ciąg bitów i odczytując je w odwrotnej kolejności. Jeśli bit jest 1, bierze prawą gałąź. Jeśli bit jest 0, bierze lewą gałąź. Zatrzymuje się, gdy osiągnie ostatni bit (który musi być a 1).

index' (Tree l m r) 1 = m
index' (Tree l m r) n = case n `divMod` 2 of
                          (n', 0) -> index' l n'
                          (n', 1) -> index' r n'

Tak jak natsjest zdefiniowane dla, indextak index nats n == njest zawsze prawdziwe, nats'jest zdefiniowane dla index'.

nats' = Tree l 1 r
  where
    l = fmap (\n -> n*2)     nats'
    r = fmap (\n -> n*2 + 1) nats'
    nats' = Tree l 1 r

Teraz natsi indexsą po prostu nats'i index'ale z wartościami przesuniętymi o 1:

index t n = index' t (n+1)
nats = fmap (\n -> n-1) nats'
Pitarou
źródło
Dzięki. Zapamiętuję funkcję wielowymiarową i to naprawdę pomogło mi dowiedzieć się, co naprawdę robią indeks i nats.
Kittsil
8

Jak stwierdzono w odpowiedzi Edwarda Kmetta, aby przyspieszyć działanie, trzeba buforować kosztowne obliczenia i mieć do nich szybki dostęp.

Aby funkcja nie była monadyczna, rozwiązanie polegające na budowaniu nieskończonego leniwego drzewa z odpowiednim sposobem jego indeksowania (jak pokazano w poprzednich postach) spełnia ten cel. Jeśli zrezygnujesz z niemonadycznego charakteru funkcji, możesz użyć standardowych kontenerów asocjacyjnych dostępnych w Haskell w połączeniu z monadami „stanowymi” (takimi jak State lub ST).

Chociaż główną wadą jest to, że otrzymujesz funkcję niemonadyczną, nie musisz już samodzielnie indeksować struktury i możesz po prostu użyć standardowych implementacji kontenerów asocjacyjnych.

Aby to zrobić, musisz najpierw ponownie napisać swoją funkcję, aby akceptowała jakąkolwiek monadę:

fm :: (Integral a, Monad m) => (a -> m a) -> a -> m a
fm _    0 = return 0
fm recf n = do
   recs <- mapM recf $ div n <$> [2, 3, 4]
   return $ max n (sum recs)

W swoich testach nadal możesz zdefiniować funkcję, która nie wykonuje zapamiętywania za pomocą Data.Function.fix, chociaż jest nieco bardziej szczegółowa:

noMemoF :: (Integral n) => n -> n
noMemoF = runIdentity . fix fm

Następnie możesz użyć State monad w połączeniu z Data.Map, aby przyspieszyć działanie:

import qualified Data.Map.Strict as MS

withMemoStMap :: (Integral n) => n -> n
withMemoStMap n = evalState (fm recF n) MS.empty
   where
      recF i = do
         v <- MS.lookup i <$> get
         case v of
            Just v' -> return v' 
            Nothing -> do
               v' <- fm recF i
               modify $ MS.insert i v'
               return v'

Po drobnych zmianach można zamiast tego dostosować kod do pracy z Data.HashMap:

import qualified Data.HashMap.Strict as HMS

withMemoStHMap :: (Integral n, Hashable n) => n -> n
withMemoStHMap n = evalState (fm recF n) HMS.empty
   where
      recF i = do
         v <- HMS.lookup i <$> get
         case v of
            Just v' -> return v' 
            Nothing -> do
               v' <- fm recF i
               modify $ HMS.insert i v'
               return v'

Zamiast trwałych struktur danych możesz także wypróbować zmienne struktury danych (takie jak Data.HashTable) w połączeniu z monadą ST:

import qualified Data.HashTable.ST.Linear as MHM

withMemoMutMap :: (Integral n, Hashable n) => n -> n
withMemoMutMap n = runST $
   do ht <- MHM.new
      recF ht n
   where
      recF ht i = do
         k <- MHM.lookup ht i
         case k of
            Just k' -> return k'
            Nothing -> do 
               k' <- fm (recF ht) i
               MHM.insert ht i k'
               return k'

W porównaniu z implementacją bez zapamiętywania, każda z tych implementacji pozwala, przy dużych nakładach, na uzyskanie wyników w mikro-sekundach, zamiast czekać kilka sekund.

Używając Criterion jako benchmarku, mogłem zauważyć, że implementacja z Data.HashMap faktycznie działała nieco lepiej (około 20%) niż ta z Data.Map i Data.HashTable, dla których czasy były bardzo podobne.

Wyniki testu porównawczego okazały się nieco zaskakujące. Moje początkowe wrażenie było takie, że HashTable przewyższy implementację HashMap, ponieważ jest zmienna. W tej ostatniej implementacji może być ukryty jakiś defekt wydajności.

Quentin
źródło
2
GHC bardzo dobrze radzi sobie z optymalizacją niezmiennych struktur. Intuicja z C nie zawsze się sprawdza.
John Tyree
3

Kilka lat później przyjrzałem się temu i zdałem sobie sprawę, że istnieje prosty sposób na zapamiętanie tego w czasie liniowym przy użyciu zipWithfunkcji pomocniczej:

dilate :: Int -> [x] -> [x]
dilate n xs = replicate n =<< xs

dilatema tę przydatną właściwość dilate n xs !! i == xs !! div i n.

Więc zakładając, że mamy f (0), to upraszcza obliczenia do

fs = f0 : zipWith max [1..] (tail $ fs#/2 .+. fs#/3 .+. fs#/4)
  where (.+.) = zipWith (+)
        infixl 6 .+.
        (#/) = flip dilate
        infixl 7 #/

Wygląda bardzo podobnie do naszego oryginalnego opisu problemu i daje liniowe rozwiązanie ( sum $ take n fszajmie O (n)).

rampion
źródło
2
więc jest to rozwiązanie generatywne (corecursive?) lub programowanie dynamiczne. Biorąc O (1) czasu na każdą wygenerowaną wartość, tak jak zwykle robi to Fibonacci. Wspaniały! Rozwiązanie EKMETT jest podobne do logarytmicznego dużego Fibonacciego, uzyskując duże liczby znacznie szybciej, pomijając wiele operacji między. Czy to prawda?
Will Ness
a może jest bliżej tej dla liczb Hamminga, z trzema wskazówkami wstecznymi do tworzonej sekwencji i różnymi prędkościami dla każdej z nich poruszającymi się wzdłuż niej. naprawdę ładny.
Will Ness
2

Kolejny dodatek do odpowiedzi Edwarda Kmetta: samodzielny przykład:

data NatTrie v = NatTrie (NatTrie v) v (NatTrie v)

memo1 arg_to_index index_to_arg f = (\n -> index nats (arg_to_index n))
  where nats = go 0 1
        go i s = NatTrie (go (i+s) s') (f (index_to_arg i)) (go (i+s') s')
          where s' = 2*s
        index (NatTrie l v r) i
          | i <  0    = f (index_to_arg i)
          | i == 0    = v
          | otherwise = case (i-1) `divMod` 2 of
             (i',0) -> index l i'
             (i',1) -> index r i'

memoNat = memo1 id id 

Użyj go w następujący sposób, aby zapamiętać funkcję z pojedynczą liczbą całkowitą arg (np. Fibonacci):

fib = memoNat f
  where f 0 = 0
        f 1 = 1
        f n = fib (n-1) + fib (n-2)

Buforowane będą tylko wartości argumentów nieujemnych.

Aby również buforować wartości argumentów ujemnych, użyj memoInt, zdefiniowane w następujący sposób:

memoInt = memo1 arg_to_index index_to_arg
  where arg_to_index n
         | n < 0     = -2*n
         | otherwise =  2*n + 1
        index_to_arg i = case i `divMod` 2 of
           (n,0) -> -n
           (n,1) ->  n

Aby buforować wartości dla funkcji z dwoma argumentami całkowitymi, użyj memoIntInt, zdefiniowanych w następujący sposób:

memoIntInt f = memoInt (\n -> memoInt (f n))
Neal Young
źródło
2

Rozwiązanie bez indeksowania, a nie oparte o Edward KMETT.

Oddzielam wspólne poddrzewa do wspólnego rodzica ( f(n/4)jest dzielone między f(n/2)i f(n/4)i f(n/6)jest dzielone między f(2)i f(3)). Zapisując je jako pojedynczą zmienną w rodzicu, obliczenie poddrzewa jest wykonywane raz.

data Tree a =
  Node {datum :: a, child2 :: Tree a, child3 :: Tree a}

f :: Int -> Int
f n = datum root
  where root = f' n Nothing Nothing


-- Pass in the arg
  -- and this node's lifted children (if any).
f' :: Integral a => a -> Maybe (Tree a) -> Maybe (Tree a)-> a
f' 0 _ _ = leaf
    where leaf = Node 0 leaf leaf
f' n m2 m3 = Node d c2 c3
  where
    d = if n < 12 then n
            else max n (d2 + d3 + d4)
    [n2,n3,n4,n6] = map (n `div`) [2,3,4,6]
    [d2,d3,d4,d6] = map datum [c2,c3,c4,c6]
    c2 = case m2 of    -- Check for a passed-in subtree before recursing.
      Just c2' -> c2'
      Nothing -> f' n2 Nothing (Just c6)
    c3 = case m3 of
      Just c3' -> c3'
      Nothing -> f' n3 (Just c6) Nothing
    c4 = child2 c2
    c6 = f' n6 Nothing Nothing

    main =
      print (f 123801)
      -- Should print 248604.

Kod nie jest łatwo rozszerzany do ogólnej funkcji zapamiętywania (przynajmniej nie wiedziałbym, jak to zrobić) i naprawdę musisz pomyśleć, w jaki sposób podproblemy się pokrywają, ale strategia powinna działać dla ogólnych parametrów wielu niecałkowitych . (Wymyśliłem to dla dwóch parametrów ciągów.)

Notatka jest odrzucana po każdym obliczeniu. (Ponownie myślałem o dwóch parametrach ciągów.)

Nie wiem, czy to jest bardziej wydajne niż inne odpowiedzi. Technicznie każde wyszukiwanie składa się z jednego lub dwóch kroków („Spójrz na swoje dziecko lub dziecko”), ale może wymagać dużo dodatkowej pamięci.

Edycja: to rozwiązanie nie jest jeszcze poprawne. Udostępnianie jest niekompletne.

Edycja: Powinno teraz poprawnie udostępniać poddzieci, ale zdałem sobie sprawę, że ten problem ma wiele nietrywialnych współdzielenia: n/2/2/2i n/3/3może być taki sam. Problem nie pasuje do mojej strategii.

leewz
źródło