SGEMM на видеокарте и CPU, серия 6

Уже становится какой-то дурной традицией переделывать тест умножения матриц каждые несколько месяцев. Однако с прошлого забега случилось несколько событий:

  • В форумах по CUDA появилась очень быстрая реализация SGEMM, сделанная Василием Волковым.
  • Вышел CUDA Toolkit 2.0 beta, использующий реализацию Волкова в библиотеке CuBLAS.
  • 4-ядерный Penryn стал не просто доступен, а появился у меня прямо под столом.

Участники забега

Hardware

  • Intel Core2Quad 9300, на частоте 3GHz (работает отлично, температура 49-59 градусов), память DDR2 PC6400 4-4-4-12
  • NVidia GeForce 8800GTX на штатной частоте
Software
  • Windows XP Professional x64
  • Intel MKL 10.0.3 (использовалась совместно с MS Visual C 2005 в 32-битном режиме)
  • Goto BLAS 1.26 (использовалась под Cygwin, 32-битный режим)
  • NVidia CUDA 2.0 beta (32-битная версия)

Методика тестирования

Методика совершенно незатейлива: берем квадратные матрицы, размер которых равен степени двойки (я знаю, что это дает наибольшие цифры), перемножаем, определяем ушедшее на умножение время, считаем гигафлопсы. Содержимое матриц - случайное. Перед забегом проверяем правильность работы, перемножая маленькие матрицы (8x8) и проверяя глазом корректность результата.

Для видеокарты жизнь чуть сложнее, там есть время потраченное на вычисления и дополнительное время, нужное для загрузки-выгрузки данных, поэтому для CUDA приводятся два значения производительности.

Результаты

Производительность умножения матриц (SGEMM), GFLOP/s
Описание теста размер задачи
  1024 2048 4096
Intel MKL 10.0.3 54.8 59.5 60.5
Goto BLAS 1.26 76.7 83.4 87.1
NVidia CUDA 2.0 beta, чистое время вычислений 202 204 205
NVidia CUDA 2.0 beta, вычисления+ввод-вывод 151 177 190

Обсуждение результатов

Сначала проговорим сами результаты:

  1. Penryn - феноменально быстрый. Надо смотреть не только на абсолютные цифры производительности (которые тоже большие), но и в количество реальных операций на такт. Предельная (теоретическая) производительность у этого процессора - 8 операций с одинарной точностью на ядро * 4 ядра * 3Ghz = 96 GFLOP/s, а мы получили 87 или больше 90%. И это все - за $350
  2. Несмотря на анонсированные в 10-й версии Intel MKL ускорения, Goto BLAS так и остается раза в полтора быстрее, чем библиотека от производителя.
  3. CUDA 2.0 ускорилась относительно первой версии в 1.7 раза на вычислениях и в 1.9 раза на вычислениях+I/O (вклад в скорость I/O могла дать и замена hardware хост-компьютера, хотя по тестам bandwidth я принципиального ускорения не вижу).

Удивительно, но произошло следущее

  • Стоимость гигафлопса на CPU упала с октября прошлого года примерно втрое (для десктопных процессоров) и достигла стоимости гигафлопса на GPU полугодовалой давности.
  • При этом, стоимость гигафлопса (только по одному тесту) на GPU тоже упала почти вдвое при том, что сам GPU остался прежним и подешевел несильно. Падение стоимости исключительно софтовое.
  • Конечно, 8800GTX - далеко не лидер по price-performance сейчас, новые карты могут на некоторых задачах и обгонять (но надо проверять).
  • Впрочем, соотношение в очередной раз изменится буквально в ближайшие недели, ходят упорные слухи, что новый чип от NVidia всех опять порвет.

Comments

Я думал, что Core2Quad способен выполнять только 4 операции за такт. Где упоминается 8 операций? Можно ссылочку?

А это просто креативный подсчет, такой же как у производителей видеокарт.

SSE делает 4 32-битных MAD за такт. А MAD считается за две операции.

