LAPACK: SRC/zla_geamv.f Source File

Go to the documentation of this file.
00001       SUBROUTINE ZLA_GEAMV ( TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA,
00002      $                       Y, INCY )
00003 *
00004 *     -- LAPACK routine (version 3.2.2)                                 --
00005 *     -- Contributed by James Demmel, Deaglan Halligan, Yozo Hida and --
00006 *     -- Jason Riedy of Univ. of California Berkeley.                 --
00007 *     -- June 2010                                                    --
00008 *
00009 *     -- LAPACK is a software package provided by Univ. of Tennessee, --
00010 *     -- Univ. of California Berkeley and NAG Ltd.                    --
00011 *
00012       IMPLICIT NONE
00013 *     ..
00014 *     .. Scalar Arguments ..
00015       DOUBLE PRECISION   ALPHA, BETA
00016       INTEGER            INCX, INCY, LDA, M, N
00017       INTEGER            TRANS
00018 *     ..
00019 *     .. Array Arguments ..
00020       COMPLEX*16         A( LDA, * ), X( * )
00021       DOUBLE PRECISION   Y( * )
00022 *     ..
00023 *
00024 *  Purpose
00025 *  =======
00026 *
00027 *  ZLA_GEAMV  performs one of the matrix-vector operations
00028 *
00029 *          y := alpha*abs(A)*abs(x) + beta*abs(y),
00030 *     or   y := alpha*abs(A)'*abs(x) + beta*abs(y),
00031 *
00032 *  where alpha and beta are scalars, x and y are vectors and A is an
00033 *  m by n matrix.
00034 *
00035 *  This function is primarily used in calculating error bounds.
00036 *  To protect against underflow during evaluation, components in
00037 *  the resulting vector are perturbed away from zero by (N+1)
00038 *  times the underflow threshold.  To prevent unnecessarily large
00039 *  errors for block-structure embedded in general matrices,
00040 *  "symbolically" zero components are not perturbed.  A zero
00041 *  entry is considered "symbolic" if all multiplications involved
00042 *  in computing that entry have at least one zero multiplicand.
00043 *
00044 *  Arguments
00045 *  ==========
00046 *
00047 *  TRANS   (input) INTEGER
00048 *           On entry, TRANS specifies the operation to be performed as
00049 *           follows:
00050 *
00051 *             BLAS_NO_TRANS      y := alpha*abs(A)*abs(x) + beta*abs(y)
00052 *             BLAS_TRANS         y := alpha*abs(A')*abs(x) + beta*abs(y)
00053 *             BLAS_CONJ_TRANS    y := alpha*abs(A')*abs(x) + beta*abs(y)
00054 *
00055 *           Unchanged on exit.
00056 *
00057 *  M       (input) INTEGER
00058 *           On entry, M specifies the number of rows of the matrix A.
00059 *           M must be at least zero.
00060 *           Unchanged on exit.
00061 *
00062 *  N       (input) INTEGER
00063 *           On entry, N specifies the number of columns of the matrix A.
00064 *           N must be at least zero.
00065 *           Unchanged on exit.
00066 *
00067 *  ALPHA  - DOUBLE PRECISION
00068 *           On entry, ALPHA specifies the scalar alpha.
00069 *           Unchanged on exit.
00070 *
00071 *  A      - COMPLEX*16         array of DIMENSION ( LDA, n )
00072 *           Before entry, the leading m by n part of the array A must
00073 *           contain the matrix of coefficients.
00074 *           Unchanged on exit.
00075 *
00076 *  LDA     (input) INTEGER
00077 *           On entry, LDA specifies the first dimension of A as declared
00078 *           in the calling (sub) program. LDA must be at least
00079 *           max( 1, m ).
00080 *           Unchanged on exit.
00081 *
00082 *  X      - COMPLEX*16         array of DIMENSION at least
00083 *           ( 1 + ( n - 1 )*abs( INCX ) ) when TRANS = 'N' or 'n'
00084 *           and at least
00085 *           ( 1 + ( m - 1 )*abs( INCX ) ) otherwise.
00086 *           Before entry, the incremented array X must contain the
00087 *           vector x.
00088 *           Unchanged on exit.
00089 *
00090 *  INCX    (input) INTEGER
00091 *           On entry, INCX specifies the increment for the elements of
00092 *           X. INCX must not be zero.
00093 *           Unchanged on exit.
00094 *
00095 *  BETA   - DOUBLE PRECISION
00096 *           On entry, BETA specifies the scalar beta. When BETA is
00097 *           supplied as zero then Y need not be set on input.
00098 *           Unchanged on exit.
00099 *
00100 *  Y       (input/output) DOUBLE PRECISION  array, dimension
00101 *           ( 1 + ( m - 1 )*abs( INCY ) ) when TRANS = 'N' or 'n'
00102 *           and at least
00103 *           ( 1 + ( n - 1 )*abs( INCY ) ) otherwise.
00104 *           Before entry with BETA non-zero, the incremented array Y
00105 *           must contain the vector y. On exit, Y is overwritten by the
00106 *           updated vector y.
00107 *
00108 *  INCY    (input) INTEGER
00109 *           On entry, INCY specifies the increment for the elements of
00110 *           Y. INCY must not be zero.
00111 *           Unchanged on exit.
00112 *
00113 *
00114 *  Level 2 Blas routine.
00115 *
00116 *  =====================================================================
00117 *
00118 *     .. Parameters ..
00119       COMPLEX*16         ONE, ZERO
00120       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
00121 *     ..
00122 *     .. Local Scalars ..
00123       LOGICAL            SYMB_ZERO
00124       DOUBLE PRECISION   TEMP, SAFE1
00125       INTEGER            I, INFO, IY, J, JX, KX, KY, LENX, LENY
00126       COMPLEX*16         CDUM
00127 *     ..
00128 *     .. External Subroutines ..
00129       EXTERNAL           XERBLA, DLAMCH
00130       DOUBLE PRECISION   DLAMCH
00131 *     ..
00132 *     .. External Functions ..
00133       EXTERNAL           ILATRANS
00134       INTEGER            ILATRANS
00135 *     ..
00136 *     .. Intrinsic Functions ..
00137       INTRINSIC          MAX, ABS, REAL, DIMAG, SIGN
00138 *     ..
00139 *     .. Statement Functions ..
00140       DOUBLE PRECISION   CABS1
00141 *     ..
00142 *     .. Statement Function Definitions ..
00143       CABS1( CDUM ) = ABS( DBLE( CDUM ) ) + ABS( DIMAG( CDUM ) )
00144 *     ..
00145 *     .. Executable Statements ..
00146 *
00147 *     Test the input parameters.
00148 *
00149       INFO = 0
00150       IF     ( .NOT.( ( TRANS.EQ.ILATRANS( 'N' ) )
00151      $           .OR. ( TRANS.EQ.ILATRANS( 'T' ) )
00152      $           .OR. ( TRANS.EQ.ILATRANS( 'C' ) ) ) ) THEN
00153          INFO = 1
00154       ELSE IF( M.LT.0 )THEN
00155          INFO = 2
00156       ELSE IF( N.LT.0 )THEN
00157          INFO = 3
00158       ELSE IF( LDA.LT.MAX( 1, M ) )THEN
00159          INFO = 6
00160       ELSE IF( INCX.EQ.0 )THEN
00161          INFO = 8
00162       ELSE IF( INCY.EQ.0 )THEN
00163          INFO = 11
00164       END IF
00165       IF( INFO.NE.0 )THEN
00166          CALL XERBLA( 'ZLA_GEAMV ', INFO )
00167          RETURN
00168       END IF
00169 *
00170 *     Quick return if possible.
00171 *
00172       IF( ( M.EQ.0 ).OR.( N.EQ.0 ).OR.
00173      $    ( ( ALPHA.EQ.ZERO ).AND.( BETA.EQ.ONE ) ) )
00174      $   RETURN
00175 *
00176 *     Set  LENX  and  LENY, the lengths of the vectors x and y, and set
00177 *     up the start points in  X  and  Y.
00178 *
00179       IF( TRANS.EQ.ILATRANS( 'N' ) )THEN
00180          LENX = N
00181          LENY = M
00182       ELSE
00183          LENX = M
00184          LENY = N
00185       END IF
00186       IF( INCX.GT.0 )THEN
00187          KX = 1
00188       ELSE
00189          KX = 1 - ( LENX - 1 )*INCX
00190       END IF
00191       IF( INCY.GT.0 )THEN
00192          KY = 1
00193       ELSE
00194          KY = 1 - ( LENY - 1 )*INCY
00195       END IF
00196 *
00197 *     Set SAFE1 essentially to be the underflow threshold times the
00198 *     number of additions in each row.
00199 *
00200       SAFE1 = DLAMCH( 'Safe minimum' )
00201       SAFE1 = (N+1)*SAFE1
00202 *
00203 *     Form  y := alpha*abs(A)*abs(x) + beta*abs(y).
00204 *
00205 *     The O(M*N) SYMB_ZERO tests could be replaced by O(N) queries to
00206 *     the inexact flag.  Still doesn't help change the iteration order
00207 *     to per-column.
00208 *
00209       IY = KY
00210       IF ( INCX.EQ.1 ) THEN
00211          IF( TRANS.EQ.ILATRANS( 'N' ) )THEN
00212             DO I = 1, LENY
00213                IF ( BETA .EQ. 0.0D+0 ) THEN
00214                   SYMB_ZERO = .TRUE.
00215                   Y( IY ) = 0.0D+0
00216                ELSE IF ( Y( IY ) .EQ. 0.0D+0 ) THEN
00217                   SYMB_ZERO = .TRUE.
00218                ELSE
00219                   SYMB_ZERO = .FALSE.
00220                   Y( IY ) = BETA * ABS( Y( IY ) )
00221                END IF
00222                IF ( ALPHA .NE. 0.0D+0 ) THEN
00223                   DO J = 1, LENX
00224                      TEMP = CABS1( A( I, J ) )
00225                      SYMB_ZERO = SYMB_ZERO .AND.
00226      $                    ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
00227 
00228                      Y( IY ) = Y( IY ) + ALPHA*CABS1( X( J ) )*TEMP
00229                   END DO
00230                END IF
00231 
00232                IF ( .NOT.SYMB_ZERO ) Y( IY ) =
00233      $              Y( IY ) + SIGN( SAFE1, Y( IY ) )
00234 
00235                IY = IY + INCY
00236             END DO
00237          ELSE
00238             DO I = 1, LENY
00239                IF ( BETA .EQ. 0.0D+0 ) THEN
00240                   SYMB_ZERO = .TRUE.
00241                   Y( IY ) = 0.0D+0
00242                ELSE IF ( Y( IY ) .EQ. 0.0D+0 ) THEN
00243                   SYMB_ZERO = .TRUE.
00244                ELSE
00245                   SYMB_ZERO = .FALSE.
00246                   Y( IY ) = BETA * ABS( Y( IY ) )
00247                END IF
00248                IF ( ALPHA .NE. 0.0D+0 ) THEN
00249                   DO J = 1, LENX
00250                      TEMP = CABS1( A( J, I ) )
00251                      SYMB_ZERO = SYMB_ZERO .AND.
00252      $                    ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
00253 
00254                      Y( IY ) = Y( IY ) + ALPHA*CABS1( X( J ) )*TEMP
00255                   END DO
00256                END IF
00257 
00258                IF ( .NOT.SYMB_ZERO ) Y( IY ) =
00259      $              Y( IY ) + SIGN( SAFE1, Y( IY ) )
00260 
00261                IY = IY + INCY
00262             END DO
00263          END IF
00264       ELSE
00265          IF( TRANS.EQ.ILATRANS( 'N' ) )THEN
00266             DO I = 1, LENY
00267                IF ( BETA .EQ. 0.0D+0 ) THEN
00268                   SYMB_ZERO = .TRUE.
00269                   Y( IY ) = 0.0D+0
00270                ELSE IF ( Y( IY ) .EQ. 0.0D+0 ) THEN
00271                   SYMB_ZERO = .TRUE.
00272                ELSE
00273                   SYMB_ZERO = .FALSE.
00274                   Y( IY ) = BETA * ABS( Y( IY ) )
00275                END IF
00276                IF ( ALPHA .NE. 0.0D+0 ) THEN
00277                   JX = KX
00278                   DO J = 1, LENX
00279                      TEMP = CABS1( A( I, J ) )
00280                      SYMB_ZERO = SYMB_ZERO .AND.
00281      $                    ( X( JX ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
00282 
00283                      Y( IY ) = Y( IY ) + ALPHA*CABS1( X( JX ) )*TEMP
00284                      JX = JX + INCX
00285                   END DO
00286                END IF
00287 
00288                IF ( .NOT.SYMB_ZERO ) Y( IY ) =
00289      $              Y( IY ) + SIGN( SAFE1, Y( IY ) )
00290 
00291                IY = IY + INCY
00292             END DO
00293          ELSE
00294             DO I = 1, LENY
00295                IF ( BETA .EQ. 0.0D+0 ) THEN
00296                   SYMB_ZERO = .TRUE.
00297                   Y( IY ) = 0.0D+0
00298                ELSE IF ( Y( IY ) .EQ. 0.0D+0 ) THEN
00299                   SYMB_ZERO = .TRUE.
00300                ELSE
00301                   SYMB_ZERO = .FALSE.
00302                   Y( IY ) = BETA * ABS( Y( IY ) )
00303                END IF
00304                IF ( ALPHA .NE. 0.0D+0 ) THEN
00305                   JX = KX
00306                   DO J = 1, LENX
00307                      TEMP = CABS1( A( J, I ) )
00308                      SYMB_ZERO = SYMB_ZERO .AND.
00309      $                    ( X( JX ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
00310 
00311                      Y( IY ) = Y( IY ) + ALPHA*CABS1( X( JX ) )*TEMP
00312                      JX = JX + INCX
00313                   END DO
00314                END IF
00315 
00316                IF ( .NOT.SYMB_ZERO ) Y( IY ) =
00317      $              Y( IY ) + SIGN( SAFE1, Y( IY ) )
00318 
00319                IY = IY + INCY
00320             END DO
00321          END IF
00322 
00323       END IF
00324 *
00325       RETURN
00326 *
00327 *     End of ZLA_GEAMV
00328 *
00329       END