Использование жесткого диска в матричных вычислениях Текст научной статьи по специальности «Математика»

Np — число операций умножения в алгоритме с полиномиальной сложностью,

Rp — нижняя оценка числа операций умножения в алгоритме с экспоненциальной сложностью,

Nt — время вычисления матричной степени в алгоритме с полиномиальной сложностью,

Rt — время вычисления матричной степени при использовании алгоритма с экспоненциальной сложностью,

NRt — время вычисления матричной степени при использовании алгоритма с полиномиальной сложностью,

NRp — оценка числа операций умножения при использовании алгоритма с полиномиальной сложностью.

Литература

1. Robuk V.N. A constructive formula for function of matrix. Alternative to the Lagrange-Silvester formula. Nuclear Instruments and Methods in Physics Research A 534 (2004) 319-323.

2. Ланкастер П. Теория матриц. М.: Наука, 1973, с. 280.

3. Математическая энциклопедия/ Гл. ред. И.М. Виноградов. — М.: Сов. Энциклопедия, 1979. Т.5.

Поступила в редакцию 15 декабря 2007г.

Использование жесткого диска в матричных вычислениях

© М.С.Зуев

Вычисления с использованием жесткого диска реализованы, например, в известном пакете численных методов POOCLAPACK (Parallel out-of-core Linear algebra package, см. [5]).

В данном пакете реализованы два способа хранения матриц на диске, один из которых повторяет способ хранения матриц в ОЗУ по строкам или столбцам, а второй предполагает хранение матрицы по блокам. Во втором случае требуется использование блочных матричных алгоритмов.

В современной литературе блочные и блочно-рекурсивные матричные методы рассматриваются как одни из самых эффективных. В отличие от стандартных строчно-ориентированных методов, такие методы более эффективно используют процессорный кэш, и поэтому позволяют достигать производительности процессора, близкой к пиковой (см. [2-4]). Например, в работе [2] указывается достижение 92 % производительности процессора IBM Power 3 с блочной процедурой разложения Холецкого и 76 % с соответствующей поэлементной процедурой.

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

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

Для исследования эффективности таких алгоритмов были проведены вычислительные эксперименты на однопроцессорной машине. Сравнивались алгоритмы умножения и обращения матриц (одностороннее обращение, см. [1]) размера 4096 х 4096 в конечном поле по 28-битному простому модулю. Как и в POOCLAPACK, матрица разбивается на блоки фиксированного размера, которые записываются в отдельные файлы. Матрице сопоставляется папка на диске. Если разбить матрицу на четыре равных блока, то возможно их хранение в отдельных папках.

Для определения лучшего алгоритма умножения были рассмотрены четыре алгоритма: стандартный блочный алгоритм (алгоритм 0), стандартный блочный алгоритм с записью блоков на

Рис. 1. Сравнение скорости выполнения алгоритма 3 с различными параметрами

диск (алгоритм 1). алгоритм Штрассена с листом переменного размера (алгоритм 2) и он же с записью блоков на диск (алгоритм 3).

Под стандартным алгоритмом будем понимать алгоритм умножения матричных строк на матричные столбцы. Рассмотрены три реализации этого алгоритма: поэлементная, т.е. вычисление произведения по одному элементу, блочная, т.е. вычисление произведения по блокам фиксированного размера и блочно-рекурсивная, т.е. вычисление произведения по блокам, порядок которых равен половине порядка матрицы, с помощью рекурсивной формулы. Самым быстрым оказывается второй алгоритм, поэтому эксперименты были проведены с ним. Он обозначен как алгоритм 0. Этому-алгоритму передается один параметр размер блока.

Алгоритм 1 стандартный блочный алгоритм с использованием жесткого диска. Он имеет ту же вычислительную структуру, что и алгоритм 0: блочные строки фиксированного размера умножаются на блочные столбцы. Каждый блок хранится в отдельном файле. В текущий момент вычислений в памяти хранятся только текущие блоки операндов и блок результата. Произведение блоков вычисляется с помощью стандартного поэлементного алгоритма. Алгоритму передается один параметр размер блока, хранящегося в одном файле.

Алгоритм 2 алгоритм Штрассена для матриц, элементами которых являются блоки фиксированного размера. Эти блоки, в свою очередь, умножаются с помощью стандартного поэлементного алгоритма. Алгоритму передается один параметр размер блока.

Алгоритм 3 имеет вычислительную структуру, аналогичную алгоритму 2. Блоки матриц хранятся в файлах. Их произведение вычисляется с помощью алгоритма 1. Алгоритму передаются 2 параметра размер блока матрицы, хранящегося в одном файле, и размер блоков матриц, произведение которых вычисляется с помощью алгоритма 1.

Результаты вычислительных экспериментов с алгоритмом 3 приведены на рис. 1. Вычисления производились с машиной на базе АШоп-3800. 1 Гб ОЗУ. ОС ^¥іпХР. для записи на диск использовался пустой ХТГБ-раздел. На рисунке приведены только самые быстрые алгоритмы, со временем выполнения не более 2000 с. По оси отложен размер блоков матриц, для которых выполняется умножение с помощью алгоритма 1. по оси время выполнения программы в секундах. Сплошная линия соответствует времени выполнения программ с блоком матриц, хранящихся в одном файле, равным 64. далее в порядке учащеиия пунктирных линий 128. 256 и 512. Лучшая скорость получена при храпении блоков порядка 128 в одном файле и вычислении произведения блоков порядка 2048 с помощью алгоритма 1 (см. рис. 1). Алгоритм с параметрами (128, •) участвовал в следующем эксперименте в качестве алгоритма 3.