свежий результат: 980 Gflop/s в умножении матриц на ATI HD4870, http://forum.beyond3d.com/showthread.php?t=54842

Все-таки 880 (980 - это теоретический расчет), но безумно круто в любом случае.

"This is 11% faster in theory. It achieves 980 Gflop/s in practice" - http://forum.beyond3d.com/showpost.php?p=1325769&postcount=27

Во понаписали, и не прочесть за раз.

OpenCL для ATI будет? Или может уже есть, просто я пропустил?

Вроде обещают в третьем квартале - "The GPU version of OpenCL to be bundled with the ATI Stream SDK is scheduled to be available later in Q3 of 2009", http://developer.amd.com/gpu/ATIStreamSDKBetaProgram/assets/ATI_Stream_S...

А между тем ATI Stream SDK 2.0 Beta2:
* First public beta release of ATI Stream SDK with OpenCL™ CPU support.
* The OpenCL™ 1.0 conformance logs from this release have been submitted to the Khronos OpenCL™ Working Group.

http://developer.amd.com/GPU/ATISTREAMSDKBETAPROGRAM/Pages/default.aspx

Надо, пожалуй, купить ATI-шную карточку, ибо есть тут меня идея-фикус. И на маке уже (вроде) все есть или с минуты на минуту будет....

по моему опыту, SGEMM в Intel MKL (какая-то из 10-х версий) работает гораздо быстрее в 64-битном режиме. Разница примерно такая: 24 Gflop/s при запуске в 32-битных виндах и 38 Gflop/s в 64-битных на 2.67 GHz Core2 Duo E6700. То есть 64-битный Intel MKL тоже может работать под 90% пика, как и Гото. Наблюдал аналогичную картину и на четырёх ядерных процессорах.

Одно объяснений этому явлению что я слышал --- в 2 раза больший размер доступного файла регистров в 64-битном режиме по сравнению с 32-битным. Соответственно, можно использовать большие блоки в регистрах чтобы разрешить нехватку скорости кеша.

Ну вот буду (на праздниках?) бенчмаркать 280-GTX, заодно и это проверю, спасибо за идею.

Не плохо бы еще добавить результаты для ATI карт.
К сожалению, новых 48*0 карт у меня нет, но вот результат, полученный мною 3 года назад, для старенькой ATI1950XT, с использованием CTM (sgemm):
1024x1024 - 77.5 Gflops
2048x2048 - 80.4 Gflops

Интересно, сколько может выдать на гора ATI4870 ?

На X1900XT в умножении матриц получалось до 120 Gflop/s. В прошлогодней публикации от AMD цитируется 200 Gflop/s на AMD HD 2900. На AMD FireStream 9250 заявлено 174 Gflop/s в FFT.

А разве три года назад был CTM ?

Но это придирки, конечно. Давайте результаты свежих карт, опубликуем. У меня карт от АТИ нет

Давайте у меня 3870 согласен ее погонять говорите че качать и как запускать

СTM родился в 2005. Только всех кому его давали заставили подписать NDA. Дошли до 5 беты, а потом ATI всех кинула, заявив, что дальше будет только brook+ и об CTM можно забыть.
ххх...... Канадская компания ....ххх..
Я вот думаю, NDA уже кончилась или нет? Могу я уже продемонстрировать ATIшный алгоритм перемножения матриц ( весьма крутой) народу, например на cude, или у меня будут юр.последствия? Ребята с http://sourceforge.net/projects/amdctm вроде написали, что для CTM NDA закончилась. А то сижу на исходниках CTM BLAS как собака на сене, обидно.

а что в нём крутого? Имеется в виду именно BLAS или просто умножение матриц? ATI публиковало детали реализации умножения матриц на CTM (тогда CTM назывался DVPM) в 2006. Потом, насколько я знаю, CTM не просто умер, а был развит и переименован в CAL, который сейчас свободно доступен для скачивания. CAL включает исходники умножения матриц и главу в программинг-гиде, расписывающую на три страницы как они работают.

