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