Результаты сравнения всех алгоритмов приведены на рис. 2. По оси отложен размер листа дерева алгоритма, по оси время выполнения программы в секундах. Для алгоритмов 0-2 размер листа дерева алгоритма принимаем равным размеру блоков, произведение которых вычисляется с помощью стандартного поэлементного алгоритма, а для алгоритма 3 алгоритма 1. На рисунке приведены только самые быстрые алгоритмы, со временем выполнения не более 2600 с. Сплошная линия соответствует временам выполнения алгоритмов 0. остальные в порядке учащения пунктир-

2344 -14113

596

I ? г.;'

V ' ...

«—I--------I—

32. 256. 512.

1024.

—Ь X 2048.

Рис. 2. Выбор лучшего алгоритма умножения

пых линий алгоритмам 1-3.

Результаты экспериментов 1-2:

• Лучшим алгоритмом оказался алгоритм 2, вторым - алгоритм 3, далее - алгоритмы 1 и 0. Алгоритмы с использованием жесткого диска не оказались самыми медленными, как ожидалось.

•

поднялся 9787 с.. а алгоритм Штрассена с размером листа 1 (обычный) 6544 с. Блочные и блочно-рекурсивные алгоритмы с использованием блоков на нижнем уровне и алгоритмы с использованием жесткого диска существенно быстрее алгоритмов в общепринятой записи.

•

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

•

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

Сравнивались алгоритмы одностороннего обращения матриц: алгоритм, использующий блочный стандартный алгоритм матричного умножения (алгоритм 0) и алгоритм, использующий запись блоков матриц на диск (алгоритм 1). Блоки матриц, считываемые с диска, умножаются с помощью стандартного алгоритма. Алгоритму 0 передается один параметр размер блока, используемого блочным алгоритмом умножения. Алгоритму 1 передается один параметр размер блока, хранящегося в одном файле.

На рис. 3 приведены результаты их сравнения. По оси отложен размер блока алгоритма умножения. по оси время выполнения программы. Сплошная линия соответствует времени выполнения алгоритма 0. пунктирная времени выполнения алгоритма 1.

Заметим, что алгоритм 0 с параметром 512 работает в два раза быстрее алгоритма 1 с параметром 256. а он. в свою очередь, работает в два раза быстрее алгоритма, использующего стандартное поэлементное умножение матричных блоков (этот алгоритм выполнялся 2635 секунд).

Результаты экспериментов 1-3 позволяют сделать вывод об эффективности алгоритмов, использующих хранение матриц в файлах. Такие алгоритмы эффективно используют оперативную память и кэш процессора, поэтому сравнимы с блочными алгоритмами по эффективности. Поэтому такие алгоритмы могут представлять интерес для реальных высокопроизводительных вычислений

4------1-----1-----------1-----------------------h x

256 512 1024 2048

Рис. 3. Выбор лучшего алгоритма одностороннего обращения

Литература

1. Малашонок Г. И.. Зуев М. С. Обобщенный алгоритм вычисления обратной матрицы // XI Державинские чтения. Тез. докл. Тамбов: Изд-во ТГУ им. Г.Р.Державина. 2006. С. 48-50.

2. Gustavson F. G. High-performance linear algebra algorithms using new generalized data structures for matrices // IBM Journal for Research and development. January 2003. V. 41. Issue 1. P. 31-55.

3. Gustavson F. G. Recursion leads to automatic variable blocking for dense linear-algebra algorithms // IBM Journal for Research and development. November 1997. V. 41. Issue 6. P. 737-756.

4. Lorton K. P.. Wise D. S. Analyzing Block Locality in Morton-Order and Morton-Hybrid Matrices // Proc. MEDEA Workshop’06. New York: ACM Press', 2006. P. 5-12.

5. Reiley W. C., van de Geijn R. A. POOCLAPACK: Parallel Out-of-Core Linear Algebra Package. Technical report CS-TR-99-33. Austin, TX, USA: University of Texas at Austin. 1999.

Поступила в редакцию 15 декабря 2007г.

О вычислении ядра оператора, действующего в модуле

© Г.И. Малашонок

Пусть Я - область с единицей, Яп - модуль над Я размерности и, М € Япхп - оператор, действующий в Я”. Считается фиксированным базис модуля и оператор представлен матрицей. Требуется найти ядро оператора. Предлагается алгоритм, состоящий из двух этапов. Сначала находится обобщенная присоединенная матрица А, такая, что матрица S = АМ имеет ступенчатый вид. Затем по матрице £ находится ядро оператора М. Введем некоторые обозначения.

Пусть полугруппа Рп образована матрицами, у которых число единичных элементов совпадает с рангом матрицы, а остальные элементы равны нулю. Полугруппа Бп образована диагональными матрицами порядка и с элементами 0 и 1 на диагонали, |Еп| =2” и единичная матрица I - единица в Бп и в Рп. Па Бп введем частичный по рядок: I < 7 ^ 7 — I € Юп и чертой будем обозначать инволюцию на Пп : I = I — I. Для каждой матрицы Е € Рп определим диагональные матрицы: 1е = ЕЕТ € Юп, = ЕтЕ € Юп. Для матрицы Е матрица 1е является левым аннулятором, а

матрица 1е - правым аннулятором.

Обозначим множество ступенчатых матриц порядка и: Бп = |ЗЕ € Рп, Эй € Д\0 : ^ =

1е¿Е = Б1е}. Другими словами, для каждой ступенчатой матрицы £ € Бп найдется матрица