Avoid out-of-bounds accesses in complex EIG tests

fixes Reference-LAPACK#333
martin-frbg · Apr 27, 2019 · ea2a102 · ea2a102
1 parent 12d8253
commit ea2a102
Show file tree

Hide file tree

Showing 4 changed files with 65 additions and 78 deletions.
diff --git a/TESTING/EIG/chet21.f b/TESTING/EIG/chet21.f
@@ -29,9 +29,8 @@
 *>
 *> CHET21 generally checks a decomposition of the form
 *>
-*>    A = U S U**H
-*>
-*> where **H means conjugate transpose, A is hermitian, U is unitary, and
+*>    A = U S UC>
+*> where * means conjugate transpose, A is hermitian, U is unitary, and
 *> S is diagonal (if KBAND=0) or (real) symmetric tridiagonal (if
 *> KBAND=1).
 *>
@@ -43,19 +42,18 @@
 *>
 *> Specifically, if ITYPE=1, then:
 *>
-*>    RESULT(1) = | A - U S U**H | / ( |A| n ulp ) and
-*>    RESULT(2) = | I - U U**H | / ( n ulp )
+*>    RESULT(1) = | A - U S U* | / ( |A| n ulp ) *andC>    RESULT(2) = | I - UU* | / ( n ulp )
 *>
 *> If ITYPE=2, then:
 *>
-*>    RESULT(1) = | A - V S V**H | / ( |A| n ulp )
+*>    RESULT(1) = | A - V S V* | / ( |A| n ulp )
 *>
 *> If ITYPE=3, then:
 *>
-*>    RESULT(1) = | I - U V**H | / ( n ulp )
+*>    RESULT(1) = | I - UV* | / ( n ulp )
 *>
 *> For ITYPE > 1, the transformation U is expressed as a product
-*> V = H(1)...H(n-2),  where H(j) = I  -  tau(j) v(j) v(j)**H and each
+*> V = H(1)...H(n-2),  where H(j) = I  -  tau(j) v(j) v(j)C> and each
 *> vector v(j) has its first j elements 0 and the remaining n-j elements
 *> stored in V(j+1:n,j).
 *> \endverbatim
@@ -68,15 +66,14 @@
 *>          ITYPE is INTEGER
 *>          Specifies the type of tests to be performed.
 *>          1: U expressed as a dense unitary matrix:
-*>             RESULT(1) = | A - U S U**H | / ( |A| n ulp ) and
-*>             RESULT(2) = | I - U U**H | / ( n ulp )
+*>             RESULT(1) = | A - U S U* | / ( |A| n ulp )   *andC>             RESULT(2) = | I - UU* | / ( n ulp )
 *>
 *>          2: U expressed as a product V of Housholder transformations:
-*>             RESULT(1) = | A - V S V**H | / ( |A| n ulp )
+*>             RESULT(1) = | A - V S V* | / ( |A| n ulp )
 *>
 *>          3: U expressed both as a dense unitary matrix and
 *>             as a product of Housholder transformations:
-*>             RESULT(1) = | I - U V**H | / ( n ulp )
+*>             RESULT(1) = | I - UV* | / ( n ulp )
 *> \endverbatim
 *>
 *> \param[in] UPLO
@@ -174,7 +171,7 @@
 *> \verbatim
 *>          TAU is COMPLEX array, dimension (N)
 *>          If ITYPE >= 2, then TAU(j) is the scalar factor of
-*>          v(j) v(j)**H in the Householder transformation H(j) of
+*>          v(j) v(j)* in the Householder transformation H(j) of
 *>          the product  U = H(1)...H(n-2)
 *>          If ITYPE < 2, then TAU is not referenced.
 *> \endverbatim
