Haskellで最長増加部分列（LIS）IOArrayやData.Setで試行錯誤

今度もAOJの練習問題から
最長増加部分列の問題に挑む。
最長増加部分列 | 動的計画法 | Aizu Online Judge
 
コードを書く練習と割り切っているので
考え方は素直にググる。

thさんという方のブログに
分かりやすい考え方が載っていたので
これをHaskellで実装する事にした。

;; 最後の要素から順番に、各要素を先頭とするLIS長を求めていく。
;; 計算を具体的に手でやってみると、次のような流れになる（リストに作用させる関数をfとおく）。
(f 15) => 1
(f 7 15) => 2 ;7 < 15 なので (f 15) + 1 => 1 + 1
(f 11 7 15) => 2 ;11 < 15 なので (f 15) + 1 => 1 + 1
(f 3 11 7 15) => 3 ;3より大きい11,7,15で始まるLISはそれぞれ長さが2,2,1。maxの2に1加えて3
(f 13 3 11 7 15) => 2 ;13 < 15 なので (f 15) + 1 => 1 + 1
(f 5 13 3 11 7 15) => 3 ;5より大きい13か11か7で始まるLISが最長で2。それに1加えて3
(f 9 5 13 3 11 7 15) => 3 ;9より大きい13か11で始まるLISが最長で2。それに1加えて3
(f 1 9 5 13 3 11 7 15) => 4 ;1より大きい9や5で始まるLISが最長で3。それに1加えて4
(f 14 1 9 5 13 3 11 7 15) => 2 ;14 < 15 なので (f 15) + 1 => 1 + 1
(f 6 14 1 9 5 13 3 11 7 15) => 4 ;6より大きい数のうち、9で始まるLISが最長で3。それに1加えて4
(f 10 6 14 ...略) => 3 ;10より大きい11,13のLISが最長で2。それに1加えて3
(f 2 10 6 ...略) => 5 ;2より大きい6のLISが最長で4。それに1加えて5
(f 12 2 10 ...略) => 3 ;12より大きい13, 14のLISが最長で4。それに1加えて3
(f 4 12 2 ...略) => 5 ;5より大きい数のうち、6で始まるLISが最長で4。それに1加えて5
(f 8 4 12 ...略) => 4 ;8より大きい12,10,9のLISが最長で3。それに1加えて4
(f 0 8 4 ...略) => 6 ;0より大きい数のうち、2で始まるLISが最長で5。それに1加えて6

Technical Memorandum: 最長増加部分列の長さ in Elisp

Arrayで書いてみる

最初は数字と最長がいくつなのかの2つの数をArrayに収める事を考えてみた。

