лекции (MPI), страница 8

длявысокопроизводительных систем",24.11.201522MPI_CART_SUB примерlastrow = nrow - 1;MPI_Gather(aij, 1, MPI_INT, alocal, 1, MPI_INT, lastrow, commcol);MPI_Barrier(MPI_COMM_WORLD);printf("%6d|%6d|%6d %6d|%6d|%8d|",Iam,id2D,coords2D[0],coords2D[1],colID,coords1D[0]);for (i=0; i<=lastrow; i++) {printf("%6d ",alocal[i]);}printf("\n");MPI_Finalize();/* let MPI finish up ... */}Лекции спецкурс "Пар. прогр. длявысокопроизводительных систем",24.11.201523MPI_CART_GETПараметрыsubgrid_commMPI_CommInputКоммуникаторndimsintInputЧисло измеренийdimsint *OutputМассив длин по каждому изизмеренийperiodsint *OutputПериодичность понаправлениямcoordsint *OutputКоординаты вызывающегопроцесса в решетке (массивразмера ndims )Лекции спецкурс "Пар.

прогр. длявысокопроизводительных систем",24.11.201524MPI_CART_GETИспользуется для получения информации о параметрахдекартовой топологии для заданного коммуникатораint MPI_Cart_sub( MPI_Comm comm, int maxdims,int *dims, int *periods, int *coords)Лекции спецкурс "Пар. прогр. длявысокопроизводительных систем",24.11.201525MPI_CART_GET пример/* create Cartesian topology for processes */dims[0] = nrow;dims[1] = mcol;MPI_Cart_create(MPI_COMM_WORLD, ndim, dims, period,reorder, &grid_comm);MPI_Comm_rank(grid_comm, &me);MPI_Cart_coords(grid_comm, me, ndim, coords);/* create row subgrids */belongs[0] = 1;belongs[1] = 0;MPI_Cart_sub(grid_comm, belongs, &row_comm);/* Retrieve subgrid dimensions and other info */MPI_Cartdim_get(row_comm, &mdims);MPI_Cart_get(row_comm, mdims, dims, period, row_coords);Лекции спецкурс "Пар. прогр.

длявысокопроизводительных систем",24.11.201526MPI_CART_SHIFTПолучение номеров процессов, находящихся назаданном расстоянии displacement от текущего процессав заданном измерении direction декартовой решетки.int MPI_Cart_shift( MPI_Comm comm, intdirection, int displ, int *source, int *dest )Лекции спецкурс "Пар. прогр. длявысокопроизводительных систем",24.11.201527MPI_CART_SHIFTПараметрыcommMPI_CommInputКоммуникаторdirectionintInputРазмерность, по которой будетпроизводиться сдвигdisplintInputВеличина и направлениесдвига (<0; >0; or 0)sourceint *OutputПроцесс- источникdestint *OutputПроцесс- получательЛекции спецкурс "Пар.

прогр. длявысокопроизводительных систем",24.11.201528MPI_CART_SHIFTДля периодических измерений осуществляетсяциклический сдвиг, для непериодических – линейныйсдвиг.Для n-мерной декартовой решетки значение directionдолжно быть в пределах от 0 до n-1.Значения source и dest можно использовать,например, для обмена функцией MPI_Sendrecv.Не является коллективной операцией!Лекции спецкурс "Пар. прогр. длявысокопроизводительных систем",24.11.201529Пример декартовой решетки(send&recv, mesh)MPI_Request reqs[8];MPI_Status stats[8];MPI_Comm cartcomm;int dims[2],periods[2], reorder, coords[2];MPI_Init(&argc,&argv);MPI_Comm_size(MPI_COMM_WORLD, &numtasks);if (numtasks == SIZE) {// define dims, periods, reorderMPI_Cart_create(MPI_COMM_WORLD, 2, dims, periods, reorder, &cartcomm);MPI_Comm_rank(cartcomm, &rank);MPI_Cart_coords(cartcomm, rank, 2, coords);MPI_Cart_shift(cartcomm, 1, 1, &nbrs[UP], &nbrs[DOWN]);MPI_Cart_shift(cartcomm, 0, 1, &nbrs[LEFT], &nbrs[RIGHT]);outbuf = rank;Лекции спецкурс "Пар.

прогр. длявысокопроизводительных систем",24.11.201530Пример декартовой решетки(send&recv, mesh)for (i=0; i<4; i++) {dest = nbrs[i];source = nbrs[i];MPI_Isend(&outbuf, 1, MPI_INT, dest, tag, MPI_COMM_WORLD, &reqs[i]);MPI_Irecv(&inbuf[i], 1, MPI_INT, source, tag, MPI_COMM_WORLD,&reqs[i+4]);}MPI_Waitall(8, reqs, stats);printf("rank= %d coords= %d %d neighbors(u,d,l,r)= %d %d %d %dinbuf(u,d,l,r)= %d %d %d %d\n",rank,coords[0],coords[1],nbrs[UP],nbrs[DOWN],nbrs[LEFT],inbuf[UP],inbuf[DOWN],inbuf[LEFT],inbuf[RIGHT]);}elseprintf("Must specify %d tasks. Terminating.\n",SIZE);MPI_Finalize(); }Лекции спецкурс "Пар. прогр. длявысокопроизводительных систем",24.11.201531Блочные параллельные алгоритмыматричного умноженияЛекции спецкурс "Пар.

