Программирование для многопроцессорных систем в стандарте MPI - Шпаковский Г.И., Серикова Н.В

7.2.2. Способы создания коммуникаторов

Из клеточного алгоритма следует, что было бы удобным трактовать каждую строку и каждый столбец процессов как новый коммуникатор. Это можно организовать несколькими способами.

Способ 1. Для иллюстрации создадим коммуникатор, который составляет группа из процессов первой строки виртуальной сетки. Предположим, что MPI_COMM_WORLD состоит из p процессов, где q 2 = p. Предположим, что φ(x) = (x/q; x mod q). Тогда первая строка процессов состоит из процессов с номерами 0, 1,..., q-1.

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

MPI_Group MPI_GROUP_WORLD; MPI_Group first_row_group; MPI_Comm first_row_comm;

int row_size;

int* process_ranks;

/* Создание списка процессов в новом коммуникаторе */ process_ranks = (int*) malloc(q*sizeof(int));

for (proc = 0; proc < q; proc++) process_ranks[proc] = proc;

/* Получение группы MPI_COMM_WORLD */ MPI_Comm_group(MPI_COMM_WORLD, &MPI_GROUP_WORLD);

/* Создание новой группы */ MPI_Group_incl(MPI_GROUP_WORLD, q, process_ranks, &first_row_group);

/* Создание нового коммуникатора */ MPI_Comm_create(MPI_COMM_WORLD, first_row_group,&first_row_comm);

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

MPI_COMM_WORLD; создание группы – MPI_Group_incl; созда-

ние коммуникатора – MPI_Comm_create. Результатом является коммуникатор первой строки comm. Теперь процессы в первой строке comm могут совершать коллективные операции связи. Например, процесс 0 (в группе первой строки) может передавать A00 другим процессам своей группы. Фрагмент кода программы для выполнения коллективной операции представлен ниже.

181

int my_rank_in_first_row; float* A_00;

/* my_rank есть номер процесса в MPI_GROUP_WORLD */

if (my_rank < q)

{

MPI_Comm_rank( first_row_comm, &my_rank_in_first_row);

/* Выделяем память для A_00, order = n_bar */ A_00 = (float*) malloc (n_bar*n_bar*sizeof(float));

if (my_rank_in_first_row == 0)

{ … } /* Инициализация A_00 */ MPI_Bcast(A_00, n_bar*n_bar, MPI_FLOAT, 0, first_row_comm);

}

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

MPI_Comm my_row_comm; int my_row;

/* my_rank есть номер процесса в MPI_COMM_WORLD. q*q = p */ my_row = my_rank/q;

MPI_Comm_split( MPI_COMM_WORLD, my_row, my_rank, &my_row_comm);

Единственный вызов MPI_Comm_split создает q новых коммуникаторов, имеющих одно и то же имя my_row_comm.

Способ 2. В клеточном алгоритме удобно идентифицировать процессы в MPI_COMM_WORLD координатами квадратной сетки, то есть для каждой строки и каждого столбца сетки следует сформировать собственный коммуникатор. Для этого необходимо связать квадратную сетку с MPI_COMM_WORLD. Чтобы сделать это, нужно определить:

1.Число размерностей сетки (равно двум).

2.Размер каждого измерения (q строк и q столбцов).

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

182

int varying_coords[2]; MPI_Comm row_comm;

varying_coords[0] = 0; varying_coords[1] = 1; MPI_Cart_sub(grid_comm, varying_coords, row_comm);

Вызов MPI_Cart_sub создает q новых коммуникаторов. Переменная Varying_coords есть массив булевых элементов. Она определяет, принадлежит ли каждое измерение новому коммуникатору. Так как создаются коммуникаторы для строк сетки, каждый новый коммуникатор состоит из процессов, полученных при фиксировании номера строки и изменения номеров столбцов. В каждом процессе новый коммуникатор возвращается в row_comm. Чтобы создать коммуникаторы для столбцов, значения для varying_coords подаются в другом порядке. Вызов функции тогда будет таким.

MPI_Comm col_comm;

varying_coords[0] = 1; varying_coords[1] = 0; MPI_Cart_sub(grid_comm, varying_coord, col_comm);

Заметим схожесть функций MPI_Cart_sub и MPI_Comm_split.

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

7.2.3. Параллельная программа для клеточного алгоритма

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

Тогда подпрограмма создания коммуникаторов для решения задачи будет иметь следующий код:

void Setup_grid(GRID_INFO_TYPE* grid)

