2  Схемы параллельных алгоритмов задач

В данном разделе приводятся примеры параллельных алгоритмов решения следующих задач: умножения матрицы на матрицу, задача Дирихле, решение систем линейных уравнений (СЛАУ) методом Гаусса и методом простой итерации. Здесь рассматривается простой вариант сеточной задачи (задача Дирихле), когда шаг сетки в пространстве вычислений одинаков и не меняется в процессе вычислений. При динамически изменяющемся шаге сетки потребовалось бы решать такую задачу параллельного программирования, как перебалансировка вычислительного пространства между компьютерами, для выравнивания вычислительной нагрузки компьютеров, а эта задача здесь не рассматривается.

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

Рассматриваемые задачи распараллеливаются крупнозернистыми методами. Для представления алгоритмов используется SPMD - модель вычислений (распараллеливание по данным). Однородное распределение данных по компьютерам – основа для хорошего баланса времени, затрачиваемого на вычисления, и времени, затрачиваемого на взаимодействия ветвей параллельной программы. При таком распределении преследуется цель: равенство объёмов распределяемых частей данных и соответствие нумерации распределяемых частей данных нумерации компьютеров в системе. Исходными данными рассматриваемых здесь алгоритмов являются матрицы, векторы и 2D (двумерное) пространство вычислений. В этих алгоритмах применяются следующие способы однородного распределения данных: горизонтальными полосами, вертикальными полосами и циклическими горизонтальными полосами. При распределении горизонтальными полосами матрица, вектор или 2D пространство "разрезается" на полосы по строкам (далее слово "разрезанная" будем писать без кавычек и матрицу, вектор или 2D пространство обозначать для краткости словом - данные). Пусть M – количество строк матрицы, количество элементов вектора или количество строк узлов 2D пространства, P – количество виртуальных компьютеров в системе, С1 = М /Р – целая часть от деления, С2 = М %Р – дробная часть. Данные разрезаются на Р полос. Первые (Р–С2) полос имеют по С1 строки, а остальные С2 полосы имеют по С1+1 строки. Полосы данных распределяются по компьютерам следующим образом. Первая полоса помещается в компьютер с номером 0, вторая полоса – в компьютер 1, и т. д. Такое распределение полос по компьютерам учитывается в параллельном алгоритме. Распределение вертикальными полосами аналогично предыдущему, только в распределении участвуют столбцы матрицы или столбцы узлов 2D пространства. И, наконец, распределение циклическими горизонтальными полосами. При таком распределении данные разрезаются на количество полос значительно большее, чем количество компьютеров. И чаще всего полоса состоит из одной строки. Первая полоса загружается в компьютер 0, вторая – в компьютер 1, и т.д., затем, Р-1-я полоса снова в компьютер 0, Р-я полоса в компьютер 1, и т.д.

Приведенные два алгоритма решения СЛАУ методом Гаусса показывают, что однородность распределения данных сама по себе еще недостаточна для эффективности алгоритма. Эффективность алгоритмов зависит еще и от способа распределения данных. Разный способ представления данных влечет, соответственно, и разную организацию алгоритмов, обрабатывающих эти данные.

 

2.1  Запуск параллельной программы

Под виртуальным компьютером понимается программно реализуемый компьютер. Виртуальный компьютер работает в режиме интерпретации его физическим процессором. В одном физическом компьютере, в общем случае, может находиться и работать одновременно виртуальных компьютеров - столько, сколько позволяет память физического компьютера. На системе МВС1000 в одном физическом компьютере создается только один виртуальный. Под виртуальной топологией здесь понимается программно реализуемая топология связей между виртуальными компьютерами на физической системе.

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

Запуск параллельной программы продемонстрируем на примере. Допустим, требуется решить задачу program.c. Алгоритм задачи распараллелен на N процессов, независимо выполняющихся и взаимодействующих друг с другом. Задана нужная для решения задачи топология связей между этими процессами: - top (например, двумерная решетка). Для решения этой задачи было бы оптимально иметь вычислительную систему из N компьютеров (Для МВС1000 должно быть N), с той же структурой связей, что и top. Далее в каждый компьютер необходимо загрузить по одному исполняемому модулю, реализующему ветвь параллельной программы, и стартовать эти модули. Ветви параллельной программы могут реализовываться копиями одной и той же программы (для МВС1000), а могут реализовываться разными программами (в общем случае). Опции и подробности загрузки нужно смотреть в соответствующих инструкциях.

