在哈斯克尔纪念?

关于如何在哈斯克尔有效地解决以下函数的任何建议,对于大数 (n > 108)

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

我在哈斯克尔见过用记忆法解决斐波那契问题的例子 数字,包括计算(惰性)所有斐波那契数字 但是在这种情况下,对于给定的 n,我们只需要 计算很少的中间结果。

谢谢

26496 次浏览

虽然不是最有效的方法,但是可以记住:

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

当请求 f !! 144时,将检查 f !! 143是否存在,但不计算其确切值。它仍然被设置为一些未知的计算结果。计算出的唯一精确值就是所需的数值。

所以最初,到底计算了多少,程序一无所知。

f = ....

当我们提出请求时,它开始做一些模式匹配:

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 : ...

现在它开始计算

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

这会递归地对 f 做出另一个要求,所以我们计算

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

现在我们可以慢慢回去了

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

也就是说程序现在知道了:

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 : ...

不断涌现:

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

也就是说程序现在知道了:

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

现在我们继续计算 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

也就是说程序现在知道了:

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

现在我们继续计算 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

也就是说程序现在知道了:

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

因此,计算是相当懒惰地进行的。程序知道 f !! 8的某个值存在,它等于 g 8,但它不知道 g 8是什么。

我们可以非常有效地做到这一点,通过建立一个结构,我们可以在次线性时间索引。

但首先,

{-# LANGUAGE BangPatterns #-}


import Data.Function (fix)

让我们定义 f,但是使用“开放递归”而不是直接调用它自己。

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)

通过使用 fix f,您可以得到一个未备忘的 f

这将使您可以通过调用(例如: fix f 123 = 144)来测试 f对于 f的小值所起的作用

我们可以通过定义:

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


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

这表现得还可以,并且用一些能够记录中间结果的东西替代了原本需要花费 O (n ^ 3)时间的东西。

但仅仅是索引就能找到 mf的制表答案仍然需要线性时间。这意味着结果如下:

*Main Data.List> faster_f 123801
248604

是可以忍受的,但是结果没有比这更好的了。我们可以做得更好!

首先,让我们定义一个无限树:

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)

然后我们将定义一种方法来索引它,这样我们就可以在 O (log n)时间内找到一个索引为 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

... 我们可能会发现一棵满是自然数的树是很方便的,这样我们就不必去摆弄那些指数了:

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

因为我们可以建立索引,所以你可以直接把一棵树转换成一个列表:

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

你可以通过验证 toList nats给你的 [0..]来检查到目前为止的工作

现在,

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


fastest_f :: Int -> Int
fastest_f = index f_tree

可以像上面的 list 一样工作,但是不需要花费线性时间来找到每个节点,而是可以用对数时间来追踪它。

结果要快得多:

*Main> fastest_f 12380192300
67652175206


*Main> fastest_f 12793129379123
120695231674999

事实上,它是如此之快,你可以通过并取代 IntInteger以上,并得到可笑的大答案几乎瞬间

*Main> fastest_f' 1230891823091823018203123
93721573993600178112200489


*Main> fastest_f' 12308918230918230182031231231293810923
11097012733777002208302545289166620866358

对于实现基于树的制表的开箱即用库,请使用 记忆:

