Инженерная методика адаптации приложения к гибридному кластеру с ускорителями на ПЛИС
Часть 4. Искусственное расширение теневой грани в локальной порции.
Обратим внимание на то, что, при разбиении локальной порции на мелкие блоки, каждый из которых умещается в сопроцессоре, блоки эти придется нарезать с перекрытиями (теневыми гранями), подобно тому, как нарезается на локальные порции распределенный массив при распараллеливании приложения. Очевидно, что размер перекрытия при нарезке мелких блоков такой же, как и при нарезке локальных порций в параллельной реализации (в данном случае - одна строка). А что будет, если нарезать мелкие блоки с перекрытием в две строки? Каждый такой мелкий блок можно будет, единожды скопировав в сопроцессор, обработать два раза. Правда, при обратном копировании обработанных данных в память, придется "выбросить" области перекрытия удвоенного размера, но, по сравнению с двойной экономией на копировании данных, это не страшно. Легко показать, что при увеличении размера перекрытия в N раз по сравнению с требованиями исходного численного метода, единожды скопированный в сопроцессор мелкий блок может быть обработан внутри сопроцессора N раз. Значит, следующий шаг преобразования нашей программы - это включение в нее нарезки локальной порции на мелкие блоки с произвольной шириной перекрытия? Нет, к этому преобразованию переходить еще рано. Для того, чтобы многократная обработка мелких блоков данных была корректной, необходимо, чтобы была корректной многократная (без обменов краями с соседями) обработка всей локальной порции. Для этого, в свою очередь, потребуется нарезать уже сами локальные порции с расширенными в N раз перекрытиями, а обмены такими расширенными краями сделать в N раз более редкими. С этого преобразования программы и начнем. Текст преобразованной таким образом программы приводится ниже.
Текст функции main():
#include <mpi.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "vsda.h"
#include "dimension.h"
// Change here the number of steps, the cell geometry, etc
#define STEPITER 1000
#define delx 0.5
#define dely 0.25
#define NITER 10002
#define STW 3
// end change here.
/***/
extern void itstep( int mx, int my, void *f, void *newf, void *r,
double rdx2, double rdy2, double beta, int niters );
int main( int argc, char **argv )
{
int i, j, k, n, mx1, mx2, my1, my_number, n_of_nodes,
totalmx, partmx, leftmx, mx, my, niter, stw;
FILE *fp;
double howlong;
double rdx2, rdy2, beta;
VSDA vsda;
long *distr_out;
MPI_Request rq[5];
/***/
MPI_Init( &argc, &argv );
MPI_Comm_rank( MPI_COMM_WORLD, &my_number );
MPI_Comm_size( MPI_COMM_WORLD, &n_of_nodes );
if ( argc != 3 )
{
if ( !my_number )
{
fprintf( stderr, "Usage: %s <nrows> <ncolumns>\n", argv[0] );
fflush( stderr );
}
return (-1);
}
totalmx = mx = (int)atol( argv[1] );
my = (int)atol( argv[2] );
if ( my < 1 )
{
if ( !my_number )
{
fprintf( stderr, "Number of columns (%d) should be positive\n", my );
fflush( stderr );
}
return (-1);
}
distr_out = (long*)malloc( n_of_nodes*sizeof( *distr_out ) );
if ( !distr_out )
{
fprintf( stderr, "No memory\n" );
fflush( stderr );
return( -1 );
}
totalmx += 2;
my += 2;
stw = STW;
niter = NITER;
/* Make STW odd, NITER divisible by STW: */
if ( !(stw%2) ) stw--;
niter -= (niter%stw);
/* Here we know the array width, so make the arrays themselves: */
{
DIM2( double, f, my );
DIM2( double, newf, my );
DIM2( double, r, my );
void *tmpf;
/* Make the distributed array, discover the local portion length: */
VSDA_Create( MPI_COMM_WORLD, VSDA_BLOCKED, 0,
(long)totalmx, (long)sizeof(*f), (long)stw, NULL, distr_out, &vsda );
mx = distr_out[my_number];
if ( n_of_nodes > 1 )
{
if ( my_number > 0 ) mx += stw;
if ( my_number < (n_of_nodes-1) ) mx += stw;
}
f = (typeof(f))malloc( mx*sizeof(*f) );
newf = (typeof(newf))malloc( mx*sizeof(*newf) );
r = (typeof(r))malloc( mx*sizeof(*r) );
if ( (!f) || (!newf) || (!r) )
{
fprintf( stderr, "Cannot allocate, exiting\n" );
exit( -1 );
}
rdx2 = 1./delx/delx;
rdy2 = 1./dely/dely;
beta = 1.0/(2.0*(rdx2+rdy2));
if ( !my_number )
{
printf( "Solving heat conduction task on %d by %d grid by %d processors\n",
totalmx-2, my-2, n_of_nodes );
fflush( stdout );
}
for ( i = 0; i < mx; i++ )
{
for ( j = 0; j < my; j++ )
{
if ( ((i==0)&&(my_number==0)) || (j==0) ||
((i==(mx-1))&&(my_number==(n_of_nodes-1))) || (j==(my-1)) )
newf[i][j] = f[i][j] = 1.0;
else
newf[i][j] = f[i][j] = 0.0;
r[i][j] = 0.0;
}
}
/* Iteration loop: */
howlong = MPI_Wtime();
for ( n = 0; n < niter/stw; n++ )
{
if ( !my_number )
{
if ( !(n%STEPITER) )
printf( "Iteration %d\n", n );
}
/* Step of calculation starts here: */
itstep( mx, my, f, newf, r, rdx2, rdy2, beta, stw );
/* Do all the transfers: */
VSDA_Update_begin( &vsda, rq, newf );
VSDA_Update_end( rq );
/* swap the halves: */
tmpf = newf;
newf = f;
f = (typeof(f))tmpf;
}
if ( !my_number )
{
printf( "Elapsed time: %f sec\n", (float)(MPI_Wtime()-howlong) );
fp = fopen( "output_par_2.dat", "w" );
fclose( fp );
}
for ( j = 0; j < n_of_nodes; j++ )
{
MPI_Barrier( MPI_COMM_WORLD );
if ( j == my_number )
{
int istart, isize;
istart = 1;
isize = mx-1;
if ( n_of_nodes > 1 )
{
if ( my_number > 0 ) istart = stw;
if ( my_number < (n_of_nodes-1) ) isize = mx-stw;
}
fp = fopen( "output_par_2.dat", "a" );
for ( i = istart; i < isize; i++ )
{
fwrite( &(f[i][1]), my-2, sizeof(f[0][0]), fp );
}
fclose( fp );
}
}
}
MPI_Finalize();
return 0;
}
Обратим внимание на то, что не только усложнять, но и вообще как-либо изменять ту часть логики программы, которая организует обмены данными между MPI - процессами, нам не пришлось, что свидетельствует о действительной пригодности библиотеки VSDA для упрощения разработки задач этого класса.
Текст функции itstep():
#include "dimension.h"
void itstep( int mx, int my, void *pf, void *pnewf,
void *pr, double rdx2, double rdy2, double beta, int niter )
{
DIM2( double, f, my ) = (typeof(f))pf;
DIM2( double, newf, my ) = (typeof(newf))pnewf;
DIM2( double, r, my ) = (typeof(r))pr;
int i, j, k, mx1, my1;
void *tmpf;
/***/
mx1 = mx - 1;
my1 = my - 1;
for ( k = 0; k < niter; k++ )
{
for ( i = 1; i < mx1; i++ )
{
for ( j = 1; j < my1; j++ )
{
newf[i][j] = ((f[i-1][j]+f[i+1][j])*rdx2+(f[i][j-1]+f[i][j+1])*rdy2-r[i][j])*beta;
}
}
tmpf = newf;
newf = f;
f = (typeof(f))tmpf;
}
}
Приведем необходимые пояснения. Коэффициент расширения теневых граней задан константой STW (STencil Width, ширина шаблона) в тексте функции main(). Чтобы избежать еще большего усложнения программы, нам пришлось наложить довольно противоестественные ограничения на исходные данные. Мы требуем, чтобы число итераций было кратно STW, и чтобы само STW было нечетным. Последнее требование связано с тем, что и в функции main(), и в вызываемой из нее функции itstep(), независимо друг от друга меняются местами указатели f и newf. Для корректной работы программы необходимо, чтобы они менялись местами согласованно, то есть чтобы после вызова itstep() результаты счета были в том массиве, на которую в вызывающей программе указывает newf. А для этого необходимо и достаточно, чтобы число замен местами f и newf, которое равно STW, было нечетным.
◄ Часть 3 Часть 5 ►