Диссертация (Оптимизация размещения массивов в общей памяти), страница 10

11: Разбиение блока Block A матрицы A , блока Block B матрицы B и блокарезультата Block C матрицы С для размещения их в оперативной памяти.Данные внутри умножаемых блоков матрицы A хранятся вертикальнымиполосами ширины L. Внутри полосы матрицы A элементы хранятся блокамиразмера 4.  L . Будем называть такие блоки блоками второго уровня. Далее наРисунке 12 продемонстрировано размещение элементов, лежащих внутри блоковматрицы A.65Рис. 12: Размещение блоков матрицы A в оперативной памятиЭлементы, лежащие внутри блоков матрицы B и матрицы C, хранятся построкам. Высота блока второго уровня матрицы A ( Block A ) и ширина блокаматрицы B ( Block C ) выбирается, исходя из количества и размеров векторныхрегистров процессора. Так, процессоры Intel Sandy Bridge/Ivy Bridge/Haswellвключают 16 регистров размера 256 бит, каждый из которых может хранить 4числа типа double.

Для такой архитектуры эффективно разбивать матрицу A наблоки размера 4.  L , а матрицу B разбивать на блоки размера L  8 . Приумножении таких блоков восемь AVX - регистров резервируются для накоплениярезультата, в то время как оставшиеся восемь регистров используются длязагрузки данных из памяти.

Подробнее использование данной методики описано в[26].Параметр L выбирается, исходя из размера кеш-памяти L1, а именно - Lдолжно быть таким, чтобы умножаемый блок второго уровня матрицы A,умножаемый блок второго уровня матрицы B, а также блок матрицы С срезультатами вычислений помещались в кеш-памяти L1:4 * 8  4  L  8  L   sizeof (double)  sizeof ( L1Cache)66В этом случае элементы блока матрицы B не будут вытесняться междуитерациями.Умножение блоков 4  L и L  8 реализовано в отдельной процедуре, котораядля достижения максимальной производительности запрограммирована на языкеассемблера.

Поэтому при переносе на новую архитектуру, программнуюреализацию алгоритма умножения блоков потребуется каждый раз переписывать,что является непростой задачей. Так, для архитектуры Ivy Bridge размерпрограммной реализации алгоритма умножения блоков у автора составил порядка800 строк.Для внутреннего цикла применяется преобразование раскрутки. Раскруткацикла (loop unrolling) - преобразование циклов, заключающееся в многократномдублировании тела циклов с уменьшением числа итераций [99]. Такжезадействованыпрограммнаяпредвыборкаданных(softwareprefetching),переупорядочивание инструкций, выравнивание данных.Следует отметить, что для низкоуровневой оптимизации программнойреализации умножения блоков не требуется напрямую писать код на языкеассемблера.

В современных компиляторах языка C++ (например, ICC, GCC,Clang) существует возможность вставки в программную реализацию на языкеC++ ассемблерных команд с использованием встроенных функций (intrinsics).Также с помощью специальных параметров компиляции (например, параметр “–march”компилятораGCC)имеетсявозможностьуказатькомпиляторугенерировать векторные инструкции для конкретной архитектуры.Автор провел эксперименты, в которых сравнивалась программнаяреализация алгоритма умножения блоков, реализованная полностью на языкеассемблера, а также программная реализация алгоритма умножения блоков,оптимизированная компилятором ICC с использованием встроенных функций.Результаты экспериментов показали, что в случае использования ручнойоптимизациидостигаетсядополнительноеускорениев6%.Анализсгенерированного компилятором кода показал, что компилятор неоптимальнораспределяет векторные регистры, что увеличивает количество обращений к67оперативной памяти.

