Haskell Repa трафарет хаки

Эта проблема

Я пытаюсь понять, как работает Repa, и я пытался описать пример кода "Blur" из пакета Repa examples. Код использует stencil2 Quasi Quote:

[stencil2|   2  4  5  4  2
             4  9 12  9  4
             5 12 15 12  5
             4  9 12  9  4
             2  4  5  4  2 |]

Что просто TemplateHaskell фрагмент, который генерирует функцию:

makeStencil2 5 5 coeffs where
     {-# INLINE[~0] coeffs #-}
     coeffs = \ ix -> case ix of
                      Z :. -2 :. -2 -> Just 2
                      Z :. -2 :. -1 -> Just 4
                      Z :. -2 :. 0 -> Just 5
                      Z :. -2 :. 1 -> Just 4
                      Z :. -2 :. 2 -> Just 2
                      [...]
                      _ -> Nothing

Можно использовать TH, но я бы хотел сохранить значения в массиве Repa, поэтому я изменил код для использования вместо него массива Repa, но мой код работает в 2 раза медленнее по сравнению с исходным.

Некоторые модные заметки

Я заметил, что авторы Repa используют жестко закодированную матрицу значений 7 на 7 для получения коэффициентов: http://hackage.haskell.org/package/repa-3.2.3.3/docs/src/Data-Array-Repa-Stencil-Dim2.html (см.: template7x7)

Вопросы

  1. Я хочу спросить вас, почему он не оптимизирован как оригинальный и как мы можем это исправить? Я хочу написать "свернутую" функцию, которая позволит мне выполнять свертывание трафарета (массива репа) над изображением.
  2. Нужно ли нам использовать такие жестко закодированные матрицы, чтобы GHC оптимизировал код? Неужели нет способа создать быстрый код на Haskell без использования таких "хаков"?

Код

Оригинальная функция размытия:

blur    :: Monad m => Int -> Array U DIM2 Double -> m (Array U DIM2 Double)
blur !iterations arrInit
 = go iterations arrInit
 where  go !0 !arr = return arr
        go !n !arr  
         = do   arr'    <- computeP
                         $ A.smap (/ 159)
                         $ forStencil2 BoundClamp arr
                           [stencil2|   2  4  5  4  2
                                        4  9 12  9  4
                                        5 12 15 12  5
                                        4  9 12  9  4
                                        2  4  5  4  2 |]
                go (n-1) arr'

моя функция размытия:

blur !iterations arrInit = go iterations arrInit
    where 
          stencilx7 = fromListUnboxed (Z :. 7 :. 7) 
                    [  0,  0,  0,  0,  0,  0, 0
                    ,  0,  2,  4,  5,  4,  2, 0
                    ,  0,  4,  9, 12,  9,  4, 0
                    ,  0,  5, 12, 15, 12,  5, 0
                    ,  0,  4,  9, 12,  9,  4, 0
                    ,  0,  2,  4,  5,  4,  2, 0
                    ,  0,  0,  0,  0,  0,  0, 0
                    ] :: Array U DIM2 Int
          magicf (Z :. x :. y) = Just $ fromIntegral $ unsafeIndex stencilx7 (Z:. (x+3) :. (y+3))
          go !0 !arr = return arr
          go !n !arr  
           = do   
                  arr'    <- computeP
                           $ A.smap (/ 159)
                           $ A.forStencil2 BoundClamp arr 
                            $ makeStencil2 5 5 magicf
                  go (n-1) arr'

Остальной код:

{-# LANGUAGE PackageImports, BangPatterns, TemplateHaskell, QuasiQuotes #-}
{-# OPTIONS -Wall -fno-warn-missing-signatures -fno-warn-incomplete-patterns #-}

import Data.List
import Control.Monad
import System.Environment
import Data.Word
import Data.Array.Repa.IO.BMP
import Data.Array.Repa.IO.Timing
import Data.Array.Repa                          as A
import qualified Data.Array.Repa.Repr.Unboxed   as U
import Data.Array.Repa.Stencil                  as A
import Data.Array.Repa.Stencil.Dim2             as A
import Prelude                                  as P

main 
 = do   args    <- getArgs
        case args of
         [iterations, fileIn, fileOut]  -> run (read iterations) fileIn fileOut
         _                              -> usage

usage   = putStr $ unlines
        [ "repa-blur <iterations::Int> <fileIn.bmp> <fileOut.bmp>" ]


-- | Perform the blur.
run :: Int -> FilePath -> FilePath -> IO ()
run iterations fileIn fileOut
 = do   arrRGB  <- liftM (either (error . show) id) 
                $  readImageFromBMP fileIn

        arrRGB `deepSeqArray` return ()
        let (arrRed, arrGreen, arrBlue) = U.unzip3 arrRGB
        let comps                       = [arrRed, arrGreen, arrBlue]

        (comps', tElapsed)
         <- time $ P.mapM (process iterations) comps

        putStr $ prettyTime tElapsed

        let [arrRed', arrGreen', arrBlue'] = comps'
        writeImageToBMP fileOut
                (U.zip3 arrRed' arrGreen' arrBlue')


process :: Monad m => Int -> Array U DIM2 Word8 -> m (Array U DIM2 Word8)
process iterations 
        = promote >=> blur iterations >=> demote
{-# NOINLINE process #-}


promote :: Monad m => Array U DIM2 Word8 -> m (Array U DIM2 Double)
promote arr
 = computeP $ A.map ffs arr
 where  {-# INLINE ffs #-}
        ffs     :: Word8 -> Double
        ffs x   =  fromIntegral (fromIntegral x :: Int)
{-# NOINLINE promote #-}


demote  :: Monad m => Array U DIM2 Double -> m (Array U DIM2 Word8)
demote arr
 = computeP $ A.map ffs arr

 where  {-# INLINE ffs #-}
        ffs     :: Double -> Word8
        ffs x   =  fromIntegral (truncate x :: Int)

Компилировать с: ghc -O2 -threaded -fllvm -fforce-recomp Main.hs -ddump-splices

1 ответ

  1. Считывание коэффициентов свертки из массива теоретически не может быть таким же быстрым, как пайка констант прямо в скомпилированном коде, потому что последний подход ничего не стоит на уровне машины.

  2. Нет, GHC может сваривать статические трафареты произвольного размера. Смотрите мою реализацию статических сверток с fixed-vector с лямбды:

    [dim2St| 1   2   1
             0   0   0
            -1  -2  -1 |]
    -->
    Dim2Stencil
     n3
     n3
     (VecList
        [VecList
           [\ acc a -> return (acc + a),
            \ acc a -> (return $ (acc + (2 * a))),
            \ acc a -> return (acc + a)],
         VecList
           [\ acc _ -> return acc,
            \ acc _ -> return acc,
            \ acc _ -> return acc],
         VecList
           [\ acc a -> return (acc - a),
            \ acc a -> (return $ (acc + (-2 * a))),
            \ acc a -> return (acc - a)]])
     (\ acc a reduce -> reduce acc a)
     (return 0)
    
Другие вопросы по тегам