latest.c

/*
   latest.c

   test program for the lalib library

   usage:
              latest rows cols size
   
       rows, cols: dimensions of processor grid
       size:       dimension N of global vector

*/

#include <stdio.h>
#include <math.h>
#include "lalib.h"


void main(int argc, char *argv[]) {
  int             rows;       /* dims of processor grid */
  int             cols;
  int             size;       /* size of vectors */
  int             blocklength; 
  int             me;         /* rank in WORLD */
  Lalib_grid *    grid;       /* 2d processor grid */

  void test_dot(Lalib_grid *grid, int blocklen, int size);
  
  MPI_Init(&argc, &argv);  /* enroll in MPI */
  MPI_Comm_rank(MPI_COMM_WORLD, &me);
  
  /* get arguments */
  if (argc != 5) {
    fprintf(stderr,"Usage: %s <rows> <cols> <vector size> <block size>\n",
	    argv[0]);
    MPI_Abort(MPI_COMM_WORLD, -1);
    exit(1);
  }
   
  rows        = atoi(argv[1]);
  cols        = atoi(argv[2]);
  size        = atoi(argv[3]);
  blocklength = atoi(argv[4]);
   
  /* initialize grid */
  grid = lalib_build_grid(MPI_COMM_WORLD, rows, cols);

  /* test scalar product */
  test_dot(grid, blocklength, size);

  MPI_Finalize();
}

void test_dot(Lalib_grid *grid, int blocklen, int size) {
  Lalib_vector *  a;
  Lalib_vector *  b;
  double          sum_s;
  double          sum_no;
  double          c;
  int             i;
  int             i_glob;
  
  /* build a new col vector a(i) = (i+1) / sqrt(size) */
  c = 1.0 / sqrt(size);
  a = lalib_build_vector(size, LALIB_COL, blocklen, grid);
  for (i=0; i<a->local_length; i++) {
    lalib_vector_global_index(a,  grid->myrow,  i, &i_glob);
    a->data[i] =  c * (i_glob + 1);
  }

  /* build a new col vector b(i) = 1 / (i+1) * sqrt(size) */
  b = lalib_build_vector(size, LALIB_COL, blocklen, grid);
  for (i=0; i<b->local_length; i++) {
    lalib_vector_global_index(a,  grid->myrow,  i, &i_glob);
    b->data[i] = c / (i_glob + 1);
  }

  /* compute scalar product */
  lalib_dot_stride(a, b, &sum_s);
  lalib_dot_nostride(a, b, &sum_no);

  /* print them */
  printf("  proc i_glob  i   a[i]   b[i]\n");
  printf("==============================\n");
  for (i=0; i<a->local_length; i++) {
    lalib_vector_global_index(a,  grid->myrow,  i, &i_glob);
    printf("(%2d,%2d) %2d    %2d %5.2lf  %5.2lf\n", 
	   grid->myrow, grid->mycol, i_glob, i, 
	   a->data[i], b->data[i]);
  }
  printf("\n scalar product a*b: %5.2lf (stride) %5.2lf (nostride)\n",
	 sum_s, sum_no);
  printf("====================================================\n");
}