Отладка узких мест производительности для алгоритма наиболее длинных общих подпоследовательностей

Я пишу длинный общий алгоритм подпоследовательности в Haskell с использованием векторной библиотеки и государственной монады (для инкапсуляции самой императивной и изменчивой природы Miller O (NP) алгоритм). Я уже написал это в C для какого-то проекта, в котором я нуждался, и теперь я пишу его в Haskell, чтобы изучить, как написать эти императивные алгоритмы grid-walk с хорошей производительностью, которые соответствуют C. Версия, которую я написал с незапакованными векторами примерно в 4 раза медленнее, чем версия C для одних и тех же входов (и скомпилирована с правильными флагами оптимизации), я использовал как системное таксовое время, так и методы Criterion для проверки относительных временных измерений между версиями Haskell и C и одинаковыми типами данных, как больших и крошечные входы). Я пытался выяснить, где могут быть проблемы с производительностью, и по достоинству оценит обратную связь - возможно, есть известная проблема с производительностью, с которой я, возможно, столкнулся здесь, особенно в векторной библиотеке, которую я сильно использую здесь.

В моем коде у меня есть одна функция, называемая gridWalk, которая вызывается чаще всего, а также выполняет большую часть работы. Вероятно, замедление производительности будет, но я не могу понять, что это может быть. Полный код Haskell здесь. Фрагменты кода ниже:

import Data.Vector.Unboxed.Mutable as MU
import Data.Vector.Unboxed as U hiding (mapM_)
import Control.Monad.ST as ST
import Control.Monad.Primitive (PrimState)
import Control.Monad (when) 
import Data.STRef (newSTRef, modifySTRef, readSTRef)
import Data.Int


type MVI1 s  = MVector (PrimState (ST s)) Int

cmp :: U.Vector Int32 -> U.Vector Int32 -> Int -> Int -> Int
cmp a b i j = go 0 i j
               where
                 n = U.length a
                 m = U.length b
                 go !len !i !j| (i<n) && (j<m) && ((unsafeIndex a i) == (unsafeIndex b j)) = go (len+1) (i+1) (j+1)
                                    | otherwise = len

