Участник:Blizn/Хранение ненулевых элементов разреженной матрицы. Умножение разреженной матрицы на вектор.: различия между версиями
Blizn (обсуждение | вклад) |
Blizn (обсуждение | вклад) |
||
Строка 190: | Строка 190: | ||
Параллельная реализация алгоритма заключается в следующем: | Параллельная реализация алгоритма заключается в следующем: | ||
− | # | + | # Каждый процессор считывает входные данные из файла, зависящего от его номера. Данные имеют вид: вектор <math>b</math> и набор определённых строк матрицы <math>A</math> в формате '''RR(C)U''' в зависимости от номера процессора и их числа; |
− | + | # Процессоры проводят вычисления значений вектора <math>c</math> в зависимости от имеющихся данных; | |
− | + | # Все процессоры отправляют результаты вычисления процессору, выбранному главным, и главный процессор выводит результат вычислений. | |
− | # | ||
− | # Все процессоры отправляют результаты вычисления | ||
− | Проведём исследование масштабируемости параллельной реализации алгоритма умножения разреженной матрицы на вектор. Исследование | + | Проведём исследование масштабируемости параллельной реализации алгоритма умножения разреженной матрицы на вектор. Исследование проводилось на суперкомпьютере "Ломоносов"<ref name="Lom">Воеводин Вл., Жуматий С., Соболев С., Антонов А., Брызгалов П., Никитенко Д., Стефанов К., Воеводин Вад. Практика суперкомпьютера «Ломоносов» // Открытые системы, 2012, N 7, С. 36-39.</ref> [http://parallel.ru/cluster Суперкомпьютерного комплекса Московского университета]. |
== Динамические характеристики и эффективность реализации алгоритма == | == Динамические характеристики и эффективность реализации алгоритма == |
Версия 23:11, 29 ноября 2016
Умножение разреженной матрицы на вектор | |
Последовательный алгоритм | |
Последовательная сложность | O(l) |
Объём входных данных | 2l+m+n+1 |
Объём выходных данных | n |
Параллельный алгоритм | |
Высота ярусно-параллельной формы | O(m) |
Ширина ярусно-параллельной формы | O(n) |
Выполнила: И.В. Близнякова (611 группа).
Содержание
- 1 Свойства и структура алгоритма
- 1.1 Общее описание алгоритма[1]
- 1.2 Математическое описание алгоритма
- 1.3 Вычислительное ядро алгоритма
- 1.4 Макроструктура алгоритма
- 1.5 Схема реализации последовательного алгоритма
- 1.6 Последовательная сложность алгоритма
- 1.7 Информационный граф
- 1.8 Ресурс параллелизма алгоритма
- 1.9 Входные и выходные данные алгоритма
- 1.10 Свойства алгоритма
- 2 Программная реализация алгоритма
- 2.1 Особенности реализации последовательного алгоритма
- 2.2 Локальность данных и вычислений
- 2.3 Возможные способы и особенности параллельной реализации алгоритма
- 2.4 Масштабируемость алгоритма и его реализации
- 2.5 Динамические характеристики и эффективность реализации алгоритма
- 2.6 Выводы для классов архитектур
- 2.7 Существующие реализации алгоритма
- 3 Литература
1 Свойства и структура алгоритма
1.1 Общее описание алгоритма[1]
Разрежённая матрица — это матрица с преимущественно нулевыми элементами. В противном случае, если бо́льшая часть элементов матрицы ненулевые, матрица считается плотной.
Среди специалистов нет единства в определении того, какое именно количество ненулевых элементов делает матрицу разрежённой. Разные авторы предлагают различные варианты. Для матрицы порядка n число ненулевых элементов:
- есть O(n). Такое определение подходит разве что для теоретического анализа асимптотических свойств матричных алгоритмов;
- в каждой строке не превышает 10 в типичном случае;
- ограничено n^{1+\gamma}, где \gamma \lt 1.
- таково, что для данного алгоритма и вычислительной системы имеет смысл извлекать выгоду из наличия в ней нулей.
1.1.1 Хранение разреженной матрицы
1.1.1.1 Формат RR(C)O
Рассмотрим сначала формат RR(C)O. Сокращенное название данного формата происходит от английского словосочетания "Row - wise Representation Complete and Ordered" (строчное представление, полное и упорядоченное). В данном формате вместо одного двумерного массива, используются три одномерных. Значения ненулевых элементов матрицы и соответствующие им столбцовые индексы хранятся в этом формате по строкам в двух массивах AN и JA. Массив указателей IA, используется для ссылки на компоненты массивов AN и JA, с которых начинается описание очередной строки. Последняя компонента массива IA содержит указатель первой свободной компоненты в массивах AN и JA, т.е. равна числу ненулевых элементов матрицы, увеличенному на единицу. Здесь уместно привести пример.
Рассмотрим матрицу A:
- \begin{pmatrix} 0 & 0 & 0 & 2 & 0 \\ 1 & 0 & 3 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 \\ 0 & 6 & 0 & 0 & 0 \\ 0 & 0 & 4 & 0 & 0 \\ \end{pmatrix},
тогда ее представление в формате RR(C)O будет иметь вид:
IA = [ 1 2 4 4 5 6 ] JA = [ 4 1 3 2 3 ] AN = [ 2 1 3 6 4 ]
Т.е. массив AN содержит все не нулевые элементы исходной матрицы A, массив JA номер столбца в котором находится соответствующий элемент из AN и наконец массив IA содержит номер с которого начинается описание элементов в массивах JA и AN. Таким образом информация об элементах 2-ой строки матрицы хранится в элементах с IA[2] = 2 по IA[3] - 1 = 3 включительно массивов JA и AN. Можно обратить внимание, что IA[3] = IA[4] = 4, а это означает, что 3-я строка матрицы A нулевая.
В общем случае описание r-й строки матрицы A хранится в компонентах с IA[r] до IA[r + 1] - 1 включительно массивов AN и JA. Если IA[r + 1] = IA[r], то это означает, что r-я строка нулевая. Количество элементов в массиве IA на единицу больше, чем число строк исходной матрицы, а количество элементов в массивах JA и AN равно числу ненулевых элементов исходной матрицы.
Данный способ представления называют полным, поскольку представлена вся матрица A, упорядоченным, поскольку элементы каждой строки матрицы A хранятся в соответствии с возрастанием столбцовых индексов, и строчным, поскольку информация о матрице A указывается по строкам.
Массивы IA и JA представляют портрет (структуру) матрицы A, задаваемый как множество списков смежности ассоциированного с A графа. Если алгоритм, реализующий какую-либо операцию над разреженными матрицами, разбит на этапы символической обработки, на котором определяется портрет результирующей матрицы, и численной обработки, на котором определяются значения элементов результирующей матрицы, то массивы IA и JA заполняются на первом этапе, а массив AN - на втором.
1.1.1.2 Формат RR(C)U
Рассмотрим теперь формат RR(C)U.
Сокращенное название данного формата происходит от английского словосочетания "Row - wise Representation Complete and Unordered" (строчное представление, полное, но неупорядоченное). Формат RR(C)U отличается от RR(C)O тем, что в данном случае соблюдается упорядоченность строк, но внутри каждой строки элементы исходных матриц могут храниться в произвольном порядке. Для матрицы A нашего примера вполне можно было бы использовать и строчное представление, полное, но неупорядоченное такое:
IA = [ 1 2 4 4 5 6 ] JA = [ 4 3 1 2 3 ] AN = [ 2 1 3 6 4 ]
Такие неупорядоченные представления могут быть очень удобны в практических вычислениях. Результаты большинства матричных операций получаются неупорядоченными, а их упорядочение стоило бы значительных затрат машинного времени. В то же время, за немногими исключениями, алгоритмы для разреженных матриц не требуют, чтобы их представления были упорядоченными.
1.1.1.3 Замечания
Несколько замечаний по поводу рассмотренных форматов представления:
- Очевидно, что представление матрицы в формате RR(C)O так же является и представлением в формате RR(C)U, но не наоборот.
- Из представления матрицы в формате RR(C) нельзя получить информацию о точном количестве столбцов исходной матрицы.
- Целесообразно (в вопросе экономии памяти) использовать представления RR(C) в случае, если матрица содержит значительное число нулевых элементов.
1.1.2 Умножение разреженной матрицы на вектор
Важным приложением этих алгоритмов является вычисление векторов Ланцоша, необходимое при итерационном решении линейных уравнений методом сопряженных градиентов, а также при вычислении собственных значений и собственных векторов матрицы. Достоинство этих процедур, с вычислительной точки зрения, состоит в том, что единственная требуемая матричная операция - это повторное умножение матрицы на последовательность заполненных векторов; сама матрица не меняется.
Мы рассмотрим умножение разреженной матрицы общего вида, хранимой в форме RR(C)U посредством массивов IA, JA, AN на заполненный вектор-столбец.
1.2 Математическое описание алгоритма
Исходные данные: разреженная матрица общего вида A с элементами a_{ij} (i = 1,...,n и j = 1,...,m). Заполненный вектор-столбец b с элементами b_{j} (j =1,...,m).
Вычисляемые данные: заполненный вектор-столбец c с элементами c_{i} (i = 1,...,n).
Формулы метода:
c_{i} = \sum\limits_{k = 1}^{l_{i}} a_{i,j=j(k)}b_{j=j(k)},
где l_{i} - количество ненулевых элементов строки i матрицы A, j(k) - индекс k-го ненулевого элемента матрицы A.
1.3 Вычислительное ядро алгоритма
Вычислительное ядро последовательной версии умножения разреженной матрицы на вектор можно составить из множественных (всего их n) вычислений скалярных произведений строк матрицы:
c_{i} = \sum\limits_{k = 1}^{l_{i}} a_{i,j=j(k)}b_{j=j(k)}.
1.4 Макроструктура алгоритма
Как записано и в описании ядра алгоритма, основную часть метода составляют множественные (всего n) вычисления сумм:
c_{i} = \sum\limits_{k = 1}^{l_{i}} a_{i,j=j(k)}b_{j=j(k)},
которые могут вычисляться в произвольном порядке.
1.5 Схема реализации последовательного алгоритма
Далее предполагаем, что разреженная матрица общего вида A хранится в форме RR(C)U посредством массивов IA, JA, AN. Последовательность исполнения метода следующая:
Выполнять для i от 1 до n
- c_{i} = 0
- IAA = IA[i]
- IAB = IA[i + 1] - 1
- c_{i} = \sum\limits_{k = IAA}^{IAB} AN[k]b_{JA[k]}.
После этого (если i \le n) происходит переход к шагу 1 с бо́льшим i.
1.6 Последовательная сложность алгоритма
Для умножения разреженной матрицы общего вида, хранимой в форме RR(C)U, размером n \times m на заполненный вектор m \times 1 в последовательном (наиболее быстром) варианте требуется:
- l сложений,
- l умножений.
Умножения и сложения составляют основную часть алгоритма.
При классификации по последовательной сложности, таким образом, алгоритм умножения разреженной матрицы на вектор относится к алгоритмам O(l).
1.7 Информационный граф
Опишем граф алгоритма[2][3][4] как аналитически, так и в виде рисунка.
Граф алгоритма состоит из двух групп вершин, расположенных в целочисленных узлах двух областей одной размерности.
Первая группа вершин расположена в двумерной области, соответствующая ей операция вычисляет функцию a \cdot b. Естественно введённые координаты области таковы:
- i — меняется в диапазоне от 1 до n, принимая все целочисленные значения;
- j — меняется в диапазоне от 1 до k, принимая все целочисленные значения,
где k = k(i) = IA[i+1]-IA[i].
Аргументы операции следующие:
- a - элемент входных данных, а именно AN[IA_{i}+j-1].
- b - элемент входных данных, а именно b[IA_{i}+j-1].
Результат срабатывания операции является промежуточным данным алгоритма.
Вторая группа вершин расположена в двумерной области, соответствующая ей операция a + b. Естественно введённые координаты области таковы:
- i - меняется в диапазоне от 1 до n, принимая все целочисленные значения;
- j - меняется в диапазоне от 1 до k-1, принимая все целочисленные значения,
где k = k(i) = IA[i+1]-IA[i].
Аргументы операции следующие:
- a:
- при j = 1 - результат срабатывания операции, соответствующей вершине из первой группы, с координатами (i, j);
- при j \gt 1 - результат срабатывания операции, соответствующей вершине из второй группы, с координатами (i, j - 1);
- b — результат срабатывания операции, соответствующей вершине из первой группы, с координатами (i, j + 1).
Результат срабатывания операции:
- при j \lt k - 1 является промежуточным данным алгоритма;
- при j = k - 1 является выходным данным c_{i}.
Описанный граф можно посмотреть на рис.1. Здесь вершины первой группы обозначены красным цветом и отмечены знаком умножения, вершины второй - зелёным цветом и знаком сложения. Вершины, соответствующие входным данным обозначены белым цветом и выходным - синим.
1.8 Ресурс параллелизма алгоритма
Для умножения разреженной матрицы общего вида, хранимой в форме RR(C)U, размерности n \times m на заполненный вектор m \times 1 в параллельном варианте требуется последовательно выполнить следующие ярусы:
- Не более чем m сложений и умножений (n вычислений в каждом из ярусов)
При классификации по высоте ЯПФ, таким образом, алгоритм умножения разреженной матрицы на вектор относится к алгоритмам со сложностью O(m). При классификации по ширине ЯПФ его сложность будет O(n).
1.9 Входные и выходные данные алгоритма
Входные данные: разреженная матрица общего вида A размерности n \times m, хранимая в форме RR(C)U посредством массивов IA, JA, AN. Заполненный вектор-столбец b с элементами b_{j} размерности m \times 1.
Объём входных данных: 2l+m+n+1, где l - количество ненулевых элементов в матрице A.
Выходные данные: заполненный вектор-столбец c с элементами c_{i} (i = 1,...,n).
Объём выходных данных: n.
1.10 Свойства алгоритма
Соотношение последовательной и параллельной сложности в случае неограниченных ресурсов, как хорошо видно, является линейной.
При этом вычислительная мощность алгоритма, как отношение числа операций к суммарному объему входных и выходных данных – константа.
Пусть l - количество ненулевых элементов матрицы A размерности n \times m. Пусть x - объем памяти, используемый для хранения значения элемента матрицы, y - объём памяти, используемый для хранения номера столбца или строки. В таком случае для хранения матрицы в стандартном представлении нам потребуется объем памяти, равный x \cdot n \cdot m, для хранения в формате RR(C) - y(n + 1) + (x + y)l. Таким образом хранение в формате RR(C) не является эффективным для матриц, в которых l \gt \frac{xnm - y(n+1)}{x+y}.
2 Программная реализация алгоритма
2.1 Особенности реализации последовательного алгоритма
2.2 Локальность данных и вычислений
2.3 Возможные способы и особенности параллельной реализации алгоритма
2.4 Масштабируемость алгоритма и его реализации
Параллельная реализация алгоритма заключается в следующем:
- Каждый процессор считывает входные данные из файла, зависящего от его номера. Данные имеют вид: вектор b и набор определённых строк матрицы A в формате RR(C)U в зависимости от номера процессора и их числа;
- Процессоры проводят вычисления значений вектора c в зависимости от имеющихся данных;
- Все процессоры отправляют результаты вычисления процессору, выбранному главным, и главный процессор выводит результат вычислений.
Проведём исследование масштабируемости параллельной реализации алгоритма умножения разреженной матрицы на вектор. Исследование проводилось на суперкомпьютере "Ломоносов"[5] Суперкомпьютерного комплекса Московского университета.
2.5 Динамические характеристики и эффективность реализации алгоритма
2.6 Выводы для классов архитектур
2.7 Существующие реализации алгоритма
Существующие реализации:
SciPy [1],
MatLab [2],
Intel MKL [3].
3 Литература
- ↑ С.Писсанецки. Технология разреженных матриц. - М.: Мир, 1988.
- ↑ Воеводин В.В. Математические основы параллельных вычислений// М.: Изд. Моск. ун-та, 1991. 345 с.
- ↑ Воеводин В.В., Воеводин Вл.В. Параллельные вычисления. – СПб.: БХВ - Петербург, 2002. – 608 с.
- ↑ Фролов А.В.. Принципы построения и описание языка Сигма. Препринт ОВМ АН N 236. М.: ОВМ АН СССР, 1989.
- ↑ Воеводин Вл., Жуматий С., Соболев С., Антонов А., Брызгалов П., Никитенко Д., Стефанов К., Воеводин Вад. Практика суперкомпьютера «Ломоносов» // Открытые системы, 2012, N 7, С. 36-39.