@@ -297,7 +294,7 @@ SUBROUTINE CHET21( ITYPE, UPLO, N, KBAND, A, LDA, D, E, U, LDU, V,
 *
       IF( ITYPE.EQ.1 ) THEN
 *
-*        ITYPE=1: error = A - U S U**H
+*        ITYPE=1: error = A - U S U*
 *
          CALL CLASET( 'Full', N, N, CZERO, CZERO, WORK, N )
          CALL CLACPY( CUPLO, N, N, A, LDA, WORK, N )
@@ -307,7 +304,8 @@ SUBROUTINE CHET21( ITYPE, UPLO, N, KBAND, A, LDA, D, E, U, LDU, V,
    10    CONTINUE
 *
          IF( N.GT.1 .AND. KBAND.EQ.1 ) THEN
-            DO 20 J = 1, N - 1
+CMK            DO 20 J = 1, N - 1
+            DO 20 J = 2, N - 1
                CALL CHER2( CUPLO, N, -CMPLX( E( J ) ), U( 1, J ), 1,
      $                     U( 1, J-1 ), 1, WORK, N )
    20       CONTINUE
@@ -316,7 +314,7 @@ SUBROUTINE CHET21( ITYPE, UPLO, N, KBAND, A, LDA, D, E, U, LDU, V,
 *
       ELSE IF( ITYPE.EQ.2 ) THEN
 *
-*        ITYPE=2: error = V S V**H - A
+*        ITYPE=2: error = V S V* - A
 *
          CALL CLASET( 'Full', N, N, CZERO, CZERO, WORK, N )
 *
@@ -373,7 +371,7 @@ SUBROUTINE CHET21( ITYPE, UPLO, N, KBAND, A, LDA, D, E, U, LDU, V,
 *
       ELSE IF( ITYPE.EQ.3 ) THEN
 *
-*        ITYPE=3: error = U V**H - I
+*        ITYPE=3: error = U V* - I
 *
          IF( N.LT.2 )
      $      RETURN
@@ -409,7 +407,7 @@ SUBROUTINE CHET21( ITYPE, UPLO, N, KBAND, A, LDA, D, E, U, LDU, V,
 *
 *     Do Test 2
 *
-*     Compute  U U**H - I
+*     Compute  UU* - I
 *
       IF( ITYPE.EQ.1 ) THEN
          CALL CGEMM( 'N', 'C', N, N, N, CONE, U, LDU, U, LDU, CZERO,

diff --git a/TESTING/EIG/chpt21.f b/TESTING/EIG/chpt21.f
@@ -29,9 +29,8 @@
 *>
 *> CHPT21  generally checks a decomposition of the form
 *>
-*>         A = U S U**H
-*>
-*> where **H means conjugate transpose, A is hermitian, U is
+*>         A = U S UC>
+*> where * means conjugate transpose, A is hermitian, U is
 *> unitary, and S is diagonal (if KBAND=0) or (real) symmetric
 *> tridiagonal (if KBAND=1).  If ITYPE=1, then U is represented as
 *> a dense matrix, otherwise the U is expressed as a product of
@@ -42,16 +41,15 @@
 *>
 *> Specifically, if ITYPE=1, then:
 *>
-*>         RESULT(1) = | A - U S U**H | / ( |A| n ulp ) and
-*>         RESULT(2) = | I - U U**H | / ( n ulp )
+*>         RESULT(1) = | A - U S U* | / ( |A| n ulp ) *andC>         RESULT(2) = | I - UU* | / ( n ulp )
 *>
 *> If ITYPE=2, then:
 *>
-*>         RESULT(1) = | A - V S V**H | / ( |A| n ulp )
+*>         RESULT(1) = | A - V S V* | / ( |A| n ulp )
 *>
 *> If ITYPE=3, then:
 *>
-*>         RESULT(1) = | I - U V**H | / ( n ulp )
+*>         RESULT(1) = | I - UV* | / ( n ulp )
 *>
 *> Packed storage means that, for example, if UPLO='U', then the columns
 *> of the upper triangle of A are stored one after another, so that
@@ -72,16 +70,14 @@
 *>
 *>    If UPLO='U', then  V = H(n-1)...H(1),  where
 *>
-*>        H(j) = I  -  tau(j) v(j) v(j)**H
-*>
+*>        H(j) = I  -  tau(j) v(j) v(j)C>
 *>    and the first j-1 elements of v(j) are stored in V(1:j-1,j+1),
 *>    (i.e., VP( j*(j+1)/2 + 1 : j*(j+1)/2 + j-1 ) ),
 *>    the j-th element is 1, and the last n-j elements are 0.
 *>
 *>    If UPLO='L', then  V = H(1)...H(n-1),  where
 *>
-*>        H(j) = I  -  tau(j) v(j) v(j)**H
-*>
+*>        H(j) = I  -  tau(j) v(j) v(j)C>
 *>    and the first j elements of v(j) are 0, the (j+1)-st is 1, and the
 *>    (j+2)-nd through n-th elements are stored in V(j+2:n,j) (i.e.,
 *>    in VP( (2*n-j)*(j-1)/2 + j+2 : (2*n-j)*(j-1)/2 + n ) .)
@@ -95,15 +91,14 @@
 *>          ITYPE is INTEGER
 *>          Specifies the type of tests to be performed.
 *>          1: U expressed as a dense unitary matrix:
-*>             RESULT(1) = | A - U S U**H | / ( |A| n ulp ) and
-*>             RESULT(2) = | I - U U**H | / ( n ulp )
+*>             RESULT(1) = | A - U S U* | / ( |A| n ulp )   *andC>             RESULT(2) = | I - UU* | / ( n ulp )
 *>
 *>          2: U expressed as a product V of Housholder transformations:
-*>             RESULT(1) = | A - V S V**H | / ( |A| n ulp )
+*>             RESULT(1) = | A - V S V* | / ( |A| n ulp )
 *>
 *>          3: U expressed both as a dense unitary matrix and
 *>             as a product of Housholder transformations:
-*>             RESULT(1) = | I - U V**H | / ( n ulp )
+*>             RESULT(1) = | I - UV* | / ( n ulp )
 *> \endverbatim
 *>
 *> \param[in] UPLO
@@ -186,7 +181,7 @@
 *> \verbatim
 *>          TAU is COMPLEX array, dimension (N)
 *>          If ITYPE >= 2, then TAU(j) is the scalar factor of
-*>          v(j) v(j)**H in the Householder transformation H(j) of
+*>          v(j) v(j)* in the Householder transformation H(j) of
 *>          the product  U = H(1)...H(n-2)
 *>          If ITYPE < 2, then TAU is not referenced.
 *> \endverbatim
@@ -318,7 +313,7 @@ SUBROUTINE CHPT21( ITYPE, UPLO, N, KBAND, AP, D, E, U, LDU, VP,
 *
       IF( ITYPE.EQ.1 ) THEN
 *
-*        ITYPE=1: error = A - U S U**H
+*        ITYPE=1: error = A - U S U*
 *
          CALL CLASET( 'Full', N, N, CZERO, CZERO, WORK, N )
          CALL CCOPY( LAP, AP, 1, WORK, 1 )
@@ -328,7 +323,7 @@ SUBROUTINE CHPT21( ITYPE, UPLO, N, KBAND, AP, D, E, U, LDU, VP,
    10    CONTINUE
 *
          IF( N.GT.1 .AND. KBAND.EQ.1 ) THEN
-            DO 20 J = 1, N - 1
+            DO 20 J = 2, N - 1
                CALL CHPR2( CUPLO, N, -CMPLX( E( J ) ), U( 1, J ), 1,
      $                     U( 1, J-1 ), 1, WORK )
    20       CONTINUE
@@ -337,7 +332,7 @@ SUBROUTINE CHPT21( ITYPE, UPLO, N, KBAND, AP, D, E, U, LDU, VP,
 *
       ELSE IF( ITYPE.EQ.2 ) THEN
 *
-*        ITYPE=2: error = V S V**H - A
+*        ITYPE=2: error = V S V* - A
 *
          CALL CLASET( 'Full', N, N, CZERO, CZERO, WORK, N )
 *
@@ -405,7 +400,7 @@ SUBROUTINE CHPT21( ITYPE, UPLO, N, KBAND, AP, D, E, U, LDU, VP,
 *
       ELSE IF( ITYPE.EQ.3 ) THEN
 *
-*        ITYPE=3: error = U V**H - I
+*        ITYPE=3: error = U V* - I
 *
          IF( N.LT.2 )
      $      RETURN
@@ -436,7 +431,7 @@ SUBROUTINE CHPT21( ITYPE, UPLO, N, KBAND, AP, D, E, U, LDU, VP,
 *
 *     Do Test 2
 *
-*     Compute  U U**H - I
+*     Compute  UU* - I
 *
       IF( ITYPE.EQ.1 ) THEN
          CALL CGEMM( 'N', 'C', N, N, N, CONE, U, LDU, U, LDU, CZERO,

diff --git a/TESTING/EIG/zhet21.f b/TESTING/EIG/zhet21.f
@@ -29,9 +29,8 @@
 *>
 *> ZHET21 generally checks a decomposition of the form
 *>
-*>    A = U S U**H
-*>
-*> where **H means conjugate transpose, A is hermitian, U is unitary, and
+*>    A = U S UC>
+*> where * means conjugate transpose, A is hermitian, U is unitary, and
 *> S is diagonal (if KBAND=0) or (real) symmetric tridiagonal (if
 *> KBAND=1).
 *>
@@ -43,19 +42,18 @@
 *>
 *> Specifically, if ITYPE=1, then:
 *>
-*>    RESULT(1) = | A - U S U**H | / ( |A| n ulp ) and
-*>    RESULT(2) = | I - U U**H | / ( n ulp )
+*>    RESULT(1) = | A - U S U* | / ( |A| n ulp ) *andC>    RESULT(2) = | I - UU* | / ( n ulp )
 *>
 *> If ITYPE=2, then:
 *>
-*>    RESULT(1) = | A - V S V**H | / ( |A| n ulp )
+*>    RESULT(1) = | A - V S V* | / ( |A| n ulp )
 *>
 *> If ITYPE=3, then:
 *>
-*>    RESULT(1) = | I - U V**H | / ( n ulp )
+*>    RESULT(1) = | I - UV* | / ( n ulp )
 *>
 *> For ITYPE > 1, the transformation U is expressed as a product
-*> V = H(1)...H(n-2),  where H(j) = I  -  tau(j) v(j) v(j)**H and each
+*> V = H(1)...H(n-2),  where H(j) = I  -  tau(j) v(j) v(j)C> and each
 *> vector v(j) has its first j elements 0 and the remaining n-j elements
 *> stored in V(j+1:n,j).
 *> \endverbatim
@@ -68,15 +66,14 @@
 *>          ITYPE is INTEGER
 *>          Specifies the type of tests to be performed.
 *>          1: U expressed as a dense unitary matrix:
-*>             RESULT(1) = | A - U S U**H | / ( |A| n ulp ) and
-*>             RESULT(2) = | I - U U**H | / ( n ulp )
+*>             RESULT(1) = | A - U S U* | / ( |A| n ulp )   *andC>             RESULT(2) = | I - UU* | / ( n ulp )
 *>
 *>          2: U expressed as a product V of Housholder transformations:
-*>             RESULT(1) = | A - V S V**H | / ( |A| n ulp )
+*>             RESULT(1) = | A - V S V* | / ( |A| n ulp )
 *>
 *>          3: U expressed both as a dense unitary matrix and
 *>             as a product of Housholder transformations:
-*>             RESULT(1) = | I - U V**H | / ( n ulp )
+*>             RESULT(1) = | I - UV* | / ( n ulp )
 *> \endverbatim
 *>
 *> \param[in] UPLO
@@ -174,7 +171,7 @@
 *> \verbatim
 *>          TAU is COMPLEX*16 array, dimension (N)
 *>          If ITYPE >= 2, then TAU(j) is the scalar factor of
-*>          v(j) v(j)**H in the Householder transformation H(j) of
+*>          v(j) v(j)* in the Householder transformation H(j) of
 *>          the product  U = H(1)...H(n-2)
 *>          If ITYPE < 2, then TAU is not referenced.
 *> \endverbatim
@@ -297,7 +294,7 @@ SUBROUTINE ZHET21( ITYPE, UPLO, N, KBAND, A, LDA, D, E, U, LDU, V,
 *
       IF( ITYPE.EQ.1 ) THEN
 *
-*        ITYPE=1: error = A - U S U**H
+*        ITYPE=1: error = A - U S U*
 *
          CALL ZLASET( 'Full', N, N, CZERO, CZERO, WORK, N )
          CALL ZLACPY( CUPLO, N, N, A, LDA, WORK, N )
@@ -307,7 +304,8 @@ SUBROUTINE ZHET21( ITYPE, UPLO, N, KBAND, A, LDA, D, E, U, LDU, V,
    10    CONTINUE
 *
          IF( N.GT.1 .AND. KBAND.EQ.1 ) THEN
-            DO 20 J = 1, N - 1
+CMK            DO 20 J = 1, N - 1
+            DO 20 J = 2, N - 1
                CALL ZHER2( CUPLO, N, -DCMPLX( E( J ) ), U( 1, J ), 1,
      $                     U( 1, J-1 ), 1, WORK, N )
    20       CONTINUE
@@ -316,7 +314,7 @@ SUBROUTINE ZHET21( ITYPE, UPLO, N, KBAND, A, LDA, D, E, U, LDU, V,
 *
       ELSE IF( ITYPE.EQ.2 ) THEN
 *
-*        ITYPE=2: error = V S V**H - A
+*        ITYPE=2: error = V S V* - A
 *
          CALL ZLASET( 'Full', N, N, CZERO, CZERO, WORK, N )
 *
@@ -373,7 +371,7 @@ SUBROUTINE ZHET21( ITYPE, UPLO, N, KBAND, A, LDA, D, E, U, LDU, V,
 *
       ELSE IF( ITYPE.EQ.3 ) THEN
 *
-*        ITYPE=3: error = U V**H - I
+*        ITYPE=3: error = U V* - I
 *
          IF( N.LT.2 )
      $      RETURN
@@ -409,7 +407,7 @@ SUBROUTINE ZHET21( ITYPE, UPLO, N, KBAND, A, LDA, D, E, U, LDU, V,
 *
 *     Do Test 2
 *
-*     Compute  U U**H - I
+*     Compute  UU* - I
 *
       IF( ITYPE.EQ.1 ) THEN
          CALL ZGEMM( 'N', 'C', N, N, N, CONE, U, LDU, U, LDU, CZERO,