прогр. длявысокопроизводительных систем",24.11.201532Матричное умножениеСуществует множество вариантов решения этой задачи намногопроцессорных системах.Алгоритм решения существенным образом зависит от того,производится или нет распределение матриц по процессорам,и какая топология процессоров при этом используется.Лекции спецкурс "Пар. прогр. длявысокопроизводительных систем",24.11.201533Распределение матриц по процессамКаждая из трех матриц (A,B и С) может быть распределенаодним из способов:копии матриц находятся в каждом процессе;распределена по столбцам на одномерную сетку;распределена по строкам на одномерную сетку;распределена на двумерную или трехмернуюпроцессную сетку.Могут использоваться и различные комбинации.

Всезависим от решаемой задачи.Лекции спецкурс "Пар. прогр. длявысокопроизводительных систем",24.11.201534Распределение матриц по процессамВозможное распределение матриц:- матрица A и C – 1D распределение (ленточное, построкам)- Матрица В – 1D (ленточное, по столбцам)Лекции спецкурс "Пар. прогр. длявысокопроизводительных систем",24.11.201535Умножение матрицУмножение матрицС= A * B, гдеA - n x l матрица and B l x m матрицаci,j (0 <= i < n, 0 <= j < m) вычисляются :Лекции спецкурс "Пар.

прогр. длявысокопроизводительных систем",24.11.201536Умножение матриц: C = A x BЛекции спецкурс "Пар. прогр. длявысокопроизводительных систем",24.11.201537Умножение матрицы на векторc=AxbЛекции спецкурс "Пар. прогр. длявысокопроизводительных систем",24.11.201538Блочное умножение матрицЛекции спецкурс "Пар. прогр. длявысокопроизводительных систем",24.11.201539Последовательный алгоритмfor (i = 0; i < n; i++)for (j = 0; j < n; j++) {c[i][j] = 0;for (k = 0; k < n; k++)c[i][j] = c[i][j] + a[i][k] * b[k][j];}n3 операций умножения и n3 операций сложенияСложность алгоритма O(n3).Легко распараллеливаетсяЛекции спецкурс "Пар. прогр. длявысокопроизводительных систем",24.11.201540Блочный алгоритмДелим матрицы на подматрицыПусть матрицы разделены на s^2 подматриц.Каждый блок содержит n/s x n/s элементов.Обозначим Ap,q - блок матрицы Afor (p = 0; p < s; p++)for (q = 0; q < s; q++) {Cp,q = 0;/* обнуление блоков */for (r = 0; r < m; r++)/*блочное умножение &*/Cp,q = Cp,q + Ap,r * Br,q; /*и сложение блоков*/}СтрокаCp,q = Cp,q + Ap,r * Br,q;означает умножение блоков Ap,r и Br,q, используя матричное умножениеЛекции спецкурс "Пар.

прогр. длявысокопроизводительных систем",24.11.201541Блочные алгоримы матричногоумноженияАлгоритм ФоксаАлгоритм КеннонаЛекции спецкурс "Пар. прогр. длявысокопроизводительных систем",24.11.201542Параллельный алгоритм/* distribute matrices by rows. */void Parallel_matrix_mult(MATRIX_T A, MATRIX_T B, MATRIX_T C, int n){for each column of B {Allgather(column);Compute dot product of my row of A with column;}/* can distribute matrices by blocks of rows.

Also B could be distributed by* columns*/Parallel Processing43Идея алгоритма КеннонаАлгоритм Кеннона определяет порядоксуммирования членов во внутреннем циклеC(i,j) = s-1 A(i, (i  j  k) mod s) * B( (i  j  k) mod s, j )k 0таким образом, чтобы в каждом процессе на каждомшаге алгоритма находился один из блоков матриц Aи B. Предусматривается первоначальноераспределение блоков матриц таким образом, чтобыминимизировать обмены блоками в процессевыполнения алгоритма.Лекции спецкурс "Пар. прогр. длявысокопроизводительных систем",24.11.201544Алгоритм Кеннона матричного умноженияfor all (i=0 to s-1)// начальное распределение блоков матрицы AЦиклический сдвиг влево строки i матрицы A на jтак, чтобы на место A(i,j) была записана подматрица A(i,(i+j) mod s) endfor for all (i=0 to s-1) // начальное распределение блоков матрицы BЦиклический сдвиг вверх столбца j матрицы B на jтак, чтобы на место B(i,j) была записана подматрица B((i+j) mod s,j)endfor for k=0 to s-1for all (i=0 to s-1, j=0 to s-1)C(i,j) = C(i,j) + A(i,j)*B(i,j)Лекции спецкурс "Пар.

прогр. длявысокопроизводительных систем",24.11.201545Алгоритм Кеннона матричного умноженияЦиклический сдвиг влево каждой строки матрицы A на 1так,чтобы на место A(i,j) была записана подматрица A(i,(j+1) mod s)Циклический сдвиг вверх каждого столбца матрицы B на 1так,чтобы на место B(i,j) была записана подматрица B((i+1) mod s,j)end forend forЛекции спецкурс "Пар.