{ int old_rank, dimensions[2], periods[2], coordinates[2], varying_coords[2]; /* определение параметров сетки */

MPI_Comm_size(MPI_COMM_WORLD, &(grid->p)); MPI_Comm_rank(MPI_COMM_WORLD, &old_rank); grid->q = (int) sqrt((double) grid->p);

dimensions[0] = dimensions[1] = grid->q; periods[0] = periods[1] = 1;

/*создание коммуникатора с картезианской топологии размерности 2 */ MPI_Cart_create(MPI_COMM_WORLD, 2, dimensions, periods,1, &(grid->comm)); MPI_Comm_rank(grid->comm, &(grid->my_rank));

MPI_Cart_coords(grid->comm, grid->my_rank, 2,coordinates); grid->my_row = coordinates[0]; grid->my_col = coordinates[1];

/* создание коммуникатора строки */ varying_coords[0] = 0; varying_coords[1] = 1; MPI_Cart_sub(grid->comm, varying_coords,&(grid->row_comm));

/* создание коммуникатора столбца */ varying_coords[0] = 1; varying_coords[1] = 0; MPI_Cart_sub(grid->comm, varying_coords,&(grid->col_comm));

}

Ниже представлен текст обобщающей программы клеточного умножения матриц. LOCAL_MATRIX_TYPE – тип элементов умножаемых матриц, DERIVED_LOCAL_MATRIX – тип перемножаемых матриц.

void Mult_Matr (int n, GRID_INFO_TYPE* grid, LOCAL_MATRIX_TYPE* local_A, local_B, local_C)

source = (grid->my_row + 1) % grid->q;

dest = (grid->my_row + grid->q – 1) % grid->q;

/* выделение памяти для пересылки блоков A */ temp_A = Local_matrix_allocate(n_bar);

185

for (step = 0; step < grid->q; step++)

{

bcast_root = (grid->my_row + step) % grid->q; if (bcast_root == grid->my_col) {

MPI_Bcast(local_A,1,DERIVED_LOCAL_MATRIX,bcast_root, grid->row_comm);

/* умножение подматриц */ Local_matrix_multiply(local_A, local_B,local_C);

}

else { MPI_Bcast(temp_A,1,DERIVED_LOCAL_MATRIX,bcast_root,

grid->row_comm);

/* умножение подматриц */ Local_matrix_multiply(temp_A, local_B, local_C);

} /* пересылка подматриц В*/ MPI_Send(local_B, 1, DERIVED_LOCAL_MATRIX, dest, tag,

grid->col_comm); MPI_Recv(local_B,1,DERIVED_LOCAL_MATRIX,source,tag,

grid->col_comm, &status);

}

}

КОНТРОЛЬНЫЕ ВОПРОСЫ И ЗАДАНИЯ К ГЛАВЕ 7

Контрольные вопросы к 7.1

1.Объясните понятие “самопланирующий алгоритм”.

2.Можно ли в главном процессе программы умножения матрицы на вектор заменить коллективную операцию MPI_Bcast на обычный межпроцессный обмен MPI_Send/MPI_Recv? Как изменится код подчиненного процесса?

3.Какую функцию выполняют параметры status, tag при передаче данных в программе умножения матрицы на вектор?

4.Как изменится код программы умножения матрицы на вектор, если необходимо умножить вектор на матрицу?

5.Можно ли в программе умножения матриц использовать другие виды коммуникаций?

Контрольные вопросы к 7.2

1.В чем преимущества и недостатки двух способов создания коммуникаторов для алгоритма клеточного умножения?

2.Как определить тип DERIVED_LOCAL_MATRIX средствами MPI?

3.Чему равны значения source, tag в операциях обмена процедуры Mult_Matr?

4.Каково минимальное количество процессов для клеточного алгоритма?

5.Проведите анализ эффективности для алгоритма клеточного умножения.

186

Задания для самостоятельной работы

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

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

7.3.Напишите программу, реализующую алгоритм умножения матриц, при котором происходит равномерное распределение частей матрицы А по процессам,

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

7.4.Напишите программу, реализующую клеточный алгоритм умножения

матриц.

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

Глава 8. РЕШЕНИЕ ДИФФЕРЕНЦИАЛЬНЫХ УРАВНЕНИЙ В ЧАСТНЫХ ПРОИЗВОДНЫХ

Решение дифференциальных уравнений в частных производных является ядром многих приложений. Решение задачи на конечноразностной вычислительной сетке методом Якоби, рассмотренное в этом параграфе, дает возможность продемонстрировать использование «виртуальной топологии», которая делает более удобным размещение процессов на вычислительной среде, а также позволяет исследовать эффективность различных способов коммуникаций для решения одной и той же задачи. Метод Якоби был выбран благодаря простоте вычислительной части [4, 21].

8.1. ЗАДАЧА ПУАССОНА