Программа предварительно компилируется:

mpicc  [  ] -o program.exe program.c

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

Здесь рассматривается команда запуска параллельной программы на системе МВС1000. Предполагается, что ветви параллельной программы реализуются копиями одной и той же программы. Необходимое количество физических компьютеров и виртуальных компьютеров задаются пользователем в командной строке:

mpirun  -np N program.exe

N = {1,2,3,…} - указывает количество виртуальных компьютеров, необходимых для решения рассматриваемой программы с именем - program.exe. По этой команде система MPI создает (в оперативной памяти системы из N физических компьютеров) N виртуальных компьютеров, объединенных виртуальными каналами связи со структурой полный граф. И этой группе виртуальных компьютеров присваивается стандартное системное имя  MPI_COMM_WORLD. После чего пользовательская программа program.exe загружается в память каждого из созданных виртуальных компьютеров и стартует.

Отображение виртуальных компьютеров и структуры их связи на конкретную физическую систему осуществляется системой MPI автоматически, т.е. пользователю не нужно переделывать свою программу для разных физических систем (с другими компьютерами и другой архитектурой). (Рассматриваемая версия MPI не позволяет пользователю осуществлять это отображение, либо осуществлять пересылку виртуальных компьютеров в другие физические компьютеры, т.е. не позволяет перераспределять виртуальные компьютеры по физическим компьютерам).

2.2  Умножение матрицы на матрицу

Умножение матрицы на вектор и матрицы на матрицу являются базовыми макрооперациями для многих задач линейной алгебры, например итерационных методов решения систем линейных уравнений и т. п. Поэтому приведенные алгоритмы можно рассматривать как фрагменты в алгоритмах этих методов. В этой секции приведено три алгоритма умножения матрицы на матрицу. Разнообразие вариантов алгоритмов проистекает от разнообразия вычислительных систем и размеров задач. Рассматриваются и разные варианты загрузки данных в систему: загрузка данных через один компьютер; и загрузка данных непосредственно каждым компьютером с дисковой памяти. Если загрузка данных осуществляется через один компьютер, то данные считываются этим компьютером с дисковой памяти, разрезаются и части рассылаются по остальным компьютерам. Но данные могут быть подготовлены и заранее, т.е. заранее разрезаны по частям и каждая часть записана на диск в виде отдельного файла со своим именем; затем каждый компьютер непосредственно считывает с диска, предназначенный для него файл.

2.2.1  Алгоритм 1

Заданы две исходные матрицы A и B. Вычисляется произведение C = А х B, где А - матрица n₁ х n₂, и B - матрица n₂ х n₃. Матрица результатов C имеет размер n₁ х n₃. Исходные матрицы предварительно разрезаны на полосы, полосы записаны на дисковую память отдельными файлами со своими именами и доступны всем компьютерам. Матрица результатов возвращается в нулевой процесс.

Реализация алгоритма выполняется на кольце из p₁ компьютеров. Матрицы разрезаны как показано на рисунке 2.1: матрица А разрезана на p₁ горизонтальных полос, матрица B разрезана на p₁ вертикальных полос, и матрица результата C разрезана на p₁ полосы. Здесь предполагается, что в память каждого компьютера загружается и может находиться только одна полоса матрицы А и одна полоса матрицы B.

 

 

 

 

 

 

 

 

 

                 A               B             C

 

Рис. 2.1  Разрезание данных для параллельного алгоритма произведения двух матриц при вычислении на кольце компьютеров. Выделенные полосы расположены в одном компьютере.

 

 

Поскольку по условию в компьютерах находится по одной полосе матриц, то полосы матрицы B (либо полосы матрицы A) необходимо "прокрутить" по кольцу компьютеров мимо полос матрицы A (матрицы B). Каждый сдвиг полос вдоль кольца и соответствующее умножение представлено на рисунке 2.2 в виде отдельного шага. На каждом таком шаге вычисляется только часть полосы. Процесс i вычисляет на j-м шаге произведение i-й горизонтальной полосы матрицы A и j-й вертикальной полосы матрицы B, произведение получено в подматрице (i,j) матрицы C. Текст программы, реализующий алгоритм, приведен в главе 9.

