Возьмем самую быструю реализацию задачи из предыдущей статьи.
import qualified Data.List as L
import qualified Data.List.Ordered as LO
import Data.Ord (comparing)
import Control.Arrow
alphabet :: [Char]
alphabet = "ACGT"
subSeqs :: Int -> String -> [String]
subSeqs n = takeWhile ((==n) . length) . map (take n) . L.tails
nearSeqs :: Int -> String -> [String]
nearSeqs n = LO.nubSort . nearSeqs' n
where nearSeqs' 0 s = [s]
nearSeqs' n s =
concatMap (\(a, b) -> map (a++) b) $
concatMap (\m -> [(z ++ [x], nearSeqs' (n - 1) $ tail z') |
x <- alphabet, let (z, z') = splitAt m s])
[0 .. length s - 1]
grpSeqs :: Int -> [String] -> [[String]]
grpSeqs n = L.sortBy (flip $ comparing length) . L.group . L.sort . seqs
where seqs = concatMap (\(s, m) -> concat . replicate m . nearSeqs n $ s) .
map (head &&& length) . L.group . L.sort
main = do
let seqs = ["ACGTTGCATGTCGCATGATGCATGAGAGCT",
"CACAGTAGGCGCCGGCACACACAGCCCCGGGCCCCGGGCCGCCCCGGGCCGGC\
\GGCCGCCGGCGCCGGCACACCGGCACAGCCGTACCGGCACAGTAGTACCGGCC\
\GGCCGGCACACCGGCACACCGGGTACACACCGGGGCGCACACACAGGCGGGCG\
\CCGGGCCCCGGGCCGTACCGGGCCGCCGGCGGCCCACAGGCGCCGGCACAGTA\
\CCGGCACACACAGTAGCCCACACACAGGCGGGCGGTAGCCGGCGCACACACAC\
\ACAGTAGGCGCACAGCCGCCCACACACACCGGCCGGCCGGCACAGGCGGGCGG\
\GCGCACACACACCGGCACAGTAGTAGGCGGCCGGCGCACAGCC"]
mapM_ print' $ take 10 $ grpSeqs 2 $ subSeqs 10 $ seqs !! 1
where print' = print . (head &&& length)
Я все же решил использовать Data.List.Ordered.nubSort для удаления дубликатов в nearSeqs, поэтому программа должна работать еще чуточку быстрее. В среднем новая версия работает порядка 0.51 секунды, то есть действительно быстрее на ∼0.03 секунды.
Отлично. Давайте отпрофилируем нашу программу и посмотрим, что можно сделать, чтобы она работала еще быстрее. Компилируем с поддержкой профилирования.
ghc --make -rtsopts -prof -auto-all -caf-all -fforce-recomp gen
Запускаем программу с записью данных профилирования в файл.
./gen +RTS -p
Смотрим результаты профилирования в файле gen.prof.
Mon Nov 10 16:38 2014 Time and Allocation Profiling Report (Final)
gen +RTS -p -RTS
total time = 0.77 secs (766 ticks @ 1000 us, 1 processor)
total alloc = 593,276,296 bytes (excludes profiling overheads)
COST CENTRE MODULE %time %alloc
grpSeqs Main 52.0 24.6
nearSeqs Main 20.6 10.7
nearSeqs.nearSeqs'.\ Main 8.4 19.9
nearSeqs.nearSeqs'.\ Main 7.8 19.5
nearSeqs.nearSeqs' Main 4.6 9.9
nearSeqs.nearSeqs'.\.(...) Main 4.3 12.3
grpSeqs.seqs Main 0.9 1.5
grpSeqs.seqs.\ Main 0.4 1.5
individual inherited
COST CENTRE MODULE no. entries %time %alloc %time %alloc
MAIN MAIN 51 0 0.0 0.0 100.0 100.0
CAF:main :Main 101 0 0.0 0.0 0.0 0.0
main Main 116 0 0.0 0.0 0.0 0.0
main.print' Main 118 0 0.0 0.0 0.0 0.0
CAF:main Main 100 0 0.0 0.0 100.0 100.0
main Main 102 1 0.0 0.0 100.0 100.0
main.print' Main 117 1 0.0 0.0 0.0 0.0
main.seqs Main 106 1 0.0 0.0 0.0 0.0
subSeqs Main 105 1 0.0 0.0 0.0 0.0
grpSeqs Main 103 1 52.0 24.6 100.0 99.9
grpSeqs.seqs Main 104 1 0.9 1.5 48.0 75.3
grpSeqs.seqs.\ Main 107 254 0.4 1.5 47.1 73.8
nearSeqs Main 108 254 20.6 10.7 46.7 72.3
nearSeqs.nearSeqs' Main 109 193294 4.6 9.9 26.1 61.7
nearSeqs.nearSeqs'.\ Main 112 193040 8.4 19.9 8.4 19.9
nearSeqs.nearSeqs'.\ Main 110 48260 7.8 19.5 13.2 31.8
nearSeqs.nearSeqs'.\.z Main 115 192024 0.3 0.0 0.3 0.0
nearSeqs.nearSeqs'.\.(...) Main 114 193040 4.3 12.3 4.3 12.3
nearSeqs.nearSeqs'.\.z' Main 113 183951 0.8 0.0 0.8 0.0
CAF:$dShow_r1F9 Main 99 0 0.0 0.0 0.0 0.0
CAF:$dOrd_r1F8 Main 98 0 0.0 0.0 0.0 0.0
CAF:$dEq_r1D9 Main 97 0 0.0 0.0 0.0 0.0
CAF:alphabet_ryE Main 96 0 0.0 0.0 0.0 0.0
alphabet Main 111 1 0.0 0.0 0.0 0.0
CAF GHC.Conc.Signal 92 0 0.0 0.0 0.0 0.0
CAF GHC.IO.Handle.FD 84 0 0.0 0.0 0.0 0.0
CAF GHC.IO.Encoding 79 0 0.0 0.0 0.0 0.0
CAF GHC.IO.Encoding.Iconv 68 0 0.0 0.0 0.0 0.0
Видим, что абсолютным чемпионом по времязатратам является функция grpSeqs с результатом 52% от всего времени выполнения программы, за ней идет функция nearSeqs — 20.6%. Функция grpSeqs состоит сплошь из сортировок и группировок списков — эти типы вычислений не имеют значительных перспектив для распараллеливания. С другой стороны, списки, возвращаемые функцией nearSeqs для 254 уникальных подпоследовательностей исходной последовательности seqs !! 1, могут быть вычислены независимо, а это может означать потенциально ожидаемый 15% выигрыш по времени для 4-ядерной машины (20.6% - (20.6% / 4) ≈ 15%).
Перепишем функцию grpSeqs с параллельным вычислением nearSeqs. В список импортируемых модулей добавляем
import Control.Parallel.Stategies
Определяем новую функцию grpSeqsPar.
grpSeqsPar :: Int -> [String] -> [[String]]
grpSeqsPar n = L.sortBy (flip $ comparing length) . L.group . L.sort . seqs
where seqs = concat . parMap rdeepseq
(\(s, m) -> concat . replicate m . nearSeqs n $ s) .
map (head &&& length) . L.group . L.sort
Не забываем заменить вызов grpSeqs на grpSeqsPar в функции main.
mapM_ print' $ take 10 $ grpSeqsPar 2 $ subSeqs 10 $ seqs !! 1
Компилируем с поддержкой распараллеливания.
ghc --make -O2 -rtsopts -threaded -feager-blackholing -fforce-recomp gen
На моей машине установлен 4-ядерный процессор Intel Core i5. Запускаем с распараллеливанием на 4 ядра.
time ./gen +RTS -N4
("GCACACAGAC",20)
("GCGCACACAC",20)
("ACACACACAC",19)
("CCCGCACACA",19)
("CACACACACG",18)
("CGCGCACACA",18)
("CGGCACACAC",18)
("GCACACACCC",18)
("GCACACACGC",18)
("GCCCGCACAC",18)
real 0m0.521s
user 0m1.049s
sys 0m0.339s
Что-то никакого толку. Давайте запустим программу с опцией -s, которая выводит разную полезную информацию на экран.
time ./gen +RTS -N4 -s
("GCACACAGAC",20)
("GCGCACACAC",20)
("ACACACACAC",19)
("CCCGCACACA",19)
("CACACACACG",18)
("CGCGCACACA",18)
("CGGCACACAC",18)
("GCACACACCC",18)
("GCACACACGC",18)
("GCCCGCACAC",18)
316,334,968 bytes allocated in the heap
161,177,512 bytes copied during GC
18,120,968 bytes maximum residency (10 sample(s))
310,528 bytes maximum slop
53 MB total memory in use (0 MB lost due to fragmentation)
Tot time (elapsed) Avg pause Max pause
Gen 0 365 colls, 365 par 0.54s 0.14s 0.0004s 0.0147s
Gen 1 10 colls, 9 par 0.30s 0.08s 0.0084s 0.0159s
Parallel GC work balance: 38.30% (serial 0%, perfect 100%)
TASKS: 10 (1 bound, 9 peak workers (9 total), using -N4)
SPARKS: 254 (254 converted, 0 overflowed, 0 dud, 0 GC'd, 0 fizzled)
INIT time 0.00s ( 0.00s elapsed)
MUT time 0.60s ( 0.30s elapsed)
GC time 0.84s ( 0.23s elapsed)
EXIT time 0.00s ( 0.00s elapsed)
Total time 1.44s ( 0.53s elapsed)
Alloc rate 529,187,979 bytes per MUT second
Productivity 41.9% of total user, 113.3% of total elapsed
gc_alloc_block_sync: 67687
whitehole_spin: 0
gen[0].sync: 22
gen[1].sync: 10747
real 0m0.536s
user 0m1.025s
sys 0m0.417s
Продуктивность программы составила 41.9%, остальное время ушло на сборку мусора (GC time 0.84s)! Обычно с этим борются увеличением предполагаемого размера кучи (suggested heap size) для сборщика мусора (см. документацию здесь). Давайте установим эту величину в 1 гигабайт.
time ./gen +RTS -N4 -H1G -s
("GCACACAGAC",20)
("GCGCACACAC",20)
("ACACACACAC",19)
("CCCGCACACA",19)
("CACACACACG",18)
("CGCGCACACA",18)
("CGGCACACAC",18)
("GCACACACCC",18)
("GCACACACGC",18)
("GCCCGCACAC",18)
315,385,264 bytes allocated in the heap
17,826,760 bytes copied during GC
15,081,784 bytes maximum residency (2 sample(s))
118,472 bytes maximum slop
970 MB total memory in use (0 MB lost due to fragmentation)
Tot time (elapsed) Avg pause Max pause
Gen 0 1 colls, 1 par 0.03s 0.01s 0.0110s 0.0110s
Gen 1 2 colls, 1 par 0.12s 0.04s 0.0178s 0.0321s
Parallel GC work balance: 83.58% (serial 0%, perfect 100%)
TASKS: 10 (1 bound, 9 peak workers (9 total), using -N4)
SPARKS: 254 (254 converted, 0 overflowed, 0 dud, 0 GC'd, 0 fizzled)
INIT time 0.00s ( 0.00s elapsed)
MUT time 0.56s ( 0.31s elapsed)
GC time 0.15s ( 0.05s elapsed)
EXIT time 0.00s ( 0.00s elapsed)
Total time 0.72s ( 0.37s elapsed)
Alloc rate 567,298,189 bytes per MUT second
Productivity 78.5% of total user, 153.2% of total elapsed
gc_alloc_block_sync: 15012
whitehole_spin: 0
gen[0].sync: 2831
gen[1].sync: 0
real 0m0.391s
user 0m0.600s
sys 0m0.137s
А вот это уже здорово! Время выполнения программы сократилось до 0.39 секунды, продуктивность увеличилась до 78.5%. Кстати, в выводе программы есть строка, озаглавленная SPARKS: в ней находится статистика запусков параллельных вычислений, которые называются искрами (sparks). Всего было запущено 254 искры, что соответствует числу уникальных 10-меров в исходной последовательности.
А что будет, если в оригинальной, нераспараллеленной версии выделить 1 гигабайт для кучи? Может быть она тоже ускорится? Как показала практика — нет. Продуктивность программы действительно возрастает, но время выполнения практически не изменяется, как это ни странно.
И, наконец, следующая версия функции grpSeqsPar работает еще чуточку быстрее (но незначительно).
grpSeqsPar :: Int -> Int -> [String] -> [[String]]
grpSeqsPar m n = L.sortBy (flip $ comparing length) . L.group . sortPar . seqs
where sortPar = foldr LO.merge [] . parMap rdeepseq L.sort
seqs = concat' m . parMap rdeepseq
(\(s, m) -> concat . replicate m . nearSeqs n $ s) .
map (head &&& length) . L.group . L.sort
concat' m s = let l = length s
in map concat $ chunksOf
(l `div` m + if l `mod` m == 0 then 0 else 1) s
Вот только я совсем не уверен в ее оптимальности (именно в оптимальности, а не в корректности): в самом деле, вложенные стратегии rdeepseq скорее всего не самый лучший вариант. В этой версии затратная сортировка списка, возвращаемого функцией seqs, заменена на серию последовательных слияний LO.merge его по отдельности отсортированных частей. Количество частей, на которые разбивается исходный список, задается новым параметром m функции grpSeqsPar. На мой взгляд, он должен быть равен количеству ядер процессора. Соответственно, вызов в функции main принимает вид
mapM_ print' $ take 10 $ grpSeqsPar 4 2 $ subSeqs 10 $ seqs !! 1
Для использования функции chunksOf в список импортируемых модулей следует добавить
import Data.List.Split (chunksOf)
Если скомпилировать новую программу и запустить ее на выполнение, то в строке SPARKS появятся несколько десятков погасших искр (fizzled sparks). Искры гаснут в том случае, когда танки (thunks), то есть элементы для вычисления, на которые они ссылаются, уже были вычислены в некотором параллельном потоке программы. Такая ситуация может свидетельствовать о неоптимальном выборе стратегии параллелизации. Опция компиляции -feager-blackholing позволяет предотвратить многократные вычисления одних и тех же танков, поэтому само по себе появление погасших искр не должно влиять на производительность программы.
Что еще можно распараллелить в этой программе? Ну, например, функцию subSeqs. Хотя в данном случае это не даст никаких преимуществ, поскольку время выполнения этой функции по результатам профилирования составляет 0.0%. Тем не менее, для полноты картины я приведу параллельную версию subSeqs.
subSeqsPar :: Int -> Int -> String -> [String]
subSeqsPar m n s = subSeqs n s `using` parListChunk (length s `div` m) rdeepseq
Здесь я решил применить parListChunk со стратегией rdeepseq. Размер чанка вычисляется на основании значения первого параметра m функции subSeqsPar путем разделения исходной последовательности s на m равных частей: это упрощенный аналог алгоритма из функции concat’.
Update. Вложенная стратегия parMap rdeepseq для вычисления всех мутаций с помощью nearSeqs в функции grpSeqsPar действительно не нужна: внешняя стратегия из функции sortPar должна прекрасно справиться с распараллеливанием сама. По этой же самой причине subSeqsPar тоже не нужна в принципе. Таким образом, функция grpSeqsPar переписывается так:
grpSeqsPar :: Int -> Int -> [String] -> [[String]]
grpSeqsPar m n = L.sortBy (flip $ comparing length) . L.group . sortPar . seqs
where sortPar = foldr LO.merge [] . parMap rdeepseq L.sort
seqs = concat' m .
map (\(s, m) -> concat . replicate m . nearSeqs n $ s) .
map (head &&& length) . L.group . L.sort
concat' m s = let l = length s
in map concat $ chunksOf
(l `div` m + if l `mod` m == 0 then 0 else 1) s
Новая стратегия должна сократить число искр до четырех и потенциально ускорить выполнение программы.
Кстати, для визуализации выполнения параллельной программы на haskell существует прекраснейший инструмент — программа threadscope. Для того, чтобы воспользоваться threadscope, нужно скомпилировать нашу программу с поддержкой журнала событий.
ghc --make -O2 -rtsopts -threaded -eventlog -feager-blackholing -fforce-recomp gen
Давайте запустим программу с опцией -H1G, как мы это делали раньше.
time ./gen +RTS -N4 -H1G -qa -lf
("GCACACAGAC",20)
("GCGCACACAC",20)
("ACACACACAC",19)
("CCCGCACACA",19)
("CACACACACG",18)
("CGCGCACACA",18)
("CGGCACACAC",18)
("GCACACACCC",18)
("GCACACACGC",18)
("GCCCGCACAC",18)
real 0m0.282s
user 0m0.525s
sys 0m0.118s
Не обращайте внимание на новый рекорд скорости (0.282 сек.): просто я запустил программу на более быстрой машине. Зато обратите внимание на новые опции: -qa и -lf. Первая опция экспериментальная и в документации haskell говорится, что она может ускорить выполнение параллельной программы, а может и не ускорить. Для меня она действительно работает. Вторая опция приводит к записи событий в файл gen.eventlog. Давайте откроем этот файл в программе threadscope.
Вот такая картинка появится на экране. Что ж, здесь достаточно много информации как на графике в центральном окне, так и во вкладках, расположенных в нижнем окне. Верхняя часть графика — Activity — отражает степень загрузки потоков выполнения (у нас их четыре) в течении всего времени жизни программы. Четыре плотных зеленых слоя, начинающиеся на 25 мс означают, что все четыре потока работают на полную катушку: это очень хорошо. Начиная с 0.13 с программа стала работать в одном потоке: это финальная сортировка по длинам в функции grpSeqsPar (вы можете убедиться в этом, удалив L.sortBy (flip $ comparing length) . из определения grpSeqsPar и запустив программу заново).
А что же это за провал, длившийся 20 мс в самом начале программы, когда активность была нулевой и все потоки занимались сборкой мусора (оранжевые полосы)? К счастью, график является интерактивным, и мы можем узнать, что произошло в определенной точке выполнения программы, кликнув мышкой по графику в этой точке и просмотрев события во вкладке Raw events.
Ага, переполнение кучи. Опция -H1G не влияет на первый запуск сборщика мусора: для этого нужно настроить опцию -A. Давайте запустим программу с этой опцией.
time ./gen +RTS -N4 -A128M -qa -lf
("GCACACAGAC",20)
("GCGCACACAC",20)
("ACACACACAC",19)
("CCCGCACACA",19)
("CACACACACG",18)
("CGCGCACACA",18)
("CGGCACACAC",18)
("GCACACACCC",18)
("GCACACACGC",18)
("GCCCGCACAC",18)
real 0m0.270s
user 0m0.471s
sys 0m0.115s
Я подобрал размер аллокации (128M) таким образом, чтобы единожды выделенной кучи хватило на все время выполнения программы. Теперь threadscope покажет такую картинку.
На этот раз продуктивность программы составила 99% против 93% в предыдущем запуске, сборщик мусора практически не задействован. Однако время выполнения почти не изменилось. Помните, я говорил, что это странно. Так вот, никакой странности здесь оказывается нет: просто время сборки мусора эффективно заменилось временем, затраченным на аллокацию и управление гораздо бо́льшими объемами памяти. В данном случае значительное время было потрачено в фазе инициализации (INIT): это хорошо видно на графике.
А что будет, если использовать параметры аллокации кучи для сборщика мусора по умолчанию?
time ./gen +RTS -N4 -qa -lf
("GCACACAGAC",20)
("GCGCACACAC",20)
("ACACACACAC",19)
("CCCGCACACA",19)
("CACACACACG",18)
("CGCGCACACA",18)
("CGGCACACAC",18)
("GCACACACCC",18)
("GCACACACGC",18)
("GCCCGCACAC",18)
real 0m0.312s
user 0m0.877s
sys 0m0.169s
Интересный результат. Время выполнения немного выросло. Время сборщика мусора (GC time) теперь равно времени полезной нагрузки (Mutator time), и это видно по соотношению зеленого и оранжевого цветов в потоках. Загрузка просела до двух (чуть позже трех) потоков одновременно. Несмотря на это, время выполнения, как я уже сказал, выросло ненамного: видимо это связано с более эффективной стратегией аллокации и меньшим потреблением памяти.
Update 2. Функцию concat’ из grpSeqsPar можно переписать с помощью стандартной функции ceiling.
concat' m s = let l = let (/.) = (/) `on` fromIntegral
in ceiling $ length s /. m
in map concat $ chunksOf l s
Для функции on нужен
import Data.Function (on)
Определение функции (/.) я взял отсюда. Она нужна для преобразования типов (или коэрции) аргументов операции вещественного деления (/) из Int в абстрактный тип Num a => a. Без этого преобразования типов программа не скомпилируется: неявные преобразования типов в стиле C и C++ в haskell не поддерживаются.
Комментариев нет:
Отправить комментарий