Yet Another eXchange Tool 0.11.3
Loading...
Searching...
No Matches
Round Robin example on how to use yaxt

This example(rrobin.c) is the first step to understand working with yaxt. First of all include the basic header files:

#include <mpi.h>
#include <yaxt.h>

We need to initialize MPI and yaxt as follows:

int main(void) {
//init mpi
xt_mpi_call(MPI_Init(NULL, NULL), MPI_COMM_WORLD);
@ MPI_COMM_WORLD
Definition core.h:73
void xt_initialize(MPI_Comm default_comm)
Definition xt_init.c:70
#define xt_mpi_call(call, comm)
Definition xt_mpi.h:68

Find out number of processes and the local rank like in common MPI programs.

int rank, size;
MPI_Comm_rank (MPI_COMM_WORLD, &rank);
MPI_Comm_size (MPI_COMM_WORLD, &size);

In this exmaple we are going to create a source array of length 5 for each process.
Each process will have an array the is filled with unique values. It is recommendable to use an array length < 10 and a number of processes < 10 to get a readable output.
In this case each process gets values of the form "abc", where

  • a = 1
  • b = rank of the local process
  • c = index of the element within the local array

For example: 132 is the source value of on the process with rank 3 and is on the position 3 in the local array.
The goal of this little program is to rotate the source arrays in a round robin fashion. First proces owns first five elements indexed by 0, 1, 2, 3, 4 and is going to get the next five elements from the second process indexed by 5, 6, 7, 8, 9, without knowing who owns those. In fact first rank gets value array from second rank, the second one from the third and so on. The last rank gets an array from the first. We fill the source arrays with the previously defined values, and fill the destination arrays with -1 to see the change.

{
// init source and destination data array, local data 5 elements length
enum { len = 5 };

There are many ways to define, which elements are locally available (source) and which are required (destination). We could define them with an array of indices using an index vector (xt_idxvec.h), or we could define a block of elements we want to have using index stripes (xt_idxstripes.h). Using stripes we have to name the local start index, how many elements we want to have, an the stride between the elements. Here we need for the source an index stripe containing 5 elements with a stride of 1, beginnig at 0 for rank 0, at 1*len for rank 1 etc.

//print source
printf("SOURCEvalue: %d, element_index: %d, rank: %d \n", src_array[i],
i, rank);
}
// output barrier
MPI_Barrier (MPI_COMM_WORLD);
struct Xt_stripe src_stripes = {.start = (Xt_int)(rank*len), .stride = 1, .nstrides = len};
// source index list by stripes
Xt_idxlist src_idxlist = xt_idxstripes_new(&src_stripes, 1);
struct Xt_stripe dst_stripes = {.start = (Xt_int)(((rank+1)*len)%(size*len)),
.stride = 1, .nstrides = len};
// destination index list by stripes
Xt_idxlist dst_idxlist = xt_idxstripes_new(&dst_stripes, 1);
Xt_int start
Definition xt_stripe.h:55
XT_INT Xt_int
Definition xt_core.h:72
Xt_idxlist xt_idxstripes_new(struct Xt_stripe const *stripes, int num_stripes)

Now, we need the mapping of source and destination data between the processes and a redistribution object for the actual data exchange. There multiple strategies for doing the mapping, in this example all2all is used. An alternative would be dist_dir (xt_xmap_dist_dir.h).

// xmap
Xt_xmap xmap
= xt_xmap_all2all_new(src_idxlist, dst_idxlist, MPI_COMM_WORLD);
// redist_p2p
Xt_redist redist = xt_redist_p2p_new(xmap, MPI_INTEGER);
Xt_redist xt_redist_p2p_new(Xt_xmap xmap, MPI_Datatype datatype)
Xt_xmap xt_xmap_all2all_new(Xt_idxlist src_idxlist, Xt_idxlist dst_idxlist, MPI_Comm comm)

To do the main step, we need pointers of source and destination arrays. Here it is "overdressed", but shows the main charachter if you have higher number of data arrays.

// array pointer, especially necessary for data array number > 1
int* src_array_p = &src_array[0];
int* dst_array_p = &dst_array[0];
//Exchange
xt_redist_s_exchange1(redist, src_array_p, dst_array_p);
void xt_redist_s_exchange1(Xt_redist redist, const void *src_data, void *dst_data)
Definition xt_redist.c:92

To see the result:

for (int p = 0; p < len; p++)
printf("DESTvalue: %d, element_index: %d, rank: %d \n",

Once the created yaxt objects are not needed anymore they need to be deleted.

xt_idxlist_delete(dst_idxlist);
xt_idxlist_delete(src_idxlist);
}
void xt_idxlist_delete(Xt_idxlist idxlist)
Definition xt_idxlist.c:75
void xt_redist_delete(Xt_redist redist)
Definition xt_redist.c:74
void xt_xmap_delete(Xt_xmap xmap)
Definition xt_xmap.c:86

Common MPI finalisation

MPI_Finalize();
return 0;
}