Последовательные стадии вычислений иллюстрируются на рисунке 2.2.

·         1. Каждый компьютер считывает с диска соответствующую ему полосу матрицы А. Нулевая полоса должна считываться нулевым компьютером, первая полоса - первым компьютером и т.д., последняя полоса - считывается последним компьютером. На рисунке 2.2 полосы матрицы А и B пронумерованы.

·         2. Каждый компьютер считывает с диска соответствующую ему полосу матрицы B. В данном случае нулевая полоса должна считываться нулевым компьютером, первая полоса - первым компьютером и т.д., последняя полоса - считывается последним компьютером.

·         3. Вычислительный шаг 1. Каждый процесс вычисляет одну подматрицу произведения. Вертикальные полосы матрицы B сдвигаются вдоль кольца компьютеров.

·         4. Вычислительный шаг 2. Каждый процесс вычисляет одну подматрицу произведения. Вертикальные полосы матрицы B сдвигаются вдоль кольца компьютеров. И т. д.

·         5. Вычислительный шаг p1-1. Каждый процесс вычисляет одну подматрицу произведения. Вертикальные полосы матрицы B сдвигаются вдоль кольца компьютеров.

·         6. Вычислительный шаг p1. Каждый процесс вычисляет одну подматрицу произведения. Вертикальные полосы матрицы B сдвигаются вдоль кольца компьютеров.

·         7. Матрица C собирается в нулевом компьютере.

 

 

 

 

 

 

 

 

 

 

 

 

              1. scatter A                          2. scatter B

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

           7. Сбор результатов в С

 

Рис. 2.2  Стадии вычислений произведения матриц в кольце компьютеров.

 

Если "прокручивать" вертикальные полосы матрицы B, то матрица С будет распределена горизонтальными полосами, а если "прокручивать" горизонтальные полосы матрицы A, то матрица С будет распределена вертикальными полосами.

Алгоритм характерен тем, что после каждого шага вычислений осуществляется обмен данными. Пусть t_u, t_u, и t_p время операций, соответственно, умножения, сложения и пересылки одного числа в соседний компьютер. Согласно приведенным в начале секции обозначениям суммарное время операций умножений равно:

U = (n₁*n₂)*(n₃*n₂)*t_u,

суммарное время операций сложений равно:

S = (n₁*n₂)*(n₃*(n₂-1))*t_s,

суммарное время операций пересылок данных по всем компьютерам равно:

P = (n₃*n₂)*(p₁-1)*t_p.

Тогда общее время вычислений будет равно:

T = (U+S+P)/p₁

И отношение времени "вычислений без обменов" к общему времени вычислений есть величина:

K = (U+S)/(U+S+P).

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

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

2.2.2  Алгоритм 2

Вычисляется произведение C = А х B, где А - матрица n₁ х n₂, и B - матрица n₂ х n₃. Матрица результатов C имеет размер n₁ х n₃. Исходные матрицы первоначально доступны на нулевом процессе, и матрица результатов возвращена в нулевой процесс.

Параллельное выполнение алгоритма осуществляется на двумерной (2D) решетке компьютеров размером p₁ х p₂. Матрицы разрезаны, как показано на рисунке 2.3: матрица А разрезана на p₁ горизонтальных полос, матрица B разрезана на p₂ вертикальных полос, и матрица результата C разрезана на p₁ х p₂ подматрицы (или субматрицы).

 

 

 

 

 

 

 

 

                 A            B         C

 

Рис. 2.3  Разрезание данных для параллельного алгоритма произведения двух матриц при вычислении в 2D решетке компьютеров. Выделенные данные расположены в одном компьютере.

 

 

Каждый компьютер (i,j) вычисляет произведение i-й горизонтальной полосы матрицы A и j-й вертикальной полосы матрицы B, произведение получено в подматрице (i,j) матрицы C.

Последовательные стадии вычисления иллюстрируются на рисунке 2.4:

