/***************************************************************************

  Parallel Linear Algebra Package Release R2.0.2

  6 Feb 2000
  
  Copyright (c) 1997,1998,1999,2000 Robert van de Geijn and 
  The University of Texas at Austin.
  See the file README for details on the gnu license

  This program is free software; you can redistribute it and/or modify
  it under the terms of the GNU General Public License as published by
  the Free Software Foundation; either version 1, or (at your option)
  any later version.

  This program is distributed in the hope that it will be useful,
  but WITHOUT ANY WARRANTY; without even the implied warranty of
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  GNU General Public License for more details.

  You should have received a copy of the GNU General Public License
  along with this program; if not, write to the Free Software
  Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.

  Written under the direction of: 
 	Robert van de Geijn, Department of Computer Sciences,
	University of Texas at Austin  78712.    
	rvdg@cs.utexas.edu
 
  Many people contributed to the Parallel Linear Algebra Package (PLAPACK)
  including:	
  
  Gregory A. Baker, Center for Space Research,
  University of Texas at Austin  78722      baker@csr.utexas.edu

  Philip A. Alpatov, Department of Physics,
  University of Texas at Austin 78712.     philip@physics.utexas.edu

  James Overfelt, Department of Mathematics,
  University of Texas at Austin  78712.    overfelt@math.utexas.edu
  
  John Andrew Gunnels, Department of Computer Sciences,
  University of Texas at Austin 78712.     gunnels@cs.utexas.edu
  
  Greg Morrow, Department of Physics, 
  University of Texas at Austin 78712.     morrow@physics.utexas.edu

  Wesley Reiley, Department of Computer Sciences
  University of Texas at Austin 78712.     wesley@cs.utexas.edu

  Dr. Robert Van de Geijn, Department of Computer Sciences,
  University of Texas at Austin 78712.     rvdg@cs.utexas.edu

  This release (Release 2.0.2) was written primarily by
  Robert van de Geijn
 			
***********i****************************************************************/


#include "PLA.h"

int PLA_Gemv( int transa, 
	       PLA_Obj alpha, PLA_Obj A, PLA_Obj x, 
               PLA_Obj beta,  PLA_Obj y )
/*
  Purpose : Parallel matrix vector multiplication

  IN     transa      integer, PLA_NO_TRANSPOSE, PLA_TRANSPOSE, 
                              or PLA_CONJUGATE TRANSPOSE
  IN     alpha       multiscalar, scale factor for Ax
  IN     A           matrix
  IN     x           multivector 
  IN     beta        multiscalar, scale factor for y
  IN/OUT y           multivector, overwritten with 
                     y <- alpha * A  * x + beta * y  or 
                     y <- alpha * A' * x + beta * y 

  NOTE:  For details on how to implement parallel matrix-vector multiplication,
         see
	 
	 R. van de Geijn, Using PLAPACK, The MIT Press, 1997.
*/
{
  int 
    value = PLA_SUCCESS,
    owner_row, owner_col;
  PLA_Obj
    alpha_cpy = NULL, beta_cpy = NULL,
    x_dup = NULL, y_dup = NULL, zero = NULL;
  
  if ( PLA_ERROR_CHECKING )    /* Perform parameter and error checking */
    value = PLA_Gemv_enter( transa, alpha, A, x, beta, y );

  /* Create usual constants */
  PLA_Create_constants_conf_to( A, NULL, &zero, NULL );

  if ( value == PLA_SUCCESS ){
    /* If necessary, duplicate alpha and beta to all nodes */
    PLA_Obj_owner_row( alpha, &owner_row );
    PLA_Obj_owner_col( alpha, &owner_col );
    if ( owner_row != PLA_ALL_ROWS || owner_col != PLA_ALL_COLS ){
      PLA_Mscalar_create_conf_to( alpha, PLA_ALL_ROWS, PLA_ALL_COLS,
				   &alpha_cpy );
      PLA_Copy( alpha, alpha_cpy );
    }

    PLA_Obj_owner_row( beta, &owner_row );
    PLA_Obj_owner_col( beta, &owner_col );
    if ( owner_row != PLA_ALL_ROWS || owner_col != PLA_ALL_COLS ){
      PLA_Mscalar_create_conf_to( beta, PLA_ALL_ROWS, PLA_ALL_COLS,
				   &beta_cpy );
      PLA_Copy( beta, beta_cpy );
    }

    if ( transa == PLA_NO_TRANS ){
      /* Create duplicated projected multivectors to hold copy of x and
         local contributions to y */
      PLA_Pmvector_create_conf_to( A, PLA_PROJ_ONTO_ROW, PLA_ALL_ROWS, 1,
				    &x_dup );
      PLA_Pmvector_create_conf_to( A, PLA_PROJ_ONTO_COL, PLA_ALL_COLS, 1,
				    &y_dup );
      
      /* Duplicate x */
      PLA_Copy( x, x_dup );

      /* Perform local matrix-vector multiply */
      PLA_Local_gemv( transa, 
		       ( alpha_cpy == NULL ? alpha: alpha_cpy ), 
		       A, x_dup, zero, y_dup );
	
      /* Add local contributions to vector y */
      PLA_Reduce_x( PLA_SHAPE_GENERAL, y_dup, 
		     ( beta_cpy == NULL ? beta : beta_cpy ), y );
    }
    else{
      /* Create duplicated projected multivectors to hold copy of x and
         local contributions to y */
      PLA_Pmvector_create_conf_to( A, PLA_PROJ_ONTO_COL, PLA_ALL_COLS, 1,
				    &x_dup );
      PLA_Pmvector_create_conf_to( A, PLA_PROJ_ONTO_ROW, PLA_ALL_ROWS, 1,
				    &y_dup );
      
      /* Duplicate x */
      PLA_Copy( x, x_dup );

      /* Perform local matrix-vector multiply */
      PLA_Local_gemv( transa, 
		       ( alpha_cpy == NULL ? alpha: alpha_cpy ), 
		       A, x_dup, zero, y_dup );
	
      /* Add local contributions to vector y */
      PLA_Reduce_x( PLA_SHAPE_GENERAL, y_dup, 
		     ( beta_cpy == NULL ? beta : beta_cpy ), y );
    }
  }

  /* free temporary objects and views */
  PLA_Obj_free( &alpha_cpy );
  PLA_Obj_free( &beta_cpy );
  PLA_Obj_free( &x_dup );
  PLA_Obj_free( &y_dup );
  PLA_Obj_free( &zero );

  if ( PLA_ERROR_CHECKING )    /* Perform exit error checking */
    value = PLA_Gemv_exit( transa, alpha, A, x, beta, y );

  return value;
}