y_eq_Ax.F


#include "iadefines.h"
c**********************************************************************
#include "author.inc"
c*    $Id: y_eq_Ax.F,v 1.17 1996/04/24 19:15:11 turner Exp $
c*
c*    Computes y=Ax.
c*
c*    <PARAMETER LIST>
c*
c*     Input:
c*      a - matrix
c*      ia - integer vector containing info about how "a" is stored
c*        NOTE: see description of ia below
c*      ja - index array for matrix
c*      x - x-vector
c*
c*     Output:
c*      y - y-vector
c*      status - return status
c*        -3  ==>  internal error
c*        -2  ==>  memory allocation failure
c*        -1  ==>  invalid argument(s)
c*         0  ==>  success
c*
#include "iadesc.inc"
c*
c*    <SUBROUTINES REQUIRED>
c*
c*     JT_y_eq_Ax_ELL   for A stored in ELLPACK-ITPACK format
c*     JT_y_eq_Ax_Full  for A stored in full conventional format
c*     JT_y_eq_Ax_COO   for A stored in coordinate format
c*     JT_y_eq_Ax_RSS   for A stored in row sequential storage
c*     JT_y_eq_Ax_CSS   for A stored in column sequential storage
c*
#include "copyright.inc"
c**********************************************************************
      subroutine JT_y_eq_Ax (a, ia, ja, x, y, status)
      implicit none
c
c ... Input:
      integer ia(_JT_no_of_storage_parameters_), ja(*)
      real a(*)
#ifdef strict_f77
      real x(*)
#else
      real x(ia(_JT_nrows_))
#endif
c
c ... Output:
      integer status
#ifdef strict_f77
      real y(*)
#else
      real y(ia(_JT_nrows_))
#endif
c
      if (ia(_JT_storage_) .eq. _JT_storage_full_) then
       call JT_y_eq_Ax_Full (ia(_JT_idim_), ia(_JT_nrows_), a, x, y, status)
      elseif (ia(_JT_storage_) .eq. _JT_storage_ELL_) then
       call JT_y_eq_Ax_ELL (ia(_JT_idim_), ia(_JT_nrows_), ia(_JT_maxnz_), a, ja, x,
     &      y, status)
       if (status .eq. -1) status = -3
      elseif (ia(_JT_storage_) .eq. _JT_storage_COO_) then
       call JT_y_eq_Ax_COO (ia(_JT_nelem_), ia(_JT_nrows_), a, ja, x,
     &      y, status)
       if (status .eq. -1) status = -3
      elseif (ia(_JT_storage_) .eq. _JT_storage_RSS_) then
       call JT_y_eq_Ax_RSS (ia(_JT_nelem_), ia(_JT_nrows_), a, ja, x,
     &      y, status)
       if (status .eq. -1) status = -3
      elseif (ia(_JT_storage_) .eq. _JT_storage_CSS_) then
       call JT_y_eq_Ax_CSS (ia(_JT_nelem_), ia(_JT_nrows_), a, ja, x,
     &      y, status)
       if (status .eq. -1) status = -3
      else
       status = -1   ! invalid value for ia(_JT_storage_)
      endif
c
      return
      end