·         1. Матрица А распределяется по горизонтальным полосам вдоль координаты (x,0).

·         2. Матрица B распределяется по вертикальным полосам вдоль координаты (0,y).

·         3. Полосы А распространяются в измерении  y.

·         4. Полосы B распространяются в измерении  х.

·         5. Каждый процесс вычисляет одну подматрицу произведения.

·         7. Матрица C собирается из (x,y) плоскости.

 

Осуществлять пересылки между компьютерами во время вычислений не нужно, т. к. все полосы матрицы А пересекаются со всеми полосами матрицы B в памяти компьютеров системы.

 

 

              2.scatter B     координаты

 

 

 

 

 

 

 

 

                                             3.broadcast

                                               подматриц A

 

 

 

 

 

 

 

 

 

 

               4.broadcast                    5.вычисление

                 подматриц B                    произведений

                                                (подматриц в C)

 

 

 

 

 

 

 

 

 

               6.сбор

                 результатов в C

 

 

Рис. 2.4  Стадии вычисления произведения матриц в 2D параллельном алгоритме.

 

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

T = (U+S)/(p₁*p₂)

И отношение времени "вычислений без обменов" к общему времени вычислений есть величина:

K = (U+S)/(U+S)=1.

Здесь значения T, K, U  и S - смотри в предыдущем разделе.

2.2.3  Алгоритм 3

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

В приведенном ниже алгоритме осуществляется отображение основных данных объемом  n₁ х n₂ + n₂ х n₃ + n₁ х n₃ на объемную сетку компьютеров размером p₁ х p₂ х p₃. Матрицы разрезаны, как показано на рисунке 2.5: матрица А разрезана на p₁ х p₂ субматрицы, матрица B разрезана на p₂ х p₃ субматрицы, и матрица C разрезана на p₁ х p₃ субматрицы. Компьютер (i,j,k) вычисляет произведение субматрицы (i,j) матрицы А и субматрицы (j,k) матрицы B. Субматрица (i,k) матрицы C получается суммированием промежуточных результатов, вычисленных в компьютерах (i,j,k), j = 0,…,p2-1.

Последовательные стадии вычисления иллюстрируются на рисунке 2.6.

·         1. Субматрицы А распределяются в (x,y,0) плоскости.

·         2. Субматрицы B распределяются в (0,y,z) плоскости.

·         3. Субматрицы А распространяются в измерении z.

·         4. Субматрицы B распространяются в измерении х.

·         5. Каждый процесс вычисляет одну субматрицу.

·         6. Промежуточные результаты редуцируется в измерении y.

·         7. Матрица C собирается из (x,0,z) плоскости.

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

Этот алгоритм для больших матриц является еще более эффективным, чем предыдущий. Соотношения для общего времени вычислений в этом алгоритме будет равно:

T = (U+S)/(p₁*p₂*p₃)

А отношение времени "вычислений без обменов" к общему времени вычислений есть величина:

K = (U+S)/(U+S)=1.

Здесь значения T, K, U  и S - смотри в предыдущем разделе.

 

 

 

 

 

 

 

 

                   A            B         C

 

Рис. 2.5  Разрезание данных для параллельного алгоритма произведения двух матриц при вычислении в 3D решетке компьютеров.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

        2. scatter B                            3. broadcast

                                                   подматриц из А

 

 

 

 

 

 

 

 

 

 

 

 

 

 

        4. broadcast                            5. Вычисление

           подматриц из В                          произведений

                                                   (подматриц в С)

 

 

 

 

 

 

 

 

 

 

 

 

 

        6. reduce произведений                  7. gather C

           (суммирование произведений)             (сбор результатов)

 

Рис. 2.6  Стадии вычислений в 3D параллельном алгоритме произведения матриц.

 

 

 

 

 

 

 

 

 

 

2.3  Задача Дирихле. Явная разностная схема для уравнения Пуассона

В прямоугольной области  0 ≤ x ≤ a,  0 ≤ y ≤ b требуется найти решение уравнения:

при заданных значениях функции u на границе.

Явная разностная схема решения этого уравнения имеет вид:

где  и  - значения функций   и  в точке  разностной сетки.

Ниже представлен главный фрагмент программы итерационного решения задачи. Фрагмент последовательной программы.

 