-- function to find previous y on diagonal k for furthest point 
findYP :: MVI1 s -> Int -> Int -> ST s (Int,Int)
findYP fp k offset = do
              let k0 = k+offset-1
                  k1 = k+offset+1
              y0 <- MU.unsafeRead fp k0 >>= \x -> return $ 1+x
              y1 <- MU.unsafeRead fp k1
              if y0 > y1 then return (k0,y0)
              else return (k1,y1)
{-#INLINE findYP #-}

gridWalk :: Vector Int32 -> Vector Int32 -> MVI1 s -> Int -> (Vector Int32 -> Vector Int32 -> Int -> Int -> Int) -> ST s ()
gridWalk a b fp !k cmp = {-#SCC gridWalk #-} do
   let !offset = 1+U.length a
   (!kp,!yp) <- {-#SCC findYP #-} findYP fp k offset                          
   let xp = yp-k
       len = {-#SCC cmp #-} cmp a b xp yp
       x = xp+len
       y = yp+len

   {-#SCC "updateFP" #-} MU.unsafeWrite fp (k+offset) y  
   return ()
{-#INLINE gridWalk #-}

-- The function below executes ct times, and updates furthest point as they are found during furthest point search
findSnakes :: Vector Int32 -> Vector Int32 -> MVI1 s ->  Int -> Int -> (Vector Int32 -> Vector Int32 -> Int -> Int -> Int) -> (Int -> Int -> Int) -> ST s ()
findSnakes a b fp !k !ct cmp op = {-#SCC findSnakes #-} U.forM_ (U.fromList [0..ct-1]) (\x -> gridWalk a b fp (op k x) cmp)
{-#INLINE findSnakes #-}

Я добавил некоторые аннотации для МВЗ и выполнил профилирование с помощью определенного ввода LCS для тестирования. Вот что я получаю:

  total time  =        2.39 secs   (2394 ticks @ 1000 us, 1 processor)
  total alloc = 4,612,756,880 bytes  (excludes profiling overheads)

COST CENTRE MODULE    %time %alloc

gridWalk    Main       67.5   52.7
findSnakes  Main       23.2   27.8
cmp         Main        4.2    0.0
findYP      Main        3.5   19.4
updateFP    Main        1.6    0.0


                                                         individual     inherited
COST CENTRE    MODULE                  no.     entries  %time %alloc   %time %alloc

MAIN           MAIN                     64           0    0.0    0.0   100.0  100.0
 main          Main                    129           0    0.0    0.0     0.0    0.0
 CAF           Main                    127           0    0.0    0.0   100.0  100.0
  findSnakes   Main                    141           0    0.0    0.0     0.0    0.0
  main         Main                    128           1    0.0    0.0   100.0  100.0
   findSnakes  Main                    138           0    0.0    0.0     0.0    0.0
    gridWalk   Main                    139           0    0.0    0.0     0.0    0.0
     cmp       Main                    140           0    0.0    0.0     0.0    0.0
   while       Main                    132        4001    0.1    0.0   100.0  100.0
    findSnakes Main                    133       12000   23.2   27.8    99.9   99.9
     gridWalk  Main                    134    16004000   67.5   52.7    76.7   72.2
      cmp      Main                    137    16004000    4.2    0.0     4.2    0.0
      updateFP Main                    136    16004000    1.6    0.0     1.6    0.0
      findYP   Main                    135    16004000    3.5   19.4     3.5   19.4
   newVI1      Main                    130           1    0.0    0.0     0.0    0.0
   newVI1.\   Main                    131        8004    0.0    0.0     0.0    0.0
 CAF           GHC.Conc.Signal         112           0    0.0    0.0     0.0    0.0
 CAF           GHC.IO.Encoding         104           0    0.0    0.0     0.0    0.0
 CAF           GHC.IO.Encoding.Iconv   102           0    0.0    0.0     0.0    0.0
 CAF           GHC.IO.Handle.FD         95           0    0.0    0.0     0.0    0.0

Если я правильно интерпретирую выход профилирования (и предполагая, что из-за профилирования не слишком много искажений), gridWalk занимает большую часть времени, но основные функции cmp и findYP, которые делают тяжелый подъем в gridWalk, кажется, занимает очень мало времени в отчете профилирования. Итак, возможно, узкое место находится в оболочке forM_, которую функция findSnakes использует для вызова gridWalk? Профиль кучи тоже кажется очень чистым:

Читая ядро, ничего действительно не выпрыгивает. Я думал, что некоторые значения во внутренних циклах могут быть в коробке, но я не вижу их в ядре. Я надеюсь, что проблема производительности связана с чем-то простым, что я пропустил.

Обновление

Следуя предложению @DanielFischer, я заменил forM_ of Data.Vector.Unboxed на функцию Control.Monad в findSnakes, которая улучшила производительность с 4x до 2,5x версии C. В настоящее время версии Haskell и C размещены здесь, если вы хотите попробовать их.

Я все еще копаю ядро, чтобы увидеть, где узкие места. gridWalk чаще всего называется функцией, и для того, чтобы она хорошо выполнялась, lcsh должен сократить цикл whileM_ до хорошего итеративного внутреннего цикла проверки состояния и вложенного кода findSnakes. Я подозреваю, что в сборке это не относится к циклу whileM_, но поскольку я не очень хорошо разбираюсь в переводе ядра и поиске в GHC-функциях, связанных с именами в сборке, я думаю, что это просто вопрос терпеливо заткнуться проблема, пока я не выясню это. Между тем, если есть какие-либо указатели на исправления производительности, они будут оценены.

Другая возможность, о которой я могу думать, - это накладные расходы на проверку кучи во время вызовов функций. Как видно из отчета профилирования, gridWalk называется 16004000 раз. Предполагая, что 6 циклов для проверки кучи (я предполагаю, что это меньше, но все же допустим, что), оно составляет ~ 0,02 секунды в ящике 3,33 ГГц для 96024000 циклов.

Кроме того, некоторые номера производительности:

Haskell code (GHC 7.6.1 x86_64): Это было ~ 0,25 с до forM_.

 time ./T
1

real    0m0.150s
user    0m0.145s
sys     0m0.003s

C code (gcc 4.7.2 x86_64):

time ./test
1

real    0m0.065s
user    0m0.063s
sys     0m0.000s

Обновление 2:

Обновленный код здесь. Использование STUArray также не изменяет числа. Производительность около 1,5x на Mac OS X (x86_64,ghc7.6.1), довольно похоже на то, что @DanielFischer сообщила в Linux.

Код Haskell:

$ time ./Diff
1

real    0m0.087s
user    0m0.084s
sys 0m0.003s

C-код:

$ time ./test
1

real    0m0.056s
user    0m0.053s
sys 0m0.002s

Взгляд на cmm, вызов является хвостовым рекурсивным и превращается в цикл через llvm. Но каждая новая итерация, как представляется, выделяет новые значения, которые также вызывают проверку кучи, и поэтому может объяснить разницу в производительности. Я должен думать о том, как написать хвостовую рекурсию таким образом, чтобы никакие значения не выделялись в итерациях, избегая ошибок кучи и распределения.

Ответ 1

Вы делаете огромный удар в

U.forM_ (U.fromList [0..ct-1])

в findSnakes. Я убежден, что это не должно произойти (билет?), Но он выделяет новый Vector для перемещения каждый раз, когда вызывается findSnakes. Если вы используете

Control.Monad.forM_ [0 .. ct-1]

время работы примерно вдвое, а распределение здесь примерно на 500. (GHC оптимизирует C.M.forM_ [0 :: Int .. limit] ну, список исключается, а остальное - это в основном цикл). Вы можете сделать немного лучше, написав цикл самостоятельно.

Некоторые вещи, которые вызывают беспорядочное размножение/размер кода, не повредив производительность, - это

неиспользуемый аргумент Bool lcsh
аргумент cmp для findSnakes и gridWalk; если они никогда не вызываются с другим сравнением, чем верхний уровень cmp, этот аргумент приводит к ненужному дублированию кода.
общий тип while; специализируясь на используемом типе ST s Bool -> ST s () -> ST s (), уменьшает выделение (много), а также время работы (немного, но заметно здесь).

Общее слово о профилировании: Компиляция программы для профилирования подавляет множество оптимизаций. В частности, для библиотек, таких как Vector, bytestring или text, которые сильно используют синтез, профилирование часто приводит к ошибочным результатам.

Например, ваш исходный код производит здесь

    total time  =        3.42 secs   (3415 ticks @ 1000 us, 1 processor)
    total alloc = 4,612,756,880 bytes  (excludes profiling overheads)

COST CENTRE MODULE    %time %alloc  ticks     bytes

gridWalk    Main       63.7   52.7   2176 2432608000
findSnakes  Main       20.0   27.8    682 1281440080
cmp         Main        9.2    0.0    313        16
findYP      Main        4.2   19.4    144 896224000
updateFP    Main        2.7    0.0     91         0

Просто добавление перехвата привязки len в gridWalk ничего не меняет в непрофилирующей версии, но для версии профилирования

    total time  =        2.98 secs   (2985 ticks @ 1000 us, 1 processor)
    total alloc = 3,204,404,880 bytes  (excludes profiling overheads)

COST CENTRE MODULE    %time %alloc  ticks     bytes

gridWalk    Main       63.0   32.0   1881 1024256000
findSnakes  Main       22.2   40.0    663 1281440080
cmp         Main        7.2    0.0    214        16
findYP      Main        4.7   28.0    140 896224000
updateFP    Main        2.7    0.0     82         0

это имеет большое значение. Для версии, включая упомянутые выше изменения (и перерыв на len в gridWalk), в профилирующей версии указано

total alloc = 1,923,412,776 bytes  (excludes profiling overheads)

но непрофилирующая версия

     1,814,424 bytes allocated in the heap
        10,808 bytes copied during GC
        49,064 bytes maximum residency (2 sample(s))
        25,912 bytes maximum slop
             1 MB total memory in use (0 MB lost due to fragmentation)

                                  Tot time (elapsed)  Avg pause  Max pause
Gen  0         2 colls,     0 par    0.00s    0.00s     0.0000s    0.0000s
Gen  1         2 colls,     0 par    0.00s    0.00s     0.0001s    0.0001s

INIT    time    0.00s  (  0.00s elapsed)
MUT     time    0.12s  (  0.12s elapsed)
GC      time    0.00s  (  0.00s elapsed)
EXIT    time    0.00s  (  0.00s elapsed)
Total   time    0.12s  (  0.12s elapsed)

говорит, что он выделяется в 1000 раз меньше, чем профилирующая версия.

Для кода Vector и друзей, более надежного для определения узких мест, чем профилирование (к сожалению, гораздо более трудоемкий и сложный процесс) изучает сгенерированное ядро (или сборку, если вы владеете чтением этого).

Что касается обновления, C работает немного медленнее на моем ящике (gcc-4.7.2, -O3)

$ time ./miltest1

real    0m0.074s
user    0m0.073s
sys     0m0.001s

но Haskell о том же

$ time ./hsmiller
1

real    0m0.151s
user    0m0.149s
sys     0m0.001s

Это происходит немного быстрее при компиляции через бэкэнд LLVM:

$ time ./hsmiller1

real    0m0.131s
user    0m0.129s
sys     0m0.001s

И когда мы заменим forM_ на ручной цикл,

findSnakes a b fp !k !ct op = go 0
  where
    go x
        | x < ct    = gridWalk a b fp (op k x) >> go (x+1)
        | otherwise = return ()

он становится немного быстрее,

$ time ./hsmiller
1

real    0m0.124s
user    0m0.121s
sys     0m0.002s

соотв. через LLVM:

$ time ./hsmiller
1

real    0m0.108s
user    0m0.107s
sys     0m0.000s

В общем, сгенерированное ядро выглядит прекрасно, одно небольшое раздражение было

Main.$wa
  :: forall s.
     GHC.Prim.Int#
     -> GHC.Types.Int
     -> GHC.Prim.State# s
     -> (# GHC.Prim.State# s, Main.MVI1 s #)

и немного окольной реализации. Это фиксируется, делая newVI1 strict во втором аргументе,

newVI1 n !x = do

Поскольку это часто не вызывается, влияние на производительность, конечно, незначительно.

Мясо является ядром для lcsh, и это выглядит не так уж плохо. Единственные вещи в коробке - это Int, читаемые/записываемые в STRef, и это неизбежно. Что не очень приятно, так это то, что ядро содержит много дубликатов кода, но, по моему опыту, это редко является реальной проблемой производительности, и не все дублированные коды выживают при генерации кода.

и для того, чтобы он хорошо выполнялся, lcsh должен сократить цикл whileM_ до хорошего итеративного внутреннего цикла проверки состояния и вложенного кода findSnakes.

Вы получаете внутренний цикл, когда вы добавляете прагму INLINE в whileM_, но этот цикл не является приятным, и в этом случае он намного медленнее, чем whileM_ out-of-line (I "Не уверен, что это связано только с размером кода, но это может быть).