1*> \brief \b CHET01_ROOK 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 CHET01_ROOK( UPLO, N, A, LDA, AFAC, LDAFAC, IPIV, C, LDC, 12* RWORK, RESID ) 13* 14* .. Scalar Arguments .. 15* CHARACTER UPLO 16* INTEGER LDA, LDAFAC, LDC, N 17* REAL RESID 18* .. 19* .. Array Arguments .. 20* INTEGER IPIV( * ) 21* REAL RWORK( * ) 22* COMPLEX A( LDA, * ), AFAC( LDAFAC, * ), C( LDC, * ) 23* .. 24* 25* 26*> \par Purpose: 27* ============= 28*> 29*> \verbatim 30*> 31*> CHET01_ROOK reconstructs a complex Hermitian indefinite matrix A from its 32*> block L*D*L' or U*D*U' factorization and computes the residual 33*> norm( C - A ) / ( N * norm(A) * EPS ), 34*> where C is the reconstructed matrix, EPS is the machine epsilon, 35*> L' is the transpose of L, and U' is the transpose of U. 36*> \endverbatim 37* 38* Arguments: 39* ========== 40* 41*> \param[in] UPLO 42*> \verbatim 43*> UPLO is CHARACTER*1 44*> Specifies whether the upper or lower triangular part of the 45*> complex Hermitian matrix A is stored: 46*> = 'U': Upper triangular 47*> = 'L': Lower triangular 48*> \endverbatim 49*> 50*> \param[in] N 51*> \verbatim 52*> N is INTEGER 53*> The number of rows and columns of the matrix A. N >= 0. 54*> \endverbatim 55*> 56*> \param[in] A 57*> \verbatim 58*> A is COMPLEX array, dimension (LDA,N) 59*> The original complex Hermitian matrix A. 60*> \endverbatim 61*> 62*> \param[in] LDA 63*> \verbatim 64*> LDA is INTEGER 65*> The leading dimension of the array A. LDA >= max(1,N) 66*> \endverbatim 67*> 68*> \param[in] AFAC 69*> \verbatim 70*> AFAC is COMPLEX array, dimension (LDAFAC,N) 71*> The factored form of the matrix A. AFAC contains the block 72*> diagonal matrix D and the multipliers used to obtain the 73*> factor L or U from the block L*D*L' or U*D*U' factorization 74*> as computed by CSYTRF_ROOK. 75*> \endverbatim 76*> 77*> \param[in] LDAFAC 78*> \verbatim 79*> LDAFAC is INTEGER 80*> The leading dimension of the array AFAC. LDAFAC >= max(1,N). 81*> \endverbatim 82*> 83*> \param[in] IPIV 84*> \verbatim 85*> IPIV is INTEGER array, dimension (N) 86*> The pivot indices from CSYTRF_ROOK. 87*> \endverbatim 88*> 89*> \param[out] C 90*> \verbatim 91*> C is COMPLEX array, dimension (LDC,N) 92*> \endverbatim 93*> 94*> \param[in] LDC 95*> \verbatim 96*> LDC is INTEGER 97*> The leading dimension of the array C. LDC >= max(1,N). 98*> \endverbatim 99*> 100*> \param[out] RWORK 101*> \verbatim 102*> RWORK is REAL array, dimension (N) 103*> \endverbatim 104*> 105*> \param[out] RESID 106*> \verbatim 107*> RESID is REAL 108*> If UPLO = 'L', norm(L*D*L' - A) / ( N * norm(A) * EPS ) 109*> If UPLO = 'U', norm(U*D*U' - A) / ( N * norm(A) * 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 November 2013 121* 122*> \ingroup complex_lin 123* 124* ===================================================================== 125 SUBROUTINE CHET01_ROOK( UPLO, N, A, LDA, AFAC, LDAFAC, IPIV, C, 126 $ LDC, RWORK, RESID ) 127* 128* -- LAPACK test routine (version 3.5.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* November 2013 132* 133* .. Scalar Arguments .. 134 CHARACTER UPLO 135 INTEGER LDA, LDAFAC, LDC, N 136 REAL RESID 137* .. 138* .. Array Arguments .. 139 INTEGER IPIV( * ) 140 REAL RWORK( * ) 141 COMPLEX A( LDA, * ), AFAC( LDAFAC, * ), C( LDC, * ) 142* .. 143* 144* ===================================================================== 145* 146* .. Parameters .. 147 REAL ZERO, ONE 148 PARAMETER ( ZERO = 0.0E+0, ONE = 1.0E+0 ) 149 COMPLEX CZERO, CONE 150 PARAMETER ( CZERO = ( 0.0E+0, 0.0E+0 ), 151 $ CONE = ( 1.0E+0, 0.0E+0 ) ) 152* .. 153* .. Local Scalars .. 154 INTEGER I, INFO, J 155 REAL ANORM, EPS 156* .. 157* .. External Functions .. 158 LOGICAL LSAME 159 REAL CLANHE, SLAMCH 160 EXTERNAL LSAME, CLANHE, SLAMCH 161* .. 162* .. External Subroutines .. 163 EXTERNAL CLASET, CLAVHE_ROOK 164* .. 165* .. Intrinsic Functions .. 166 INTRINSIC AIMAG, REAL 167* .. 168* .. Executable Statements .. 169* 170* Quick exit if N = 0. 171* 172 IF( N.LE.0 ) THEN 173 RESID = ZERO 174 RETURN 175 END IF 176* 177* Determine EPS and the norm of A. 178* 179 EPS = SLAMCH( 'Epsilon' ) 180 ANORM = CLANHE( '1', UPLO, N, A, LDA, RWORK ) 181* 182* Check the imaginary parts of the diagonal elements and return with 183* an error code if any are nonzero. 184* 185 DO 10 J = 1, N 186 IF( AIMAG( AFAC( J, J ) ).NE.ZERO ) THEN 187 RESID = ONE / EPS 188 RETURN 189 END IF 190 10 CONTINUE 191* 192* Initialize C to the identity matrix. 193* 194 CALL CLASET( 'Full', N, N, CZERO, CONE, C, LDC ) 195* 196* Call CLAVHE_ROOK to form the product D * U' (or D * L' ). 197* 198 CALL CLAVHE_ROOK( UPLO, 'Conjugate', 'Non-unit', N, N, AFAC, 199 $ LDAFAC, IPIV, C, LDC, INFO ) 200* 201* Call CLAVHE_ROOK again to multiply by U (or L ). 202* 203 CALL CLAVHE_ROOK( UPLO, 'No transpose', 'Unit', N, N, AFAC, 204 $ LDAFAC, IPIV, C, LDC, INFO ) 205* 206* Compute the difference C - A . 207* 208 IF( LSAME( UPLO, 'U' ) ) THEN 209 DO 30 J = 1, N 210 DO 20 I = 1, J - 1 211 C( I, J ) = C( I, J ) - A( I, J ) 212 20 CONTINUE 213 C( J, J ) = C( J, J ) - REAL( A( J, J ) ) 214 30 CONTINUE 215 ELSE 216 DO 50 J = 1, N 217 C( J, J ) = C( J, J ) - REAL( A( J, J ) ) 218 DO 40 I = J + 1, N 219 C( I, J ) = C( I, J ) - A( I, J ) 220 40 CONTINUE 221 50 CONTINUE 222 END IF 223* 224* Compute norm( C - A ) / ( N * norm(A) * EPS ) 225* 226 RESID = CLANHE( '1', UPLO, N, C, LDC, RWORK ) 227* 228 IF( ANORM.LE.ZERO ) THEN 229 IF( RESID.NE.ZERO ) 230 $ RESID = ONE / EPS 231 ELSE 232 RESID = ( ( RESID/REAL( N ) )/ANORM ) / EPS 233 END IF 234* 235 RETURN 236* 237* End of CHET01_ROOK 238* 239 END 240