1.крутой алгоритм перемножения матриц для CTM нигде не публиковался. Публиковался стандартный алгоритм, демонстрирующий общую идею матричных операций на GPU. А это и так все знают.
2. имеется ввиду именно BLAS3
3. CTM и CAL - абсолютно разные вещи
4. в 2006 уже была бета. Называлась ,правильно, Data Parallel Virtual Machine.

странно, я что-то не видел BLAS3 в CTM, может плохо глядел. Увы, так же не слышал про >120Gflop/s на X1900, которая была как раз актуальна в то время. Реализация, которая публиковалась (http://ati.amd.com/developer/siggraph06/dpvm_e.pdf) тоже была не совсем простой, она работала через fetch-4 расширение доступное только на картах ATI.

Да крутого может и ничего, просто делать такие вещи эффективно без local store - не очень удобно.

Так можно же спросить напрямую. Хотя, конечно, на такой прямой вопрос лично я бы ответил "нет, все продолжается".

Ну или можно завести анонима в ЖЖ, выложить все там и сделать вид что оно стало общеизвестным. (но это не я посоветовал , от всего откажусь, это хакеры!)

Потоковый BLAS это на самом деле жутко интересно, я бы почитал.

Не знаю, затрагивалась ли раньше эта тема. Думаю, ни для кого не секрет, что в численных методах умножение матрицы на матрицу в чистом виде встречается крайне редко. По сути, результаты такого умножения - это те же попугаи. Поэтому было бы очень интересно увидеть результаты умножения матрицы на вектор, а также умножение матрицы на вектор в разреженном формате, например, CSR. Не знаю, позволяет ли CUDA работать с разреженными форматами. У Гото Sparce BLAS нет, в MKL есть.

умножение матрицы на матрицу в численных методах встречается сплошь и рядом, например в методе Гаусса, QR факторизации и сведении к конденсированной форме. Что значит используется и является узким местом (или одним из них) в решении систем линейных уравнений, нахождении детерминантов, методе наименьших квадратов, нахождении собственных и сингулярных значений плотных матриц.

(предыдущий анонимус - это я)
Неточно выразился. Я имел в виду, что в реальных задачах такое умножение (заполненных матриц, да и вообще матриц) очень редко встречается.

По методу Гаусса, QR и прочая. А где это вы там умножение заполненных матриц нашли? Или будем умножение на матрицу поворота реализовывать по полной программе, как в тестах? http://blog.lexa.ru/2007/01/04/o_peremnozhenii_matric_i_prochix_arxitekt...

Когда решаем задачи мат физики, умножения матрицы на матрицу нет. Только умножение матрицы на вектор. Притом - разреженной матрицы, в том же формате CSR, например. Потому и сказал о попугаях.

На практике для реализации метода Гаусса (LU и Cholesky факторизации) используются блочные алгоритмы, смотрите, например, их реализацию в библиотеке LAPACK --- процедуры sgetrf и sportf. В блочных алгоритмах большинство операций совершается в умножении матриц и других, аналогичных процедурах из BLAS3. Блочные алгоритмы потребляют гораздо меньше пропускной способности памяти, поэтому работают значительно быстрее на современных машинах при больших размерах матриц.

Задачи матфизики разные бывают. В некоторых и плотные линейные системы алгебраических уравнений надо решать, и находить решение в смысле наименьших квадратов, и собственные значения искать и SVD --- всё там есть. Разреженные матрицы тоже встречаются очень часто, для них существуют отдельные библиотеки.

new

Хорошо, разобрались, в Гауссе это нужно.

Тем не менее, очень интересен пример задачи матфизики (имеющий реальное применение!) с получающейся плотной матрицей - раз. С используемым умножением матрица*матрица - два. Слова "разреженные матрицы тоже встречаются очень часто" видятся мне не совсем корректными. Продолжаю настаивать, что при решении уравнений матфизики применительно к реальным задачам подавляющее большинство получающихся систем имеют разреженные матрицы.

Почему сразу матфизика ?

Вот равновесная химия (в смысле, поиск равновесия) - небольшие плотные матрицы, но много (и много итераций в каждой задаче).
Хотя, конечно, сейчас можно уже и кинетику считать в таких задачах, если данные есть.

Почему матфизика? Да как обычно, у кого что болит...

Химия разная бывает. Горение с помощью МКЭ считается. Кстати, тот же Ansys МКЭ пользует. А МКЭ - это разреженные матрицы, притом умножения матрицы на матрицу при решении систем нет. А Ansys - это Ansys.

Всё же интересно будет увидеть результаты умножения разреженная матрица - вектор. Если уважаемый lexa, при наличии соответствующих возможностей, такие расчёты сделает, будет ему большое спасибо. Есть такое тихое подозрение, что видюхи здесь сильно провалятся.

Я вот ещё почему об этом говорю: в третьей серии про умножении матриц была получена превосходная масштабируемость в Double precision (Intel MKL 9.0). А вот у меня на C2Q 6600, Intel MKL версии 10.0.3.020 и Intel Compiler версии 10.1.015 на разреженной матрице самый лучший результат - это увеличение скорости расчета в 1.4 раза на 4 ядрах по сравнению с 1. Ключ компиляции -O3. С ключем -O0 там ситуация заметно получше.

если интересны какие-то результаты умножения разреженная матрица - вектор на GPU, то можно посмотреть He et al. "Efficient Gather and Scatter Operations on Graphics Processors", SC07. У них в одинарной точности на GeForce 8800 получилось в 6 раз быстрее чем в двойной на CPU/Intel MKL. А в сопряжённых градиентах 1.6-1.8 Гфлоп/c на GPU.

лучшее что я слышал --- 3-5 Гфлоп/c в SpMV на GPU. Недавно Университет Мериленда объявил конкурс на самую быструю реализацию, может там что лучше появится --- http://scriptroute.cs.umd.edu/gpucompete.

есть неплохая статья по изучению производительности SpMV на разных платформах, правда GPU там нет --- Williams et al. "Optimization of Sparse Matrix-Vector Multiplication on Emerging Multicore Platforms", SC07. Там тоже цитируется плохая масштабируемость по количеству ядер как на Intel так и на AMD процессорах. Дело в том, что производительность SpMV на больших матрицах ограничена пропускной способностью памяти. Поэтому количество задействованных ядер не столь важно --- они все на одном контроллере/шине сидят. Зато хорошо масштабировалось на двух процессорах. Лучше всего SpMV работало на Cell процессоре.

Спасибо за ссылки, пошукаю.

"...одинарной точности на GeForce 8800 получилось в 6 раз быстрее чем в двойной на CPU/Intel MKL..." - это некорректное сравнение, все-таки. Single и double сравнивать?

Насчет ПСП я в курсе. Инсайдерская, так сказать, информация. В MKL изначально оптимизация так сделана, что в ПСП всё упирается. Это накладывается на то, что у процов штеуда как раз там узкое место. Думаю, в тех тестах процы AMD заметно лучше масштабируемость показали.

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

Но видеокарты имеют шанс провалиться - они заточены (несмотря на все shared-memory) на потоковое чтение из памяти и на независимые (друг от друга) вычисления. А с разреженными матрицами будет хуже, там же на каждом шаге ветвления.

Про МКЭ - если правильно элементы нумеровать, то часто удается свести к матрицам, где ненулевые элементы вдоль диагонали, что жизнь сильно упрощает.

Да, я это и имел в виду, когда говорил про то, что видюхи могут провалиться. Разреженные матрицы и параллелить напрямую - то ещё удовольствие.

Перенумерация элементов в МКЭ дело хорошее, но достаточно затратное. А если задача параболическая и сетка изменяется каждый шаг по времени?

Как свести в виду, когда ненулевые элементы вдоль диагоналей в случае с более-менее сложной сеткой надо почитать. Вспоминается мне, что такое сведение задача сложная и опять же затратная.

В задачах, которыми я занимался (лет 10 тому как) генерация сетки была, конечно, затратной, но менее затратной чем все остальное.

Но занимался я 2D и осесимметричными задачами (течение смеси сред в пористой матрице, геология), на 3D тогда ресурсов не хватало. Сейчас - хватит, но уже не занимаюсь.

Кстати, о тестах для разреженных матриц. Пошукал по ссылкам, что здесь привели. Вот что нашел: The University of Florida Sparse Matrix Collection,
http://www.cise.ufl.edu/research/sparse/matrices/

Интегральные уравнения, по-моему, приводят к плотным системам. Ещё метод Хартри-Фока, задачи рассеяния.

Где можно почитать сведение метода Хартри-Фока к интегральным уравнениям?

Дифференциальные уравнения-разряженные, а интегральные -плотные.
Все задачи вычисл. томографии сводятся к интегральным.

new

Копаю дальше, расширяю кругозор. А какие задачи приводят к интегральным уравнениям?

например, дифур по объёму иногда сводят к интегральному уравнению по поверхности. Если не ошибаюсь, в механике сплошных сред встречается. Ещё вроде в элекромагнетизме. Не моя область, плохо разбираюсь.

Так это классика. Только что там потом? При получении системы линейных уравнений будут использоваться базисные функции, которые, конечно, можно так подобрать, что будет заполненная матриц, но зачем? Берется функция-шапочка, все довольны.

это в методе конечных элементов так получается. Но там и интеграл берётся по всему объёму, и степеней свободы O(N^3). А если свести к интегралу по поверхности, то степеней свободы будет O(N^2), но и зависимость будет не локальной. Насколько я понимаю, в это уравнение будет входить свёртка с функцией Грина.

Я гораздо ближе знаком с аналогичной задачей в компьютерной графике --- уравнения рендеринга и радиосити. Есть источники света, есть непрозрачная поверхность, возможно неодносвязная, рассеивает свет по известному закону. Пространство однородное, свет не рассеивает и не поглощает. Нужно рассчитать освещённость каждой точки поверхности. Получается, что она зависит освещённости каждой другой видимой точки поверхности. Выходит интегральное уравнение для решения на поверхности с интегралами по поверхности. Как вариант, можно оперировать не с освещённостью в точке, а с распределением по полусфере входящего и исходящего излучения в точке. Причём, если среда неоднородная и рассеивающая, то возможно придётся решать уравнение переноса излучения по всему пространству, которое уже дифференциальное (точнее, интегро-дифференциальное).

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

Благодарю за консультацию. А не Игорь С. ли вас зовут?

Всегда рад. Нет, я не Игорь, я другой ;-)

GTX 280 будет примерно в 1.5 раза быстрее, чем 9800GTX. печально, что это достигается банальным увеличением количества вычислительных блоков. чип становится больше, дороже и горячее. рано или поздно этот подход упрётся в ограничения техпроцесса. неспроста стали появляться двухчиповые видюхи, значит у производителя пока нет других идей по увеличению производительности.

Повторят историю с 8800 -- сначала сделают прожорливый чип (8800GTX), потом освоят новый техпроцесс и выпустят сравнимую по производительности, но более холодную и дешёвую версию (8800GTS на G92, он же 9800GTX).

а просто GT что такое?

Тоже G92 (он на нём был изначально), но с урезанным числом исполнительных блоков.

Тоже G92 (он на этом чипе был изначально), но с урезанным числом исполнительных блоков.

По блокам там 224/128 = 1.75 раза. Частоту забыл, может и чуть ниже.
Но сам путь - наращиванием блоков (и ядер)- мне кажется куда менее тупиковым, чем наращивание частоты.

Copyright © 2008-2009 Alex Tutubalin