Задача Пуассона, то есть решение дифференциальных уравнений в частных производных, выражается следующими уравнениями:

Будем считать, что область решения является квадратной. Чтобы найти приближенное решение, определим квадратную сетку, содержащую точки (xi, yi), задаваемые как

187

Таким образом, вдоль каждого края сетки имеется n+2 точки. Следует найти аппроксимацию для u(x,y) в точках (xi, yj) выбранной сетки. Обозначим через ui,j значения u в (xi, yj), через h расстояние между точками, равное 1/(n+1). Тогда формула 8.1. для каждой из точек будет выглядеть следующим образом:

Вычисляем значения ui,j в каждой точке сетки, которые замещают предыдущие значения, используя выражение

uik,+j 1 = (uik−1, j + uik, j+1 + uik, j−1 + uik=1, j − h 2fi, j ) / 4.

Этот процесс называется итерациями Якоби и повторяется до получения решения. Фрагмент программы для этого случая таков:

integer i, j, n

double precision, u (0:n+l, 0:n+l), unew(0:n+l, 0:n+l) do 10 j = l, n

do 10 i = l, n

unew (i, j) = 0.25*(u(i-l, j) + u(i, j+l) + u(i, j-l)+u(i+l, j)) – h * h * f(i, j)

10continue

8.2.ПАРАЛЛЕЛЬНЫЕ АЛГОРИТМЫ ДЛЯ МЕТОДА ИТЕРАЦИЙ ЯКОБИ

8.2.1. Параллельный алгоритм для 1D композиции

Чтобы распараллелить последовательный алгоритм, нужно параллелизовать цикл. Для этого необходимо распределить данные (в нашем случае массивы u, unew, f) по процессам. Одна из простейших декомпозиций состоит в следующем: физическая область разделяется на слои, каждый из них обрабатывается отдельным процессом. Эта декомпозиция может быть описана следующим образом:

integer i, j, n, s, e

double precision u(0:n+l, s:e), unew(0:n+l, s:e) do 10 j = s, e

do 10 i = l, n

unew(i, j) = 0.25*(u(i-l, j)+ u(i, j+l)+ u(i, j-l) + u(i+l,j)) – h * h * f(i,j) 10 continue

188

Здесь s,e указывают значения номеров строк слоя, за которые ответственен данный процесс.

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

Процесс 2

Процесс 1

Процесс 0

Рис. 8.1. Область с теневыми точками для процесса 1

Вследствие этого размерность массивов изменится: double precision u ( 0:n+1, s-1: e+1).

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

При одном измерении можно просто использовать номер в новом коммуникаторе плюс или минус 1, чтобы найти соседа и не использовать фунцию MPI_CART_CREATE, но даже в этом случае выбор может быть не наилучшим, поскольку соседи, определенные таким

189

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

Каждый процесс посылает и получает сообщения от соседей. В декомпозиции 1D это соседи, расположенные выше и ниже. Способ определения соседей заключается в следующем: копия верхней строки одного процесса является дном теневой строки процесса, который расположен выше. Следовательно, каждый процесс является как принимающим, так и посылающим данные. Фактически данные сдвигаются от одного процесса к другому, и MPI имеет процедуру MPI_CART_SHIFT, которая может быть использована, чтобы найти соседей. Опишем процедуру MPE_DECOMP1D, которая определяет декомпозицию массивов и вызывается функцией

call MPE_DECOMP1D(n, nprocs, myrank, s, e),

где nprocs – число процессов в картезианской топологии, myrank – координаты процесса, n – размерность массива. MPE_Decomp1d вычисляет значения s и e. Например следующим образом.

Если n делится на nprocs без остатка, то

s = 1 + myrank * (n/nprocs) e = s + (n/nprocs) – 1.

Когда деление без остатка не получается, если floor(x) – целая часть х плюс 1, то наиболее простой выход:

s = 1 + myrank * floor (n/nprocs)

if (myrank .eq. nprocs – 1) then e = n else e = s + floor (n/nprocs) – 1

endif.

Для каждого процесса мы должны получить теневые данные для строки с номером s–1 от процесса ниже и данные для строки с номером e+1 от процесса выше. Ниже представлена процедура обмена данными, для которой определены следующие параметры: а – массив, nx – количество пересылаемых данных строки, s – номер первой строки массива в данном процессе, e - номер последней строки массива данного процесса, nbrbottom – номер процесса, расположенного ниже данного, а nbrtop – номер процесса, расположенного выше данного.

subroutine EXCHG1(a, nx, s, e, comm1d, nbrbottom, nbrtop ) use mpi

integer nx, s, e

double precision a(0:nx+l, s-l:e+l)

190