import Control.Applicative
import Data.Array
import Data,List (foldl')

main = do
  len <- (read :: String -> Int) <$> getLine
  xs <- map (read :: String -> Int) . lines <$> getContents
  let arr = listArray (1,len) (zip xs [0,0..])
  print $ ans $syorinaiyou arr len

ans arr = maximum $ map snd $ elems $ arr

syorinaiyou :: Array Int (Int, Int) -> Int -> Array Int (Int, Int)
syorinaiyou arr len = foldl' keisan arr (reverse [1..len])
    where
      keisan :: Array Int (Int, Int) -> Int -> Array Int (Int, Int)
      keisan arr m
        | m == len = arr // [(m,(a,1))]
        | 1==1     = fst $ foldl' seisa' (arr,1) [m..(len+1)]
        where
          (a,_) = arr ! len
          seisa' (arr,x) t
            | t > len  = (koushin,x)
            | fstarr m < fstarr t = (arr,hikaku)
            | otherwise          = (arr,x)
            where
              fstarr r = fst $ arr ! r
              sndarr r = snd $ arr ! r
              hikaku = max x ((sndarr t)+1)
              koushin = arr // [(m,(fstarr m,x))]  

（叩き台として書いていたので変数名が酷い）

早速提出する

安定のMLE。

        | m == len = arr // [(m,(a,1))]
              koushin = arr // [(m,(fstarr m,x))]  

原因は多分この辺り。
Haskellの配列は更新のコストが高すぎるので
破壊的代入をしないならMap一択なのかもしれない。

という訳でMapで書き直す。

import Control.Applicative
import qualified Data.Map.Strict as Map
import Data.List (foldl')

main = do
  len <- (read :: String -> Int) <$> getLine
  xs <- map (read :: String -> Int) . lines <$> getContents
  let arr = Map.fromList $ zip [1..len] (zip xs [0,0..])
  print $ ans $ syorinaiyou arr len

ans arr = maximum $ map (snd . snd)  $ Map.toList $ arr

syorinaiyou arr len = foldl' keisan arr (reverse [1..len])
    where
      keisan arr m
        | m == len = Map.insert m (a,1) arr
        | 1==1     = fst $ foldl' seisa' (arr,1) [m..(len+1)]
        where
          (a,_) = arr Map.! len
          seisa' (arr,x) t
            | t > len  = (koushin,x)
            | fstarr m < fstarr t = (arr,hikaku)
            | otherwise          = (arr,x)
            where
              fstarr r = fst $ arr Map.! r
              sndarr r = snd $ arr Map.! r
              hikaku = max x ((sndarr t)+1)
              koushin = Map.insert m (fstarr m,x) arr  

小さな修正で済んだ。

早速再提出。

うーむ駄目。

IOArrayを試す

テストケースを見てみると
どうも要素数が1万もあるので、更新云々の問題でもなさそうだ。
まあせっかくなのでIOArrayで再実装したパターンも試してみる。

〜コンパイルエラーと格闘する事 ２日間〜

import Control.Applicative
import Data.Array.IO
import Control.Monad (mapM_)
import Data.IORef

main = do
  len <- (read :: String -> Int) <$> getLine
  xs <- map (read :: String -> Int) . lines <$> getContents
  arr <- newListArray (1,len) (zip xs [0,0..]) :: IO (IOArray Int (Int,Int))
  syorinaiyou arr len
  x <- newIORef 0
  ans arr x
  print =<< readIORef x

ans :: IOArray Int (Int, Int) -> IORef Int -> IO ()
ans arr x = do
  xs <- getElems arr
  writeIORef x (maximum $ map snd  xs)


syorinaiyou :: IOArray Int (Int, Int) -> Int -> IO ()
syorinaiyou arr len = do mapM_ (keisan arr) (reverse [1..len])
    where
      keisan :: IOArray Int (Int, Int) -> Int -> IO ()
      keisan arr m = do
          x <- newIORef 1
          mapM_ (seisa' (arr,x)) [m..(len+1)]
        where
          seisa' :: (IOArray Int (Int, Int), IORef Int) -> Int -> IO ()
          seisa' (arr,x) t
            | t > len  = koushin
            | otherwise = do
              m'' <- readArray arr m
              let m' = fst m''
              t'' <- readArray arr t
              let t' = fst t''
              if m' < t' then hikaku
                         else return ()
            where
              hikaku = do
                x' <- readIORef x
                t'' <- readArray arr t
                let t' = snd t''
                writeIORef x (max x' (t'+1))
              koushin = do
                x' <- readIORef x
                m'' <- readArray arr m
                let m' = fst m''
                writeArray arr m (m',x')  

doだらけのコードがようやく完成。
折角なのでArrayで書いたコードをIOArrayで書き直す時に苦労したポイントをメモしておく。

  arr <- newListArray (1,len) (zip xs [0,0..]) :: IO (IOArray Int (Int,Int))  

まずここ。
型注釈「:: IO (IOArray Int (Int,Int))」がないとコンパイルエラーになる。
IOArrayなのかIOUArrayなのかコンパイラが判別出来ないからだ。

seisa' (arr,x) t
  | t > len  = koushin
  | otherwise = do
    m'' <- readArray arr m
    let m' = fst m''
    t'' <- readArray arr t
    let t' = fst t''
    if m' < t' then hikaku
     else return ()  

もとのコードと比較して一番汚く冗長になったこの部分。
初めはfstarr m < fstarr tとして
fstarr内でIOArrayから取得した値を返そうとしたのだが
非純粋？な関数から純粋な関数に値を渡す事は出来ないようで
仕方なくdoブロックを作ってブロック内で値を取り出してから
if then else で条件別に振り分けた。

後半の方は「t''」みたいな投げやりな変数名になってしまったが動いたので満足。
・・・はせずに、ブログ執筆用も兼ねて変数名で意味が分かるように修正。
こんな1文字の変数名では自分でも意味不明なので。

import Control.Applicative
import Data.Array.IO
import Control.Monad (mapM_)
import Data.IORef

main = do
  len <- (read :: String -> Int) <$> getLine
  xs <- map (read :: String -> Int) . lines <$> getContents
  arr <- newListArray (1,len) (zip xs [1,1..]) :: IO (IOArray Int (Int,Int))
  syorinaiyou arr len
  lengthLIS <- newIORef 0
  ans arr lengthLIS
  print =<< readIORef lengthLIS

ans :: IOArray Int (Int, Int) -> IORef Int -> IO ()
ans arr lengthLIS = do
  xs <- getElems arr
  writeIORef lengthLIS (maximum $ map snd xs)

syorinaiyou :: IOArray Int (Int, Int) -> Int -> IO ()
syorinaiyou arr len = do mapM_ (keisan arr) (reverse [1..len])
    where
      keisan :: IOArray Int (Int, Int) -> Int -> IO ()
      keisan arr hikaku_moto_Idx = do
          mapM_ (seisa' arr) [hikaku_moto_Idx..(len)]
        where
          seisa' :: IOArray Int (Int, Int) -> Int -> IO ()
          seisa' arr hikaku_saki_Idx = do
              (hikaku_moto,_) <- readArray arr hikaku_moto_Idx
              (hikaku_saki,_) <- readArray arr hikaku_saki_Idx
              if hikaku_moto < hikaku_saki then hikaku
                         else return ()
            where
              hikaku = do
                (_,hikaku_saki_LIS)                  <- readArray arr hikaku_saki_Idx
                (hikaku_moto,hikaku_moto_LIS_zantei) <- readArray arr hikaku_moto_Idx
                if hikaku_saki_LIS + 1 > hikaku_moto_LIS_zantei
                   then do
                     writeArray arr hikaku_moto_Idx (hikaku_moto, (hikaku_saki_LIS + 1))
                   else return ()  

だいぶ読みやすくなった。
格好つけても仕方ないので変数名はローマ字にした。
処理の流れが分かりやすくなった事で無断な部分を省く事にも成功。

早速これで提出。

メモリの削減という目的は達成されるも
残念ながらTLE

O2ビルドで強行突破をはかるも少し進んでTLE

秘密兵器Data.Set

そこでData.Setに一度
各要素を先頭とするLIS長を保存し
そこから探索をする事にした。
Data.Set

takeWhileAntitone :: (a -> Bool) -> Set a -> Set a
O(log n). Take while a predicate on the elements holds .

findMax :: Set a -> a
O(log n). The maximal element of a set.

Data.Setはリストに扱い方が似ているが
takeWhileやmaximumをO(log n)で実行できるというスグレモノだ。
早速活用してみる。

修正するべき箇所はここ

      keisan :: IOArray Int (Int, Int) -> Int -> IO ()
      keisan arr hikaku_moto_Idx = do
          mapM_ (seisa' arr) [hikaku_moto_Idx..(len)]
        where
          seisa' :: IOArray Int (Int, Int) -> Int -> IO ()
          seisa' arr hikaku_saki_Idx = do
              (hikaku_moto,_) <- readArray arr hikaku_moto_Idx
              (hikaku_saki,_) <- readArray arr hikaku_saki_Idx
              if hikaku_moto < hikaku_saki then hikaku
                         else return ()
            where
              hikaku = do
                (_,hikaku_saki_LIS) <- readArray arr hikaku_saki_Idx
                (hikaku_moto,hikaku_moto_LIS_zantei) <- readArray arr hikaku_moto_Idx
                if hikaku_saki_LIS + 1 > hikaku_moto_LIS_zantei
                   then do
                     writeArray arr hikaku_moto_Idx (hikaku_moto, (hikaku_saki_LIS + 1))
                   else return ()  

この処理というのは
比較元よりも手前かつ大きな数の中で
その要素を先頭にするLIS長が一番長いものを探し出すというもの。

lookupMax $ takeWhileAntitone (>hikaku_moto) (これまでのLIS長のData.Set)  

こんな完成イメージをもちつつソースを修正していく。

いきなり問題発生

lis.hs:6:57: error:
    Module ‘Data.Set’ does not export ‘dropWhileAntitone’  

なんと古いバージョンのData.SetにはtakeWhileやdropWhileにあたる関数がなかったのだ。
仕方ないのでSet.splitを使う。

> b = Set.fromList $ zip [1..10] [1..10]
> b
fromList [(1,1),(2,2),(3,3),(4,4),(5,5),(6,6),(7,7),(8,8),(9,9),(10,10)]
> Set.split (5,0) b
(fromList [(1,1),(2,2),(3,3),(4,4)],fromList [(5,5),(6,6),(7,7),(8,8),(9,9),(10,10)])
> snd $ Set.split (5,0) b
fromList [(5,5),(6,6),(7,7),(8,8),(9,9),(10,10)]  

これでtakeWhileAntitoneとほぼ同じような動作になる。

IOArrayは必要なくなったので
大幅に書き直したのが下のコード

{-# OPTIONS_GHC -O2 #-}
import Control.Applicative
import Data.Array.Unboxed
import Data.List (foldl')
import Data.Foldable (maximumBy)
import Data.Ord (comparing)
import qualified Data.Set as Set

main = do
  len <- (read :: String -> Int) <$> getLine
  xs <- map (read :: String -> Int) . lines <$> getContents
  let arr = listArray (1,len) xs :: UArray Int Int
  print $ snd $ maximumBy (comparing snd) $ syorinaiyou arr len

syorinaiyou arr len =
  foldl' (keisan arr) Set.empty (reverse [1..len])
    where
      keisan arr set n
        | Set.null set =  Set.insert ((arr ! n),1) set
        | otherwise    = seisa arr set n
        where
          seisa arr set n
            | Set.null set' = Set.insert (hikaku_moto,1) set
            | otherwise     = Set.insert (hikaku_moto, oldLIS + 1) set
            where
              hikaku_moto = arr ! n
              set'   = snd $ Set.split (hikaku_moto + 1 , 0) set
              oldLIS = snd $ maximumBy (comparing snd) set'  

maximumBy (comparing snd)  

この部分はData.Setに収納されたタプルのsndが最大の要素を抽出する処理。
Data.Foldableはfoldableなら何でも適用できる便利な関数が入ったライブラリ。
Data.Foldable

lookupMaxは癖があったので今回は叩き台という事もあって採用しなかった。

早速再提出してみる。

IOArrayの頃よりは進んだが、まだ突破は出来ず・・・

原因を探る

探らなくても薄々わかっている

snd $ maximumBy (comparing snd) set'

この処理がボトルネックになっていると見た。

そこでsetが長くならないように修正

seisa arr set n
    | Set.null set' = Set.insert (hikaku_moto,1) set
    | otherwise     = Set.insert (hikaku_moto, oldLIS + 1) set2
    where
      hikaku_moto = arr ! n
      (set1, set') = Set.split (hikaku_moto + 1 , 0) set
      (old,oldLIS) = maximumBy (comparing snd) set'
      set2         = Set.union set1 (snd (Set.split (old,oldLIS-1) set'))  

再提出

うーん、あと一息。

もう一度修正

 (old,oldLIS) = Set.elemAt 0 set'  

この部分は先頭＝最長なので
maximumByを使う必要はなかった。

再提出

色々と間違えていたらしい。
とりあえず継ぎ接ぎで修正して再提出。

通った！！

最後の方にはヤケクソ感が漂うコードになってしまった。

{-# OPTIONS_GHC -O2 #-}
import Control.Applicative
import Data.Array.Unboxed
import Data.List (foldl')
import Data.Foldable (maximumBy)
import Data.Ord (comparing)
import qualified Data.Set as Set

main = do
  len <- (read :: String -> Int) <$> getLine
  xs <- map (read :: String -> Int) . lines <$> getContents
  let arr = listArray (1,len) xs :: UArray Int Int
  print $ snd $ maximumBy (comparing snd) $ syorinaiyou arr len

syorinaiyou arr len =
  foldl' (keisan arr) Set.empty (reverse [1..len])
    where
      keisan arr set n
        | Set.null set = Set.insert ((arr ! n),1) set
        | otherwise    = seisa arr set n
        where
          seisa arr set n
            | Set.null set'      = Set.insert (hikaku_moto,1) set
            | Set.size set1 <= 1 =  Set.insert (hikaku_moto, oldLIS + 1) set'
            | otherwise = Set.insert (hikaku_moto, oldLIS + 1) set2
            where
              hikaku_moto = arr ! n
              (set1, set') = Set.split (hikaku_moto + 1 , 0) set
              (old,oldLIS) = Set.elemAt 0 set'
              set2         = Set.union (Set.delete (Set.elemAt ((Set.size set1)-1) set1) set1) (snd (Set.split (old,oldLIS-1) set'))  

一応解決。
最終的には破壊的代入とか要らなかった。
IOArrayの為に格闘していた時間は一体・・・