Участник:Pandvik/Ортогонализация Грама - Шмидта: различия между версиями
Pandvik (обсуждение | вклад) |
Pandvik (обсуждение | вклад) |
||
| Строка 104: | Строка 104: | ||
== Ресурс параллелизма алгоритма == | == Ресурс параллелизма алгоритма == | ||
| − | На основе графа зависимости по данным, представленном на Рис. 1 или Рис. 2 раздела [[#Информационный граф]] можно выделить следующий ресурс | + | На основе графа зависимости по данным, представленном на Рис. 1 или Рис. 2 раздела [[#Информационный граф]] можно выделить следующий ресурс для параллелизации. После вычисления <math>b_{i}</math> можно начинать параллельно вычислять все <math>proj_{b_i}a_{j}</math>. На Рис. 2 блоки, которые можно вычислять параллельно представлены столбцами. |
Кроме того, при очень больших размерностях можно также распараллелить вычисление проекций. Так как в функции проекции участвует пара векторов <math>a</math> и <math>b</math>, а сама операция проекции определяется как <math>\mathbf{proj}_b a = \frac{\left \langle a,b \right \rangle}{\left \langle b,b \right \rangle }b</math>, где <math>\left \langle a,b \right \rangle</math> - скалярное произведение векторов <math>a</math> и <math>b</math>, то с помощью использования векторных операций можно получить некоторый выигрыш в производительности. | Кроме того, при очень больших размерностях можно также распараллелить вычисление проекций. Так как в функции проекции участвует пара векторов <math>a</math> и <math>b</math>, а сама операция проекции определяется как <math>\mathbf{proj}_b a = \frac{\left \langle a,b \right \rangle}{\left \langle b,b \right \rangle }b</math>, где <math>\left \langle a,b \right \rangle</math> - скалярное произведение векторов <math>a</math> и <math>b</math>, то с помощью использования векторных операций можно получить некоторый выигрыш в производительности. | ||
| Строка 112: | Строка 112: | ||
Рассмотрим случай, когда присутствует возможность запуска бесконечного числа потоков. Попробуем оценить сложность вычисления ортогональных векторов, при условии, что параллельно будут запускаться все проекции использующие одинаковую <math>b_i</math>. В этом случае, останется [[критический путь]] зависимости по данным: <math>b_1 \leftarrow proj_{b_2}a_{2} \leftarrow b_2 \leftarrow proj_{b_3}a_{3} \leftarrow \ldots \leftarrow proj_{b_{N-1}}a_{N-1} \leftarrow b_N</math>. Следовательно для оценки сложность распараллеленной версии алгоритма достаточно почитать сложность представленной цепочки вычислений. | Рассмотрим случай, когда присутствует возможность запуска бесконечного числа потоков. Попробуем оценить сложность вычисления ортогональных векторов, при условии, что параллельно будут запускаться все проекции использующие одинаковую <math>b_i</math>. В этом случае, останется [[критический путь]] зависимости по данным: <math>b_1 \leftarrow proj_{b_2}a_{2} \leftarrow b_2 \leftarrow proj_{b_3}a_{3} \leftarrow \ldots \leftarrow proj_{b_{N-1}}a_{N-1} \leftarrow b_N</math>. Следовательно для оценки сложность распараллеленной версии алгоритма достаточно почитать сложность представленной цепочки вычислений. | ||
| − | <math>{\color{Blue}\bullet} \quad | + | <math>{\color{Blue}\bullet} \quad M \cdot 3</math> умножения в одной проекции; |
| − | <math>{\color{Blue}\bullet} \quad ( | + | <math>{\color{Blue}\bullet} \quad (M-1) \cdot 2</math> сложения; |
<math>{\color{Blue}\bullet} \quad 1</math> деление; | <math>{\color{Blue}\bullet} \quad 1</math> деление; | ||
| Строка 120: | Строка 120: | ||
В цепочке зависимости данных присутствует <math>N-1</math> проекция, значит всего: | В цепочке зависимости данных присутствует <math>N-1</math> проекция, значит всего: | ||
| − | <math>{\color{Blue}\bullet} \quad N \cdot | + | <math>{\color{Blue}\bullet} \quad N \cdot M \cdot 3</math> - умножений; |
| − | <math>{\color{Blue}\bullet} \quad N \cdot ( | + | <math>{\color{Blue}\bullet} \quad N \cdot (M-1) \cdot 2</math> - сложений; |
<math>{\color{Blue}\bullet} \quad N \cdot 1</math> - делений. | <math>{\color{Blue}\bullet} \quad N \cdot 1</math> - делений. | ||
| − | Получаем итоговую сложность <math>O(N | + | Получаем итоговую сложность <math>O(N \cdot M)</math>. |
Если дополнительно рассматривать возможность параллелизации векторных операций (скалярного произведения и умножения вектора на число), [[критический путь]] можно сократить до следующих операций: | Если дополнительно рассматривать возможность параллелизации векторных операций (скалярного произведения и умножения вектора на число), [[критический путь]] можно сократить до следующих операций: | ||
| Строка 134: | Строка 134: | ||
<math>{\color{Blue}\bullet} \quad N </math> делений; | <math>{\color{Blue}\bullet} \quad N </math> делений; | ||
| − | <math>{\color{Blue}\bullet} \quad N* | + | <math>{\color{Blue}\bullet} \quad N*logM </math> сложений; |
| − | Следовательно оценка сложности параллельного алгоритма станет равна <math>O(N | + | Следовательно оценка сложности параллельного алгоритма станет равна <math>O(N logM)</math>. |
== Входные и выходные данные алгоритма == | == Входные и выходные данные алгоритма == | ||
Версия 21:57, 15 октября 2016
Авторы описания алгоритма: Павлов Андрей, Филимонов Владимир.
Содержание
- 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].
Процесс ортогонализации Грама-Шмидта используется для квадратных невырожденных матриц, которые преобразуются, либо уже преобразованы, к верхнему(нижнему) треугольному виду.
Процесс ортогонализации Грама-Шмидта нашёл применение в оптимизации оценивания параметров моделей управления объектом[2], в протоколах безопасности[3], в обработке сигналов[4], в вычислении локальных минимумов целочисленных решёток и многом другом.
Обычно, процесс ортогонализации используется как промежуточный шаг в других алгоритмах для уменьшения количества вычислений.
1.2 Математическое описание алгоритма
Исходные данные: квадратная матрица A с линейно независимыми векторами [math]\mathbf{a}_1,...,\mathbf{a}_N | \mathbf{a}_i \in \R^{M}[/math].
Определяется оператор проекции [math]\mathbf{proj}_b a = \frac{\left \langle a,b \right \rangle}{\left \langle b,b \right \rangle }b[/math], где [math]\left \langle a,b \right \rangle[/math] - скалярное произведение векторов a и b. Данный оператор используется для проецирования вектора a коллинеарно вектору b.
[math]\mathbf{b}_1 = \mathbf{a}_1[/math]
[math]\mathbf{b}_2 = \mathbf{a}_2 - \mathbf{proj}_{\mathbf{b}_1} \mathbf{a}_2 [/math]
[math]\mathbf{b}_3 = \mathbf{a}_3 - \mathbf{proj}_{\mathbf{b}_1} \mathbf{a}_3 - \mathbf{proj}_{\mathbf{b}_2} \mathbf{a}_3 [/math]
[math]\mathbf{b}_4 = \mathbf{a}_4 - \mathbf{proj}_{\mathbf{b}_1} \mathbf{a}_4 - \mathbf{proj}_{\mathbf{b}_2} \mathbf{a}_4 - \mathbf{proj}_{\mathbf{b}_3} \mathbf{a}_4[/math]
[math]\vdots[/math]
[math]\mathbf{b}_N = \mathbf{a}_N - \sum_{j=1}^{N-1} \mathbf{proj}_{\mathbf{b}_j} \mathbf{a}_N[/math]
На основе каждого вектора [math]\mathbf{b}_j (j = 1 \cdots N)[/math] может быть получен нормированный вектор [math]\mathbf{e}_j = \frac{\mathbf{b}_j}{||\mathbf{b}_j||}[/math].
Результаты процесса ортогонализации Грама-Шмидта: [math]\mathbf{b}_1\cdots\mathbf{b}_N[/math] - система ортогональных векторов либо система ортонормированных векторов [math]\mathbf{e}_1\cdots\mathbf{e}_N[/math].
Подробное описание можно найти в [5] [6].
1.3 Вычислительное ядро алгоритма
Основная вычислительная сложность данного алгоритма заключается в вычислении суммы проекций векторов [math]\sum_{j=1}^{i-1} \mathbf{proj}_{\mathbf{b}_j} \mathbf{a}_N[/math] для каждого [math]b_i[/math].
Подробное описание общей сложности алгоритма представлено в разделе #Последовательная сложность алгоритма
1.4 Макроструктура алгоритма
Из раздела описания, представленного в разделе #Математическое описание алгоритма можно сделать вывод, что основную часть метода составляют вычисления проекций векторов и сумм этих проекций.
1.5 Схема реализации последовательного алгоритма
Пусть дано множество векторов [math]{a_1, a_2, \cdots , a_N}[/math]. Тогда:
Шаг 1: [math]\mathbf{b}_1 = \mathbf{a}_1[/math]
Для вычисления остальных [math]\mathbf{b}_i [/math] следует последовательно выполнять Шаг 2 для всех [math] i [/math] таких что при [math] 2 \le i \le N [/math].
Шаг 2: [math]\mathbf{b}_i = \mathbf{a}_i - \sum_{j=1}^{i-1} \mathbf{proj}_{\mathbf{b}_{j}} \mathbf{a}_i[/math], где для каждого [math]\mathbf{proj}_b a[/math] выполняется [math]\frac{\left \langle a,b \right \rangle}{\left \langle b,b \right \rangle}b[/math].
1.6 Последовательная сложность алгоритма
Оценим сложность последовательной реализации алгоритма.
Для вычисления одной операции проекции требуется:
[math]{\color{Blue}\bullet} \quad M \cdot 3[/math] - умножений;
[math]{\color{Blue}\bullet} \quad (M-1) \cdot 2[/math] - сложений;
[math]{\color{Blue}\bullet} \quad 1 [/math] - деление;
Общее число операций, необходимых для вычисления алгоритма следующее:
[math]{\color{Blue}\bullet} \quad \frac{N(N-1)}{2} \cdot M \cdot 3[/math] - умножений;
[math]{\color{Blue}\bullet} \quad \frac{N(N-1)}{2} \cdot (M-1) \cdot 3 [/math] - сложений и вычитаний;
[math]{\color{Blue}\bullet} \quad \frac{N(N-1)}{2} [/math] - делений.
Исходя из вышесказанного, получаем итоговую сложность последовательного алгоритма [math]O(N^2 \cdot M)[/math].
1.7 Информационный граф
Опишем информационный граф алгоритма. В разделе #Схематичное описание информационного графа будет представлена схема с сокрытием подробностей реализации функции проекции. Такое ограничение позволяет продемонстрировать основные зависимости, не перегружая картинку лишними элементами. В разделе #Уточнение информационного графа будет раскрыта внутренняя структура операции проекции.
1.7.1 Схематичное описание информационного графа
Для вычисления [math]b_i[/math] требуется найти [math]proj_{b_{j}}a_{i}[/math] для всех [math]j \in [1, i][/math]. Следовательно для полного вычисления вектора [math]b_{i}[/math] требуется знать все [math]b_{j}[/math] с меньшим индексом. Такая зависимость по данным очень не удачна для параллелизма. Однако, если разбить процесс вычисления [math]b_{i}[/math] на несколько этапов, соответствующих функциям проекции ([math]proj_{b_{j}}a_{i}[/math]), то это позволит производить некоторые предварительные вычисления для [math]b_{j}[/math] до момента, когда станут известны все предшествующие ей [math]b_{i}[/math].
На Рис. 1 изображена зависимость каждого из этапов от предыдущих вычислений. Первая строка овалов содержит все проекции, зависящие от [math]b_{1}[/math], вторая строка - все проекции, зависящие от [math]b_{2}[/math], [math]\;\cdots\;[/math] самая верхняя строка содержит проекцию, которая зависит от [math]b_{N-1}[/math].
Рис. 2 показывает зависимость по данным немного в другом формате. Каждая строка представляет собой набор данных, которые требуются для вычисления [math]b_{i}[/math] из первого столбца. Второй и последующие столбцы группируют проекции, зависящие от одной из [math]b_{i}[/math] и, с помощью стрелки, показывают эту зависимость.
1.7.2 Уточнение информационного графа
В этом разделе будет представлено уточнение информационной структуры графа изображенного на Рис. 1 с помощью раскрытия внутренней структуры операции проекции. Полученный граф несёт в себе новую информацию, которую можно использовать для улучшения алгоритма, но становится менее понятным.
На Рис. 3 изображен граф, аналогичный графу Рис. 1, но функция [math]proj[/math] изображена более подробно. Оранжевым овалом обозначаются скалярные произведения двух векторов. Голубой треугольник - принимает на вход два числа (являющиеся результатами скалярных произведений) и вычисляет их частное. Зелёным шестиугольником обозначается произведение вектора [math]a_i[/math] на число. Фиолетовый ромб складывает два числа, от которых имеется зависимость.
По сравнению с Рис. 1, более детальное описание информационной модели позволяет заметить повторяющиеся построчно скалярные произведения [math](b_i, b_i)[/math]. Это замечание позволяет уменьшить суммарное число операций, но увеличивает связность графа.
1.8 Ресурс параллелизма алгоритма
На основе графа зависимости по данным, представленном на Рис. 1 или Рис. 2 раздела #Информационный граф можно выделить следующий ресурс для параллелизации. После вычисления [math]b_{i}[/math] можно начинать параллельно вычислять все [math]proj_{b_i}a_{j}[/math]. На Рис. 2 блоки, которые можно вычислять параллельно представлены столбцами.
Кроме того, при очень больших размерностях можно также распараллелить вычисление проекций. Так как в функции проекции участвует пара векторов [math]a[/math] и [math]b[/math], а сама операция проекции определяется как [math]\mathbf{proj}_b a = \frac{\left \langle a,b \right \rangle}{\left \langle b,b \right \rangle }b[/math], где [math]\left \langle a,b \right \rangle[/math] - скалярное произведение векторов [math]a[/math] и [math]b[/math], то с помощью использования векторных операций можно получить некоторый выигрыш в производительности.
1.8.1 Оценка ускорения параллельной реализации алгоритма
Рассмотрим случай, когда присутствует возможность запуска бесконечного числа потоков. Попробуем оценить сложность вычисления ортогональных векторов, при условии, что параллельно будут запускаться все проекции использующие одинаковую [math]b_i[/math]. В этом случае, останется критический путь зависимости по данным: [math]b_1 \leftarrow proj_{b_2}a_{2} \leftarrow b_2 \leftarrow proj_{b_3}a_{3} \leftarrow \ldots \leftarrow proj_{b_{N-1}}a_{N-1} \leftarrow b_N[/math]. Следовательно для оценки сложность распараллеленной версии алгоритма достаточно почитать сложность представленной цепочки вычислений.
[math]{\color{Blue}\bullet} \quad M \cdot 3[/math] умножения в одной проекции;
[math]{\color{Blue}\bullet} \quad (M-1) \cdot 2[/math] сложения;
[math]{\color{Blue}\bullet} \quad 1[/math] деление;
В цепочке зависимости данных присутствует [math]N-1[/math] проекция, значит всего:
[math]{\color{Blue}\bullet} \quad N \cdot M \cdot 3[/math] - умножений;
[math]{\color{Blue}\bullet} \quad N \cdot (M-1) \cdot 2[/math] - сложений;
[math]{\color{Blue}\bullet} \quad N \cdot 1[/math] - делений.
Получаем итоговую сложность [math]O(N \cdot M)[/math].
Если дополнительно рассматривать возможность параллелизации векторных операций (скалярного произведения и умножения вектора на число), критический путь можно сократить до следующих операций:
[math]{\color{Blue}\bullet} \quad 3*N [/math] умножений;
[math]{\color{Blue}\bullet} \quad N [/math] делений;
[math]{\color{Blue}\bullet} \quad N*logM [/math] сложений;
Следовательно оценка сложности параллельного алгоритма станет равна [math]O(N logM)[/math].
1.9 Входные и выходные данные алгоритма
Входные данные: квадратная матрица A (элементы [math]\mathbf{a}_{ij}[/math]). Матрица может быть приведена к верхнему(нижнему) треугольному виду.
Ограничение: матрица A должна быть невырожденной.
Объём входных данных: [math]\frac{n(n+1)}{2}[/math] (после преобразования матрицы к треугольному виду достаточно хранить только ненулевые элементы).
Выходные данные: матрица B (система ортогональных векторов [math]\mathbf{b}_{1} \cdots \mathbf{b}_{N}[/math]).
Объём выходных данных: [math]\frac{n(n+1)}{2}[/math] (в силу треугольного вида матрицы).
1.10 Свойства алгоритма
Соотношение последовательной и параллельной сложности в случае неограниченных ресурсов является линейной, так как последовательная реализация имеет сложность [math]O(n^3)[/math], а параллельная реализация - [math]O(n^2)[/math]).
Отношение числа операций к суммарному объему входных и выходных данных - линейное.
В силу теоремы о существовании ортонормированного базиса в результате работы алгоритма получаем систему из, минимум, одного линейно независимого ортогонального(ортонормированного) вектора. Из-за различных способов приведения исходной матрицы к треугольному виду(в случае использования разных методов или библиотек), а также ошибок округления в самом начале алгоритма ортогонализации могут получаться различные ортогональные(ортонормированные) системы в итоге.
2 ЧАСТЬ. Программная реализация алгоритма
2.1 Особенности реализации последовательного алгоритма
2.2 Локальность данных и вычислений
2.3 Возможные способы и особенности параллельной реализации алгоритма
2.4 Масштабируемость алгоритма и его реализации
TODO: (до 15 ноября)
2.5 Динамические характеристики и эффективность реализации алгоритма
2.6 Выводы для классов архитектур
2.7 Существующие реализации алгоритма
Пример реализации базового алгоритма с использованием технологий OpenMP и MPI представлен в статье [7].
Существует реализация процесса отртогонализации Грама-Шмидта, которая является вычислительно более устойчивой[8], такая реализации имеет название Модифицированный процесс Грама-Шмидта[9].
3 Литература
- ↑ А.А. Самарский, А.В. Гулин, Численные методы
- ↑ А.Е. Карелин, А.А. Светлаков, Использование ортогонализации Грама-Шмидта для повышения экономичности многоточечных алгоритмов рекуррентного оценивания параметров моделей объектов управления
- ↑ А.В. Пискова, А. А. Менщиков, А.Г. Коробейников, Использование ортогонализации Грама-Шмидта в алгоритме приведения базиса решётки для протоколов безопасности
- ↑ МГУ ВМК, Ф.Ф. Дедус, Л.И. Куликова, А.Н. Панкратов, Р.К. Тетуев, Классические ортогональные базисы в задачах аналитического описания и обработки информационных сигналов
- ↑ А.Н. Канатников, А.П.Крищенко, Линейная Алгебра
- ↑ The University of Southern Mississippi, Department of Mathematics, Jim Lambers, Gram-Schmidt Orthogonalization
- ↑ Tirana University, Department of Mathematics, Genci Berati, Gram–Schmidt Process in Different Parallel Platforms (Control Flow versus Data Flow)
- ↑ Massachusetts Institute of Technology, Per-Olof Persson, Introduction to Numerical Methods: Gram-Schmidt Orthogonalization
- ↑ С.А. Белов, Н.Ю. Золотых, Численные Методы Линейной Алгебры