Так же компилятор не самым оптимальным образомрасставляет инструкции в раскрученном цикле, что приводит к замедлениюработы конвейера. Подробнее о порядке следования инструкций с цельюподстройки под конвейер описано в работе [80].3.3.Результаты численных экспериментовВчисленыхэкспериментахсравниваласьпроизводительностьреализованного алгоритма, а также библиотек Intel MKL, PLASMA и OpenBLAS.Тестирование проводилось на вычислительной системе с процессором Intel Corei7-3820 (3600MHz, LGA2011, L3 10240Kb) с выключенной опцией Turbo Boost.Размер кеш-памяти L1 равен 64 КБ, размер кеш-памяти L2 равен 1 МБ, размеркеш-памяти L3 равен 10 МБ.

В качестве компилятора использовался Intel C++Composer XE.Пиковая производительность вычислительной системы равнаP  m  H  T  4  2  3.6  4  115,2 ГФЛОПС , где P = 4- количество физических ядер процессора; m = 2 – количество команд, выполняемых параллельно; H = 3.6 Ггц – базовая тактовая частота; T=4 – количество чисел типа double, умещающихся в векторном регистре.Исходя из размера кеш-памяти L3, Di=512, Dk=256, Dj=8, L=256. Нижеприводятсярезультатыэкспериментовпроведениясравниваласьтрехэкспериментов.производительностьВпрограммнойкаждомизреализациипредложенного автором блочного алгоритма с алгоритмами умножения матрицвысокоуровневых библиотек OpenBLAS, MKL, PLASMA.БиблиотекаPLASMAвключаетвозможностьумноженияматриц,размещенных блочно.

В экспериментах присутствуют результаты измерений прииспользовании алгоритма, принимающего на вход стандартно размещенные68матрицы, а также - алгоритма, принимающего на вход матрицы, размещенныеблочно.Всеиспользуемыевэкспериментахпрограммныереализациидляконкретных входных размеров матриц запускались по 5 раз, и считалосьминимальное время работы. Производительность программной реализациисчиталась в гигафлопсах по следующей формуле:Performance 2  N1  N 2  N 3, где109  T T – время работы программы в секундах.На Рисунке 13 представлены результаты сравнения производительностипрограммных реализаций алгоритмов для случая, когда матрица A имеет размер512  N , а матрица B имеет размер N  512 .Рис. 13: Графики производительности программных реализаций блочногоалгоритма, а также алгоритмов пакетов MKL, OpenBLAS и PLASMA. Случай,когда размер матрицы A равен Nx512 и размер матрицы B равен 512xN.69В Таблице 7 приведены замеры времени работы программных реализацийблочного алгоритма, а также пакетов MKL, OpenBLAS и PLASMA.Таблица 7.

Таблица с замерами времени работы программных реализацийблочного алгоритма, а также пакетов MKL, OpenBLAS и PLASMA. Случай, когдаразмер матрицы A равен Nx512 и размер матрицы B равен 512xNNMKLOpenBLASPLASMA,PLASMA,Алгоритм с(среднее(среднееблочноестандартноедвойнымвремя ввремя вразмещение размещениеблочнымсек.)сек.)(среднее(среднееразмещениемвремя ввремя вданных (среднеесек.)сек.)время в сек.)320.0047940.0000830.0019200.0030480.000571640.0053070.0000940.0013950.0018810.0005211280.0060970.0002270.0019740.0028900.0009711600.0064070.0003270.0024240.0029150.0011791920.0069050.0004410.0029140.0029980.0016302240.0067820.0005800.0034920.0030280.0018442560.0080370.0007240.0041230.0040060.0022153840.0089090.0015460.0044930.0054260.0023095120.0099140.0026640.0050280.0072830.0033867680.0128810.0059340.0105490.0136710.00751910240.0172080.0105120.0138540.0215340.00999115360.6967320.0238050.0284890.0429260.02348920480.8237270.0418840.0484810.0763440.03998130720.8668120.0937760.1046350.1584290.09137740960.9299560.1647100.1840750.2632730.15510651200.8474340.2580170.2882320.4212720.24402061440.9065710.3696760.4145760.6140860.3470337071680.8613400.5038570.5617180.8199730.47572581920.8050710.6560250.7358251.0693950.61377492161.5454480.8309660.9304531.3416990.779768102401.7184351.0253411.1543021.6663720.959156112641.7054991.2403181.3899261.9997161.164584122881.7195231.4729921.6599062.3802481.382080133122.5383031.7306351.9438332.7887791.624753143362.7881912.0056032.2583463.2509291.880880153602.6933152.3037092.5825653.7215682.159584163843.4151502.6146152.9416974.2542852.446925174083.6283062.9653793.3232944.7878492.773088184323.5769903.3140983.7340345.3481683.105658194564.3053383.6962744.1629805.6891523.462284204804.4776064.0940694.6121546.5806663.826870В Таблице 8 приведены значения производительности, рассчитанныеотносительно пиковой производительности, программных реализаций блочногоалгоритма, а также пакетов MKL, OpenBLAS и PLASMA.Таблица 8.