double A[n+2][m+2], B[n][m];

...

/* Главный цикл */

while(! converged) {

/* выполнение схемы "крест" */

for(j = 1; j <= m; j++)

  for(i = 1; i <= n; i++)

    B[i-1][j-1] = 0.25*(A[i-1][j]+A[i+1][j]+A[i][j-1]+A[i][j+1]);

/* копирование результата обратно в массив A */

for(j = 1; j <= m; j++)

  for(i = 1; i <= n; i++)

    A[i][j] = B[i-1][j-1];

...

 }

 

Этот фрагмент программы описывает главный цикл итерационного процесса решения, где на каждой итерации значение в окрестности точки заменяется средним значением сумм значений ее четырех соседних точек на предыдущем временном шаге итерации (смотри рисунок 2.7).

 

 

 

 

 

 

 

 

 

 

Рис. 2.7  Вычисление точки через значения точек на предыдущем шаге итерации.

 

 

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

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

Способы разрезания данных могут быть разные, и в зависимости от способа разрезания данных будут и разные схемы взаимодействий между параллельными процессами. Распределение данных по процессорам должно быть сбалансировано. Связь между процессорами должна быть минимизирована.

Двумерный массив может быть разрезан на части как в одном измерении, так и в двух измерениях. Рисунок 2.8 поясняет эти способы разрезания данных: 1D разрезание, где матрица разрезана в одном измерении полосами, и 2D разрезание, где матрица разрезана в двух измерениях.

 

 

 

                                 

 

 

                                   1D  разрезание                               2D  разрезание

 

Рис. 2.8  Два разных способа разрезания двумерного массива на блоки.

 

 

На рисунке матрица разрезана на четыре блока; каждый блок обрабатывается в отдельном процессоре. Так как связь между процессорами осуществляется границами блоков, то объем связи минимизирован в 2D разрезании, которое имеет меньший периметр области связи. В этом разбиении каждый процессор взаимодействует, в общем случае, с четырьмя соседями, быстрее, чем два соседа в 1D разрезании. Когда отношение n/P (P число процессоров) маленькое, время связи будет зависеть от внешних коммуникаций, и первое разбиение будет лучше по взаимодействиям. Когда это отношение большое, второе разбиение лучше по взаимодействиям. Здесь используется первое разбиение.

Значение каждой точки в массиве B вычисляется через значения четырех ее соседей в массиве A. Связь между процессорами необходима границами блоков, чтобы вычислять граничные точки блоков через свои соседние точки, которые принадлежат другому процессору. Так как точки вычисляются через только свои соседние точки, вычисленные на предыдущей итерации, то для вычисления граничных точек блока необходимо присутствие копии соответствующего столбца предыдущей итерации соседнего блока, находящегося в соседнем компьютере. Следовательно, при 1D разрезании необходимо распределение блоков массива A по компьютерам с перекрытием в один столбец. Таких перекрытий столбцами для массива B не нужно; - это следует из алгоритма. На рисунке 2.9 иллюстрируется перекрытие соседних полос в один столбец. На этом рисунке римскими цифрами обозначены крайние столбцы соседних полос. Пунктирной линией обозначена граница разреза массива данных.

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

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

                начальный разрез     перекрытие полос

                 массива данных       в один столбец

 

Рис. 2.9  Перекрытие двух полос в один столбец.

 

 

На рисунке 2.10 показаны дополнительные столбцы массива A и то, как данные пересылаются после каждой итерации.

 

 

                             Процесс 0                                 Процесс 1                                Процесс 2

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Рис. 2.10  Перекрытие 1D блоков в массиве A и пересылка столбцов в конце каждой итерации.

 

 

Перекрытие одним столбцом обычно используется для разностной схемы, представляемой пяти-точечным шаблоном сетки и использующей матрицу окрестности размером 3x3. При использовании матрицы большей размерности, например 5x5, необходимо делать 2-х столбцовое перекрытие полос массивов, находящихся в смежных компьютерах. В общем случае, если рассматривается окрестность с диаметром в k точек от вычисляемой точки на предыдущем шаге вычислений, то перекрытие массивов нужно делать в k столбцов.

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

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