$ stack repl --package MemoTrie
Prelude> import Data.MemoTrie
Prelude Data.MemoTrie> :set -XLambdaCase
Prelude Data.MemoTrie> :{
Prelude Data.MemoTrie| fastest_f' :: Integer -> Integer
Prelude Data.MemoTrie| fastest_f' = memo $ \case
Prelude Data.MemoTrie|   0 -> 0
Prelude Data.MemoTrie|   n -> max n (fastest_f'(n `div` 2) + fastest_f'(n `div` 3) + fastest_f'(n `div` 4))
Prelude Data.MemoTrie| :}
Prelude Data.MemoTrie> fastest_f' 12308918230918230182031231231293810923
11097012733777002208302545289166620866358

这是对爱德华 · 克米特精彩回答的补充。

当我尝试他的代码,natsindex的定义似乎相当神秘,所以我写了一个替代版本,我觉得更容易理解。

我用 index'nats'来定义 indexnats

index' t n[1..]范围内定义。(回想一下,index t是在 [0..]范围内定义的。)它的工作原理是将 n当作一串位来搜索树,然后反过来读取这些位。如果这个位是 1,那么它取右边的分支。如果这个位是 0,它采用左侧分支。当它到达最后一个位(必须是 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'

正如 nats是为 index定义的,因此 index nats n == n总是为真,nats'也是为 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

现在,natsindex只是 nats'index',但是值偏移了1:

index t n = index' t (n+1)
nats = fmap (\n -> n-1) nats'

Edward 的答案 是一个非常棒的宝石,我已经复制了它,并提供了 memoListmemoTree组合子的实现,它们以开放递归的形式制表了一个函数。

{-# 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

对于爱德华•克梅特(Edward Kmett)的回答,还有一个补充: 一个独立的例子:

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

使用它来制表一个带有单个整数 arg 的函数(例如 fibonacci) :

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

只缓存非负参数的值。

要缓存负参数的值,请使用 memoInt,定义如下:

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

为了缓存具有两个整数参数的函数的值,使用 memoIntInt,定义如下:

memoIntInt f = memoInt (\n -> memoInt (f n))

正如 Edward Kmett 的回答所述,为了加快速度,您需要缓存昂贵的计算并能够快速访问它们。

为了保持函数不是一元的,构建无限延迟树的解决方案,以适当的方式对其进行索引(如前面的帖子所示) ,实现了这一目标。如果你放弃了这个函数的非单子性质,你可以使用哈斯克尔的标准关联容器和“ State-like”monads (比如 State 或 ST)结合使用。

虽然主要的缺点是你得到了一个非单子函数,你不必再自己索引结构,只需使用关联容器的标准实现。

为此,您首先需要重写函数来接受任何类型的 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)

对于您的测试,您仍然可以定义一个不使用 Data 进行制表的函数。Fix,尽管它有点冗长:

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

然后你可以结合使用 State monad 和 Data:

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'

只要稍加修改,您就可以调整代码以适应 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'

您也可以尝试使用可变的数据结构(如 Data)来代替持久性数据结构。HashTable)与 ST monad 的组合:

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'

与没有任何制表的实现相比,这些实现中的任何一个都允许您在巨大的输入中以微秒的速度获得结果,而不必等待几秒钟。

使用 Criterion 作为基准,我可以观察到使用 Data 的实现。HashMap 实际上比 Data 执行得稍微好一点(大约20%)。地图及数据。计时非常相似的 HashTable。

我发现基准测试的结果有点令人吃惊。我最初的感觉是,HashTable 的性能会优于 HashMap 实现,因为它是可变的。在最后一个实现中可能隐藏了一些性能缺陷。

几年后,我看到了这个,并意识到有一个简单的方法可以使用 zipWith和一个 helper 函数在线性时间内记录它:

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

dilate具有方便的属性 dilate n xs !! i == xs !! div i n

假设我们得到 f (0) ,这就简化了计算

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

看起来很像我们最初的问题描述,并给出一个线性解决方案(sum $ take n fs将采取 O (n))。

一个没有索引的解决方案,而且不是基于 Edward KMETT 的。

我将公共子树分解为公共父树(f(n/4)f(n/2)f(n/4)之间共享,f(n/6)f(2)f(3)之间共享)。通过将它们保存为父变量中的单个变量,子树的计算就完成了一次。

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.

这段代码不容易扩展到一个通用的制表函数(至少,我不知道怎么做) ,而且你真的必须考虑子问题是如何重叠的,但是 策略应该适用于一般的多个非整数参数。(我想出了两个字符串参数。)

在每次计算之后都会丢弃这个备忘录(同样,我考虑的是两个字符串参数)

我不知道这是否比其他答案更有效。从技术上讲,每个查找只有一到两个步骤(“查看您的孩子或您孩子的孩子”) ,但是可能会使用大量额外的内存。

编辑: 这个解决方案还不正确。共享是不完整的。

编辑: 现在它应该正确地共享子子子元素,但是我意识到这个问题有很多重要的共享: n/2/2/2n/3/3可能是相同的。这个问题不适合我的策略。