Comments 28
Теперь нужно сравнить с DGEMM из какой-нибудь хорошей реализации BLAS на CPU (ATLAS, MKL) или на GPU (ViennaCL) и поставить жирную точку в вопросе, «а кто ж быстрее-то?».
или сравнить с haskell реализацией.
Всё же сравнивать с нативными отклаженными годами библиотеками будет не совсем корректно. Переписывать всё на С/С++ желания и времени, к сажалению, нет(
НО! Раз наша статья получила такой резонанс(Мы ожидали меньшего), то появилась идея обощить принципы решения подобных типов задач и, например, взять к и попробовать всё это дело на СКЦ(Суперкомпьютерном центре) Политехнического Университета СПб. Дело тут уже состоит в том, чтобы убедить высокопоставленных людей, что подобного рода задач интересны и важны, тогда получим доступ к очень серьёзным ресурсам. Надеюсь я Вы меня правильно поняли, большое спасибо за комментарий к нашей работе! :)
НО! Раз наша статья получила такой резонанс(Мы ожидали меньшего), то появилась идея обощить принципы решения подобных типов задач и, например, взять к и попробовать всё это дело на СКЦ(Суперкомпьютерном центре) Политехнического Университета СПб. Дело тут уже состоит в том, чтобы убедить высокопоставленных людей, что подобного рода задач интересны и важны, тогда получим доступ к очень серьёзным ресурсам. Надеюсь я Вы меня правильно поняли, большое спасибо за комментарий к нашей работе! :)
к сожалению
Странно, что нет функционала по редактированию своих же комментариев( много опечаток
Странно, что нет функционала по редактированию своих же комментариев( много опечаток
Как раз сравнить с нативной библиотекой — более чем корректно. Если для решения типовой задачи, для которой есть прекрасная (и с «хорошей» лицензией BSD) отмычка в виде ATLAS или ViennaCL, начинают лепить велосипед, нужно убедиться, что проигрыш будет хотя бы не драматическим.
Причем, написание велосипеда — это проигрыш на время самого написания и отладки (что быстрее, написать ваше или подцепить BLAS через JNI) и на скорость вычислений (ведь и на С++ можно написать так, что все колом встанет).
Причем, написание велосипеда — это проигрыш на время самого написания и отладки (что быстрее, написать ваше или подцепить BLAS через JNI) и на скорость вычислений (ведь и на С++ можно написать так, что все колом встанет).
1) в отличие от старой статьи не указан процессор — количесто ядер? :)
2) почему при удвоении задачи (скачек степени двойки, как я понимаю?) время растет в 7 раз?
7 умножений, но непонятна функция роста.
2) почему при удвоении задачи (скачек степени двойки, как я понимаю?) время растет в 7 раз?
7 умножений, но непонятна функция роста.
1. Тестировалось все на обычных ноутбуках c i5-i7 и 6gb оперативной памяти(ноутам 3-4 года).
2. Да вы всё правильно понимаете, дело в том, что да, дело в степени двойки:)… Алгоритм работает с квадратными матрицами и из сего следует, что при размерах в 1025*1025 матрица увеличивается до 2048*2048, так же если размер матрицы 1500*1500, то матрица всё также увеличивается до размерности в 2048. Вот и происходит 7 кратный взрыв при переходе с одной степени 2 к другой. И в общем случае получается, что переумножать матрицы размером, например, в 1200*1200 и 2000*2000 — нет, просто нет.(при расширении она дополняется нулями)
Спасибо за комментарий и обратную связь!
2. Да вы всё правильно понимаете, дело в том, что да, дело в степени двойки:)… Алгоритм работает с квадратными матрицами и из сего следует, что при размерах в 1025*1025 матрица увеличивается до 2048*2048, так же если размер матрицы 1500*1500, то матрица всё также увеличивается до размерности в 2048. Вот и происходит 7 кратный взрыв при переходе с одной степени 2 к другой. И в общем случае получается, что переумножать матрицы размером, например, в 1200*1200 и 2000*2000 — нет, просто нет.(при расширении она дополняется нулями)
Спасибо за комментарий и обратную связь!
А что если дать каждому алгоритму перемножать N наборов матриц, где N кратно числу ядер на машине? Тогда можно будет сделать и сравнить параллельную версию Usual Algorithm (Transposed) с Stassen Algorithm Fork-Join. А вдруг он всё же лучше?
Что Вы понимаете под набором матриц?
В каждом тесте работа велась только с одной матрицей определенной размерности. Дальше матрицы с шагом в 100 увеличивали свою размерность. Алгоритм Штрассена выделяет подматрицы и каждую такую подматрицу можно и нужно расспараллеливать на вычисление. А вот как выделить подматрицы в обычном алгоритме и можно ли это в принципе… я лично не знаю(отношусь скептически):)
В каждом тесте работа велась только с одной матрицей определенной размерности. Дальше матрицы с шагом в 100 увеличивали свою размерность. Алгоритм Штрассена выделяет подматрицы и каждую такую подматрицу можно и нужно расспараллеливать на вычисление. А вот как выделить подматрицы в обычном алгоритме и можно ли это в принципе… я лично не знаю(отношусь скептически):)
Я понял комментарий smile616 как то, что алгоритм Штрассена будучи выполняемым на 8 ядрах начинает опережать вариант с транспонированием выполняющийся на 1 ядре. Но разница далеко не приближается к 8. И что будет интересно замерять время выполнения двух алгоритмов перемножая матрицы не один раз, а как минимум 8 раз. Так, чтобы алгоритм транспонирования был запущен в 8-ми тредах, обрабатывая каждый свою матрицу.
Прочитал метод сборки матриц и вспомнил анекдот про индейцев:
В вашей реализации идут постоянные копирования, из-за этого сильно проседает производительность.
ПС: Диаграмма обрывается прямо перед скачком, что не всем очевидно.
А на третий день Зоркий Глаз заметил, что четвёртой стены нет.Насчёт транспонирования, это позволяет использовать кэш процессора намного эффективнее.
В вашей реализации идут постоянные копирования, из-за этого сильно проседает производительность.
ПС: Диаграмма обрывается прямо перед скачком, что не всем очевидно.
Траспонирование — это не какой-то magic. Оно нужно сугубо для Java, это её особенность чтения массивов. На других языках (например: Fortran) вообще не потребуется никого траспонировать.
Предложите вариант без копирования. Вообще нативное копирование работает очень даже быстро. Тому пример ArrayList, который даже с копирование в среднем случае работает быстрей LinkedList'а.
Скачки видны до этого, а если бы мы их добавили, то сложно было бы смаштабировать результаты в читабельный график (Да и памяти у нас под Штрассена не много).
Предложите вариант без копирования. Вообще нативное копирование работает очень даже быстро. Тому пример ArrayList, который даже с копирование в среднем случае работает быстрей LinkedList'а.
Скачки видны до этого, а если бы мы их добавили, то сложно было бы смаштабировать результаты в читабельный график (Да и памяти у нас под Штрассена не много).
Траспонирование — это не какой-то magic. Оно нужно сугубо для Java, это её особенность чтения массивов. На других языках (например: Fortran) вообще не потребуется никого траспонировать.
Deosis правильно написал, что транспонирование позволяет эффективнее использовать кэш, а вы говорите про какой-то «magic». Это не особенность Java, это особенность хранения двумерного массива в линейном адрессном пространстве. И в Fortran будет то-же самое, с той лишь разницей, что транспонирование там нужно будет делать иначе, так как двумерные массивы там хранятся column-wise, а не row-wise как в большинстве языков, где есть двумерные массивы.
Примечательно, что Fortran не имеет такой проблемы. Многомерные массивы в нем хранятся по столбцам. Поэтому, кеш будет эксплуатироваться как следует и штатная реализация кубического алгоритма будет работать в разы быстрее.
В нашей статье есть ссылка на другую хабра-статью с траспонированием и там сказано, что в Фортране не потребуется делать траспонирование и обычный алгоритм будет работать в разы быстрей и кэшем всё тоже ок.
Спасибо за Ваш интерес!)
P.S. «magic» я упомянул с иронией, не более
То, что там так сказано, не значит, что так оно и есть. В той статье есть комментарий по этому поводу.
Дело в том, что column-wise ни чем не лучше row-wise, это абсолютно симметричные способы хранения матриц. Поэтому column-wise не может быть лучше или хуже для перемножения.
Для нематематических приложений, row-wise удобнее, для математических column-wise удобнее.
Мысль, которую я пытался донести, это то, что это никакая не «особенность Java». И в Fortran проблема будет точно такая же. Не во всех ЯП есть многомерные массивы, вот в С/C++ нет. Но сам факт того, что двумерный массив хранится в линейном адрессном пространстве, вызывает cache miss при перемножении, и неважно как мы выбрали хранить массивы column-wise или row-wise.
Посмотрите на тело самого внутреннего цикла:
Здесь, внутренний цикл инкрементирует переменную k. Доступ к элементам матрицы A не вызывает cache miss, так как A[i][k] и A[i][k+1] находятся «рядом», но есть проблема с B, так как элементы B[k][j] и B[k+1][j] находятся «далеко». Поэтому трансонирование матрицы B так сильно улучшает результат.
Теперь представим, что у нас Fortran, и матрицы хранятся column-wise. В этом случае элементы B[k][j] и B[k+1][j] находятся «рядом» и матрицу B транспонировать не нужно, однако A[i][k] и A[i][k+1] находятся «далеко» и теперь уже доступ к ним будет вызывать cache miss. Теперь транспонировать нужно матрицу A.
Что мы выиграли заменив row-wise на column-wise? Ничего.
P.S.
Дело в том, что column-wise ни чем не лучше row-wise, это абсолютно симметричные способы хранения матриц. Поэтому column-wise не может быть лучше или хуже для перемножения.
Для нематематических приложений, row-wise удобнее, для математических column-wise удобнее.
Мысль, которую я пытался донести, это то, что это никакая не «особенность Java». И в Fortran проблема будет точно такая же. Не во всех ЯП есть многомерные массивы, вот в С/C++ нет. Но сам факт того, что двумерный массив хранится в линейном адрессном пространстве, вызывает cache miss при перемножении, и неважно как мы выбрали хранить массивы column-wise или row-wise.
Посмотрите на тело самого внутреннего цикла:
C[i][j] += A[i][k] * B[k][j];
Здесь, внутренний цикл инкрементирует переменную k. Доступ к элементам матрицы A не вызывает cache miss, так как A[i][k] и A[i][k+1] находятся «рядом», но есть проблема с B, так как элементы B[k][j] и B[k+1][j] находятся «далеко». Поэтому трансонирование матрицы B так сильно улучшает результат.
Теперь представим, что у нас Fortran, и матрицы хранятся column-wise. В этом случае элементы B[k][j] и B[k+1][j] находятся «рядом» и матрицу B транспонировать не нужно, однако A[i][k] и A[i][k+1] находятся «далеко» и теперь уже доступ к ним будет вызывать cache miss. Теперь транспонировать нужно матрицу A.
Что мы выиграли заменив row-wise на column-wise? Ничего.
P.S.
P.S. «magic» я упомянул с иронией, не болееЯ понял вашу иронию, cache miss для меня никак не «magic», а объективная реальность, а вот загадочная «особенность Java» о которой вы говорили, как раз «magic».
Я смотрю, транспонированный алгоритм имеет скорость роста 7,5-8 за удвоение. Но маленький коэффициент перед х^3
Штрасена — ровно 7. Так что в перспективе лучше даже без распаралеливанья. Но с какого момента гарантировано лучше?
Штрасена — ровно 7. Так что в перспективе лучше даже без распаралеливанья. Но с какого момента гарантировано лучше?
Гарантий я дать не могу. Но думаю, что асимптотически расспараллеленный Штрассен будет работать быстрей с точки зрения математики. На практике же существуют ограничения, например, из-за своей рекурсивной природы алгоритму требуется много памяти, и на предельных загрузках многое будет зависить от GC (Мы же всё таки говорим о реализации под JVM), да и число ядер у нас не бесконечное число т.д. и т.п… Поэтому не нужно считать это панацеей. В нашей случае, мы использовали средненькие по производительности ноутбуки с выделением 8 ядер и оперативной памяти порядка 4-6гб под jvm. GC работал на пределе. И ясно, что проверить реализации у нас не получится на матрицах размерностью 8_000 и более в силу ограничений на железо.
log2
Integer.highestOneBit(x - 1) + 1?
log2(Collections.max(Arrays.asList(a.length, a[0].length, b[0].length)))
Integer.highestOneBit((a.length-1) | (a[0].length-1) | (b[0].length-1)) + 1 — вроде так. Битовое или как раз вам сохранит старший установленный бит. Раз бьёмся за производительность, странно создавать список, объектный массив и боксинг-анбоксинг впихивать. Это один раз делается, конечно, но всё же.
Копирований, конечно, сильно много уж. Обычно подобные алгоритмы реализуют, выделяя кусок рабочей памяти один раз и используя из него соответствующие фрагменты. Тут даже в простой реализации, видимо, дополнительной памяти потребуется столько же, сколько исходная матрица. Так как у вас всё в степенях двойки, вычисление смещений будет быстро делаться с помощью битовых сдвигов. Ну и появится локальность, конечно.
Вы правы. Так действительно будет быстрей, но на итоговую скорость работы всего алгоритма вряд ли повлияет. Зато код станет кратче.
Что касается копирования, да их тут много. Подобные алгоритмы в принципе решаются на других языка, где имеется соответсвующая гибкость работы с памятью. Мы же говорим про реализацию только на Java.
Спасибо за Ваши советы!
Что касается копирования, да их тут много. Подобные алгоритмы в принципе решаются на других языка, где имеется соответсвующая гибкость работы с памятью. Мы же говорим про реализацию только на Java.
Спасибо за Ваши советы!
Мы же говорим про реализацию только на Java
Ну это вы зря. Из того, что реализация только на Java, не значит, что её нельзя сделать быстрой. В данной задаче не нужна гибкость управления памятью. Языки вроде C/C++ могут выиграть здесь за счёт применения эффективных инструкций SSE (по факту реализация будет даже не на C/C++, а на ассемблере). Но особая гибкость памяти здесь не нужна, джавовское управление памятью никак не замедляет алгоритм, если правильно им воспользоваться.
Из того, что реализация только на Java, не значит, что её нельзя сделать быстрой
А кто сказал, что реализация работает медленно?) Просто есть противоречия в словах. То память расходуется не правильно, то в джавовское управление памятью не замедляет алгоритм. Если есть конкретные идеи, то предлагайте. А то это перерастает в какой-то холивар)
Эээ. Конкретная идея была в предыдущем комментарии. Может, я её непонятно сформулировал. Надо выделить буфер один раз и в нём жить, не выделять память много раз и не копировать. Кроме того выделить линейный массив, а не многомерный.
Я, наверно, понял Вашу идею. Вы предлагаете в task'и передавать исходную матрицу и еще необходимые индиксы для того, чтобы выделять в исходной необходимые подматрицы. Хммм, да, скорее всего при грамотной реализации это даст существенный прирост производительности, но вот читаймость и ясность кода упадет в разы.( Тоже касается и с выделением линейного массива. Думаю, лучше выбор решения оставить за тем, кому оно понадобится. Всё таки приросты производительности будут на очень больших матрицах(На практике почти не встречаются). А для матриц меньше 500*500 лучше оставить код более легкий и краткий код.
Спасибо за Ваш интерес.
P.S. кст, Вы один из тех, из-за кого наш преподаватель отменил занятие(из-за Joker), так что, возможно, если бы не Вы, то этой статьи и не было бы:))
Спасибо за Ваш интерес.
P.S. кст, Вы один из тех, из-за кого наш преподаватель отменил занятие(из-за Joker), так что, возможно, если бы не Вы, то этой статьи и не было бы:))
У вас в простейшем кубическом методе перемножения ошибка.
При перемножении матрицы A размерностью mxn и матрицы B размерностью nxk, должна получаться матрица C размерностью mxk. Ваш алгоритм будет работать корректно только для квадратных матриц.
Ваш код:
Если A.columns() > A.rows(), то будет выход за пределы индексов С (если конечно вы создали матрицу С правильно размера). И в самом внутреннем цикле, граница также неправильная.
По идее, должно быть что-то вроде:
Ну и конечно где-то выше должна быть проверка, что A.columns() == B.rows(), но положим, что проверку мы здесь опускаем.
P.S. Я понимаю, что это было взято из другой статьи, но там хватает неточностей.
При перемножении матрицы A размерностью mxn и матрицы B размерностью nxk, должна получаться матрица C размерностью mxk. Ваш алгоритм будет работать корректно только для квадратных матриц.
Ваш код:
for (int i = 0; i < A.rows(); i++) {
for (int j = 0; j < A.columns(); j++) {
for (int k = 0; k < B.columns(); k++) {
C[i][j] += A[i][k] * B[k][j];
}
}
}
Если A.columns() > A.rows(), то будет выход за пределы индексов С (если конечно вы создали матрицу С правильно размера). И в самом внутреннем цикле, граница также неправильная.
По идее, должно быть что-то вроде:
for (int i = 0; i < A.rows(); i++) {
for (int j = 0; j < B.columns(); j++) {
for (int k = 0; k < A.columns(); k++) {
C[i][j] += A[i][k] * B[k][j];
}
}
}
Ну и конечно где-то выше должна быть проверка, что A.columns() == B.rows(), но положим, что проверку мы здесь опускаем.
P.S. Я понимаю, что это было взято из другой статьи, но там хватает неточностей.
Sign up to leave a comment.
Divide&Conquer над алгоритмом Штрассена