2.4  Параллельные алгоритмы решения систем линейных алгебраических

       уравнений методом Гаусса

Требуется найти решение системы линейных алгебраических уравнений:

a₁₁x₁ + a₁₂x₂ +  … + a_1nx_n = f₁

a₂₁x₂ + a₂₂x₂ +  … + a_2nx_n = f₂

. . . . . . . . . . . . . .

a_n1x₁ + a_n2x₂ +  … + a_nnx_n = f_n

Метод Гаусса основан на последовательном исключении неизвестных.

Здесь рассматриваются два алгоритма решения СЛАУ методом Гаусса. Они связаны с разными способами представления данных (матрицы коэффициентов и правых частей) в распределенной памяти мультикомпьютера.

2.4.1  Первый алгоритм решения СЛАУ методом  Гаусса

В алгоритме, представленном в данной секции, исходная матрица коэффициентов A и вектор правых частей F разрезаны горизонтальными полосами, как показано на рисунке 2.9. Каждая полоса загружается в соответствующий компьютер: нулевая полоса – в нулевой компьютер, первая полоса – в первый компьютер, и т. д., последняя полоса – в p1 компьютер.

 

 

 

 

                  p₁

 

 

 

                            A      F

 

Рис. 2.11  Разрезание данных для параллельного алгоритма 1 решения СЛАУ методом Гаусса.

 

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

 

 

 

 

 

 

 

После прямого хода полосы матрицы A в каждом узле будут иметь вид:

 

 

Рис. 2.12  Вид полос после прямого хода в алгоритме 1 решения СЛАУ методом Гаусса.

 

Пример приведен для четырех узлов; $ – вещественные числа.

Аналогично, последовательно по компьютерам, начиная с последнего по номеру компьютера, осуществляется обратный ход.

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

2.4.2  Второй алгоритм решения СЛАУ методом Гаусса

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

 

 

 

 

 

 

 

 

 

 

                                A             F

 

Рис. 2.13  Разрезание данных для параллельного алгоритма 2 решения СЛАУ методом Гаусса.

 

Первая строка матрицы помещается в компьютер 0, вторая строка – в компьютер 1, и т. д., p1-1-я строка в узел p1 (где p1 количество узлов в системе). Затем, p1–я строка, снова помещается в узел 0, p1+1-я строка – в узел 1, и т. д.

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

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

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

 

 

Рис. 2.12  Вид полос после прямого хода в алгоритме 2 решения СЛАУ методом Гаусса.

 

 

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

2.5  Параллельный алгоритм решения СЛАУ методом простой итерации

Здесь рассматривается параллельный алгоритм решения СЛАУ методом простой итерации. Приближенные решения (итерации) системы линейных уравнений последовательно находятся по формуле:

Для решения этой задачи на параллельной системе исходную матрицу коэффициентов А разрезаем на p1 горизонтальных полосы по строкам, где p1 – количество компьютеров в системе. Аналогично, горизонтальными полосами разрезаются вектор b (правая часть) и вектора y₀ (начальное приближение), y_k (текущее приближение) и y_k+1 (следующее приближение). Полосы последовательно распределяются по соответствующим компьютерам системы, как и в описанном выше, первом алгоритме умножения матрицы на матрицу.

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

 

 

Это условие в некоторых компьютерах системы, в текущий момент, может выполняться, а в некоторых нет. Но условием завершения работы каждого компьютера является безусловное выполнение условия (2) во всех компьютерах. Таким образом, прежде чем завершить работу, при выполнении условия (2), каждый компьютер должен предварительно узнать: во всех ли компьютерах выполнилось условие (2)? И если условие (2) не выполнилось хотя бы в одном компьютере, то все компьютеры должны продолжить работу. Это обстоятельство связано с тем, что в операции умножения матрицы на вектор участвуют все компьютеры, взаимодействуя друг с другом. И цепочку этих взаимодействий прерывать нельзя, если хотя бы в одном из компьютеров не выполнится условие (2). При не выполнении условия 2, каждый процесс передает всем остальным процессам полученную итерацию своих корней. И тем самым вектор y полностью восстанавливается в каждом процессе для выполнения операции его умножения на матрицу коэффициентов на следующем шаге итерации.