Таблица со значениями производительности, рассчитаннымиотносительно пиковой производительности, программных реализаций блочногоалгоритма, а также пакетов MKL, OpenBLAS и PLASMA. Случай, когда размерматрицы A равен Nx512 и размер матрицы B равен 512xNNMKLOpenBLASPLASMA,PLASMA,Алгоритм сблочноестандартноедвойнымразмещение размещениеблочнымразмещением320.0018990.1101960.0047400.0029860.015941640.0068600.3889840.0261000.0193540.069910711280.0238860.6410020.0737770.0503960.1499541600.0355140.6958890.0938760.0780640.1930401920.0474540.7427020.1124500.1092920.2010802240.0657640.7689810.1277230.1472950.2419232560.0724840.8043940.1412910.1454170.2629993840.1471260.8478140.2917120.2415810.5677555120.2350290.8746220.4634390.3199380.6881777680.4070180.8835620.4970210.3834930.69724710240.5416420.8866360.6727980.4328390.93290715360.0301000.8809560.7361270.4885480.89283120480.0452610.8901460.7690230.4883490.93250630720.0967750.8945350.8017020.5294880.91802640960.1603630.9054140.8101650.5664500.96147551200.2749680.9031070.8084350.5531260.95491061440.3701250.9076720.8093670.5464120.96689571680.5302350.9064350.8130650.5569860.96003681920.7409570.9092990.8106870.5578140.97189592160.4885150.9085500.8114050.5627000.968204102400.5423930.9090320.8074730.5593390.971758112640.6612740.9092840.8114110.5639810.968416122880.7805520.9111910.8085860.5638810.971128133120.6205700.9101830.8103550.5648330.969498143360.6552110.9108740.8089340.5619480.971275153600.7786510.9103370.8120420.5635130.971091163840.6986790.9125980.8111280.5608680.975140174080.7424060.9083750.8105440.5626070.971363184320.8442570.9112280.8087500.5646600.972386194560.7815330.9103120.8082580.5914350.971834204800.8326480.9106520.8083580.5665490.97423572Полученные результаты свидетельствуют о том, что предложенный авторомблочный алгоритм быстрее чем алгоритмы из других библиотек.

В частности,данный алгоритм для N=20480 работает на 14% быстрее MKL, на 6% быстрееOpenBLAS, на 17% быстрее алгоритма пакета PLASMA с блочным размещениемвходных матриц и на 41% быстрее алгоритма пакета PLASMA со стандартнымразмещением входных матриц.На Рисунке 14 представлены результаты сравнения производительностипрограммных реализаций алгоритмов для случая, когда матрица A и матрица Bимеют размер N  N .Рис. 14: Графики производительности программных реализаций блочногоалгоритма, а также алгоритмов из пакетов MKL, OpenBLAS и PLASMA. Случай,когда матрица A и матрица B квадратные.В Таблице 9 приведены замеры времени работы программных реализацийблочного алгоритма, а также пакетов MKL, OpenBLAS и PLASMA.73Таблица 9.