y_eq_diagAx.F
#include "iadefines.h"
c**********************************************************************
#include "author.inc"
c* $Id: y_eq_diagAx.F,v 1.6 1996/04/24 19:15:12 turner Exp $
c*
c* Compute y = diag(A)*x.
c*
c* <PARAMETER LIST>
c*
c* Input:
c* x - vector to be multiplied
c* a - matrix
c* ia - integer vector containing info about how "a" is stored
c* NOTE: see description of ia below
c* ja - column index array for a
c*
c* Output:
c* y - destination vector
c* status - return status
c* -1 ==> invalid argument(s)
c* 0 ==> success
c*
c* <SUBROUTINES REQUIRED>
c*
c* JT_y_eq_diagAx_Full
c* JT_z_eq_yx
c*
#include "iadesc.inc"
c*
#include "copyright.inc"
c**********************************************************************
subroutine JT_y_eq_diagAx (x, a, ia, ja, y, status)
implicit none
c
c ... Input:
integer ia(_JT_no_of_storage_parameters_), ja(*)
real x(*), a(*)
c
c ... Output:
integer status
real y(*)
c
if (ia(_JT_storage_) .eq. 1) then
call JT_y_eq_diagAx_Full (ia(_JT_idim_), ia(_JT_nrows_), x, a, y, status)
elseif (ia(_JT_storage_).eq.8 .or.
& ia(_JT_storage_).eq.10 .or.
& ia(_JT_storage_).eq.11 .or.
& ia(_JT_storage_).eq.12) then
call JT_z_eq_yx (ia(_JT_nrows_), a, x, y, status)
else
status = -1
endif
c
return
end
c**********************************************************************
#include "author.inc"
c*
c* Compute y = inv(diag(A))*x.
c*
c* <PARAMETER LIST>
c*
c* Input:
c* x - vector to be multiplied
c* a - matrix
c* ia - integer vector containing info about how "a" is stored
c* NOTE: see description of ia below
c* ja - column index array for a
c*
c* Output:
c* y - destination vector
c* status - return status
c* -1 ==> invalid argument(s)
c* 0 ==> success
c*
c* <SUBROUTINES REQUIRED>
c*
c* JT_y_eq_invdiagAx_Full
c* JT_z_eq_y_div_x
c*
#include "iadesc.inc"
c*
#include "copyright.inc"
c**********************************************************************
subroutine JT_y_eq_invdiagAx (x, a, ia, ja, y, status)
implicit none
c
c ... Input:
integer ia(_JT_no_of_storage_parameters_), ja(*)
real x(*), a(*)
c
c ... Output:
integer status
real y(*)
c
if (ia(_JT_storage_) .eq. _JT_storage_full_) then
call JT_y_eq_invdiagAx_Full (ia(_JT_idim_), ia(_JT_nrows_), x, a, y, status)
elseif (ia(_JT_storage_).eq._JT_storage_ELL_ .or.
& ia(_JT_storage_).eq._JT_storage_COO_ .or.
& ia(_JT_storage_).eq._JT_storage_RSS_ .or.
& ia(_JT_storage_).eq._JT_storage_CSS_) then
call JT_z_eq_y_div_x (ia(_JT_nrows_), a, x, y, status)
else
status = -1
endif
c
return
end