1*> \brief \b ZSYT03 2* 3* =========== DOCUMENTATION =========== 4* 5* Online html documentation available at 6* http://www.netlib.org/lapack/explore-html/ 7* 8* Definition: 9* =========== 10* 11* SUBROUTINE ZSYT03( UPLO, N, A, LDA, AINV, LDAINV, WORK, LDWORK, 12* RWORK, RCOND, RESID ) 13* 14* .. Scalar Arguments .. 15* CHARACTER UPLO 16* INTEGER LDA, LDAINV, LDWORK, N 17* DOUBLE PRECISION RCOND, RESID 18* .. 19* .. Array Arguments .. 20* DOUBLE PRECISION RWORK( * ) 21* COMPLEX*16 A( LDA, * ), AINV( LDAINV, * ), 22* $ WORK( LDWORK, * ) 23* .. 24* 25* 26*> \par Purpose: 27* ============= 28*> 29*> \verbatim 30*> 31*> ZSYT03 computes the residual for a complex symmetric matrix times 32*> its inverse: 33*> norm( I - A*AINV ) / ( N * norm(A) * norm(AINV) * EPS ) 34*> where EPS is the machine epsilon. 35*> \endverbatim 36* 37* Arguments: 38* ========== 39* 40*> \param[in] UPLO 41*> \verbatim 42*> UPLO is CHARACTER*1 43*> Specifies whether the upper or lower triangular part of the 44*> complex symmetric matrix A is stored: 45*> = 'U': Upper triangular 46*> = 'L': Lower triangular 47*> \endverbatim 48*> 49*> \param[in] N 50*> \verbatim 51*> N is INTEGER 52*> The number of rows and columns of the matrix A. N >= 0. 53*> \endverbatim 54*> 55*> \param[in] A 56*> \verbatim 57*> A is COMPLEX*16 array, dimension (LDA,N) 58*> The original complex symmetric matrix A. 59*> \endverbatim 60*> 61*> \param[in] LDA 62*> \verbatim 63*> LDA is INTEGER 64*> The leading dimension of the array A. LDA >= max(1,N) 65*> \endverbatim 66*> 67*> \param[in,out] AINV 68*> \verbatim 69*> AINV is COMPLEX*16 array, dimension (LDAINV,N) 70*> On entry, the inverse of the matrix A, stored as a symmetric 71*> matrix in the same format as A. 72*> In this version, AINV is expanded into a full matrix and 73*> multiplied by A, so the opposing triangle of AINV will be 74*> changed; i.e., if the upper triangular part of AINV is 75*> stored, the lower triangular part will be used as work space. 76*> \endverbatim 77*> 78*> \param[in] LDAINV 79*> \verbatim 80*> LDAINV is INTEGER 81*> The leading dimension of the array AINV. LDAINV >= max(1,N). 82*> \endverbatim 83*> 84*> \param[out] WORK 85*> \verbatim 86*> WORK is COMPLEX*16 array, dimension (LDWORK,N) 87*> \endverbatim 88*> 89*> \param[in] LDWORK 90*> \verbatim 91*> LDWORK is INTEGER 92*> The leading dimension of the array WORK. LDWORK >= max(1,N). 93*> \endverbatim 94*> 95*> \param[out] RWORK 96*> \verbatim 97*> RWORK is DOUBLE PRECISION array, dimension (N) 98*> \endverbatim 99*> 100*> \param[out] RCOND 101*> \verbatim 102*> RCOND is DOUBLE PRECISION 103*> The reciprocal of the condition number of A, computed as 104*> RCOND = 1/ (norm(A) * norm(AINV)). 105*> \endverbatim 106*> 107*> \param[out] RESID 108*> \verbatim 109*> RESID is DOUBLE PRECISION 110*> norm(I - A*AINV) / ( N * norm(A) * norm(AINV) * EPS ) 111*> \endverbatim 112* 113* Authors: 114* ======== 115* 116*> \author Univ. of Tennessee 117*> \author Univ. of California Berkeley 118*> \author Univ. of Colorado Denver 119*> \author NAG Ltd. 120* 121*> \date December 2016 122* 123*> \ingroup complex16_lin 124* 125* ===================================================================== 126 SUBROUTINE ZSYT03( UPLO, N, A, LDA, AINV, LDAINV, WORK, LDWORK, 127 $ RWORK, RCOND, RESID ) 128* 129* -- LAPACK test routine (version 3.7.0) -- 130* -- LAPACK is a software package provided by Univ. of Tennessee, -- 131* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 132* December 2016 133* 134* .. Scalar Arguments .. 135 CHARACTER UPLO 136 INTEGER LDA, LDAINV, LDWORK, N 137 DOUBLE PRECISION RCOND, RESID 138* .. 139* .. Array Arguments .. 140 DOUBLE PRECISION RWORK( * ) 141 COMPLEX*16 A( LDA, * ), AINV( LDAINV, * ), 142 $ WORK( LDWORK, * ) 143* .. 144* 145* ===================================================================== 146* 147* 148* .. Parameters .. 149 DOUBLE PRECISION ZERO, ONE 150 PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 ) 151 COMPLEX*16 CZERO, CONE 152 PARAMETER ( CZERO = ( 0.0D+0, 0.0D+0 ), 153 $ CONE = ( 1.0D+0, 0.0D+0 ) ) 154* .. 155* .. Local Scalars .. 156 INTEGER I, J 157 DOUBLE PRECISION AINVNM, ANORM, EPS 158* .. 159* .. External Functions .. 160 LOGICAL LSAME 161 DOUBLE PRECISION DLAMCH, ZLANGE, ZLANSY 162 EXTERNAL LSAME, DLAMCH, ZLANGE, ZLANSY 163* .. 164* .. External Subroutines .. 165 EXTERNAL ZSYMM 166* .. 167* .. Intrinsic Functions .. 168 INTRINSIC DBLE 169* .. 170* .. Executable Statements .. 171* 172* Quick exit if N = 0 173* 174 IF( N.LE.0 ) THEN 175 RCOND = ONE 176 RESID = ZERO 177 RETURN 178 END IF 179* 180* Exit with RESID = 1/EPS if ANORM = 0 or AINVNM = 0. 181* 182 EPS = DLAMCH( 'Epsilon' ) 183 ANORM = ZLANSY( '1', UPLO, N, A, LDA, RWORK ) 184 AINVNM = ZLANSY( '1', UPLO, N, AINV, LDAINV, RWORK ) 185 IF( ANORM.LE.ZERO .OR. AINVNM.LE.ZERO ) THEN 186 RCOND = ZERO 187 RESID = ONE / EPS 188 RETURN 189 END IF 190 RCOND = ( ONE / ANORM ) / AINVNM 191* 192* Expand AINV into a full matrix and call ZSYMM to multiply 193* AINV on the left by A (store the result in WORK). 194* 195 IF( LSAME( UPLO, 'U' ) ) THEN 196 DO 20 J = 1, N 197 DO 10 I = 1, J - 1 198 AINV( J, I ) = AINV( I, J ) 199 10 CONTINUE 200 20 CONTINUE 201 ELSE 202 DO 40 J = 1, N 203 DO 30 I = J + 1, N 204 AINV( J, I ) = AINV( I, J ) 205 30 CONTINUE 206 40 CONTINUE 207 END IF 208 CALL ZSYMM( 'Left', UPLO, N, N, -CONE, A, LDA, AINV, LDAINV, 209 $ CZERO, WORK, LDWORK ) 210* 211* Add the identity matrix to WORK . 212* 213 DO 50 I = 1, N 214 WORK( I, I ) = WORK( I, I ) + CONE 215 50 CONTINUE 216* 217* Compute norm(I - A*AINV) / (N * norm(A) * norm(AINV) * EPS) 218* 219 RESID = ZLANGE( '1', N, N, WORK, LDWORK, RWORK ) 220* 221 RESID = ( ( RESID*RCOND ) / EPS ) / DBLE( N ) 222* 223 RETURN 224* 225* End of ZSYT03 226* 227 END 228