Assignment-04/Question-05/matvec.c
// Master process, divides up the matrix into rows for computation by slaves
Master process, divides up the matrix into rows for computation by slaves
Master process divides up the matrix into rows for computation by slaves
void master(int n_proc, int n_rows, int n_cols) {
    double *x = (double *) malloc(n_cols*sizeof(double));
    double *y = (double *) malloc(n_rows*sizeof(double));
    double *z = (double *) malloc(n_rows*sizeof(double));
    double *A = (double *) malloc(n_rows*n_cols*sizeof(double));
    double *sub = (double *) malloc(n_cols*sizeof(double));
    double ans = 0.0;
    int n_sent = 0, row = 0, sender = 0;
    MPI_Status status;
    for (int i = 0; i < n_cols; ++i) {
        x[i] = i+1;
    } 
    for (int i = 0; i < n_rows; ++i) {
        for (int j = 0; j < n_cols; ++j) {
            A[(i*n_cols)+j] = i+1;
        }
    }
    MPI_Bcast(x, n_cols, MPI_DOUBLE, MASTER, MPI_COMM_WORLD);
    for (int i = 1; i < min(n_rows+1, n_proc); ++i) {
    	memcpy((void *) sub, (void *) &A[n_sent*n_cols], n_cols * sizeof(double));
        MPI_Send(sub, n_cols, MPI_DOUBLE, i, n_sent, MPI_COMM_WORLD);
        n_sent++;
    }
    for (int i = 0; i < n_rows; ++i) {
        MPI_Recv(&ans, 1, MPI_DOUBLE, MPI_ANY_SOURCE, MPI_ANY_TAG,
                 MPI_COMM_WORLD, &status);
        sender = status.MPI_SOURCE;
        row = status.MPI_TAG;
        y[row] = ans;
        if (n_sent < n_rows) {
        	memcpy((void *) sub, (void *) &A[n_sent*n_cols], n_cols * sizeof(double));
            MPI_Send(sub, n_cols, MPI_DOUBLE, sender, n_sent, MPI_COMM_WORLD);
            n_sent++;
        }
        else {
            MPI_Send(y, n_rows, MPI_DOUBLE, sender, n_rows+1, MPI_COMM_WORLD);
        }
    }
    printf("Finished!\n");
    cblas_dgemv(CblasRowMajor, CblasNoTrans, n_rows, n_cols, 1.0, A, n_cols,
	      		x, 1, 0.0, z, 1);
    for (int i = 0; i < n_rows; ++i) {
    	y[i] -= z[i];
    }
    printf("Check: ||y-Ax|| = %f\n", cblas_dnrm2(n_rows, y, 1));
    free(x);
    free(y);
    free(z);
    free(A);
    free(sub);
}
void master ( int n_proc , int n_rows , int n_cols ) { double * x = ( double * ) malloc ( n_cols * sizeof ( double ) ) ; double * y = ( double * ) malloc ( n_rows * sizeof ( double ) ) ; double * z = ( double * ) malloc ( n_rows * sizeof ( double ) ) ; double * A = ( double * ) malloc ( n_rows * n_cols * sizeof ( double ) ) ; double * sub = ( double * ) malloc ( n_cols * sizeof ( double ) ) ; double ans = 0.0 ; int n_sent = 0 , row = 0 , sender = 0 ; MPI_Status status ; for ( int i = 0 ; i < n_cols ; ++ i ) { x [ i ] = i + 1 ; } for ( int i = 0 ; i < n_rows ; ++ i ) { for ( int j = 0 ; j < n_cols ; ++ j ) { A [ ( i * n_cols ) + j ] = i + 1 ; } } MPI_Bcast ( x , n_cols , MPI_DOUBLE , MASTER , MPI_COMM_WORLD ) ; for ( int i = 1 ; i < min ( n_rows + 1 , n_proc ) ; ++ i ) { memcpy ( ( void * ) sub , ( void * ) & A [ n_sent * n_cols ] , n_cols * sizeof ( double ) ) ; MPI_Send ( sub , n_cols , MPI_DOUBLE , i , n_sent , MPI_COMM_WORLD ) ; n_sent ++ ; } for ( int i = 0 ; i < n_rows ; ++ i ) { MPI_Recv ( & ans , 1 , MPI_DOUBLE , MPI_ANY_SOURCE , MPI_ANY_TAG , MPI_COMM_WORLD , & status ) ; sender = status . MPI_SOURCE ; row = status . MPI_TAG ; y [ row ] = ans ; if ( n_sent < n_rows ) { memcpy ( ( void * ) sub , ( void * ) & A [ n_sent * n_cols ] , n_cols * sizeof ( double ) ) ; MPI_Send ( sub , n_cols , MPI_DOUBLE , sender , n_sent , MPI_COMM_WORLD ) ; n_sent ++ ; } else { MPI_Send ( y , n_rows , MPI_DOUBLE , sender , n_rows + 1 , MPI_COMM_WORLD ) ; } } printf ( " \n " ) ; cblas_dgemv ( CblasRowMajor , CblasNoTrans , n_rows , n_cols , 1.0 , A , n_cols , x , 1 , 0.0 , z , 1 ) ; for ( int i = 0 ; i < n_rows ; ++ i ) { y [ i ] -= z [ i ] ; } printf ( " \n " , cblas_dnrm2 ( n_rows , y , 1 ) ) ; free ( x ) ; free ( y ) ; free ( z ) ; free ( A ) ; free ( sub ) ; }
Master process, divides up the matrix into rows for computation by slaves
Master process divides up the matrix into rows for computation by slaves
// 0. Allocate space for arrays as needed to compute y = A*x // 1. Define data A (n_rows x n_cols matrix) and x (n_cols x 1 vector) // 2. Broadcast vector x to all processes // 3. Send first rows to other processes // Now loop over remaining rows, receiving dot products computed // Receive a computed inner product from one of the processes // Either send another row back to sender or send final y=Ax // along with a tag to signal completion. // Inform sender that there is no more work // Calculate the norm to verify ||y-Ax|| is correct