1*> \brief \b SGTT01 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 SGTT01( N, DL, D, DU, DLF, DF, DUF, DU2, IPIV, WORK, 12* LDWORK, RWORK, RESID ) 13* 14* .. Scalar Arguments .. 15* INTEGER LDWORK, N 16* REAL RESID 17* .. 18* .. Array Arguments .. 19* INTEGER IPIV( * ) 20* REAL D( * ), DF( * ), DL( * ), DLF( * ), DU( * ), 21* $ DU2( * ), DUF( * ), RWORK( * ), 22* $ WORK( LDWORK, * ) 23* .. 24* 25* 26*> \par Purpose: 27* ============= 28*> 29*> \verbatim 30*> 31*> SGTT01 reconstructs a tridiagonal matrix A from its LU factorization 32*> and computes the residual 33*> norm(L*U - A) / ( norm(A) * EPS ), 34*> where EPS is the machine epsilon. 35*> \endverbatim 36* 37* Arguments: 38* ========== 39* 40*> \param[in] N 41*> \verbatim 42*> N is INTEGTER 43*> The order of the matrix A. N >= 0. 44*> \endverbatim 45*> 46*> \param[in] DL 47*> \verbatim 48*> DL is REAL array, dimension (N-1) 49*> The (n-1) sub-diagonal elements of A. 50*> \endverbatim 51*> 52*> \param[in] D 53*> \verbatim 54*> D is REAL array, dimension (N) 55*> The diagonal elements of A. 56*> \endverbatim 57*> 58*> \param[in] DU 59*> \verbatim 60*> DU is REAL array, dimension (N-1) 61*> The (n-1) super-diagonal elements of A. 62*> \endverbatim 63*> 64*> \param[in] DLF 65*> \verbatim 66*> DLF is REAL array, dimension (N-1) 67*> The (n-1) multipliers that define the matrix L from the 68*> LU factorization of A. 69*> \endverbatim 70*> 71*> \param[in] DF 72*> \verbatim 73*> DF is REAL array, dimension (N) 74*> The n diagonal elements of the upper triangular matrix U from 75*> the LU factorization of A. 76*> \endverbatim 77*> 78*> \param[in] DUF 79*> \verbatim 80*> DUF is REAL array, dimension (N-1) 81*> The (n-1) elements of the first super-diagonal of U. 82*> \endverbatim 83*> 84*> \param[in] DU2 85*> \verbatim 86*> DU2 is REAL array, dimension (N-2) 87*> The (n-2) elements of the second super-diagonal of U. 88*> \endverbatim 89*> 90*> \param[in] IPIV 91*> \verbatim 92*> IPIV is INTEGER array, dimension (N) 93*> The pivot indices; for 1 <= i <= n, row i of the matrix was 94*> interchanged with row IPIV(i). IPIV(i) will always be either 95*> i or i+1; IPIV(i) = i indicates a row interchange was not 96*> required. 97*> \endverbatim 98*> 99*> \param[out] WORK 100*> \verbatim 101*> WORK is REAL array, dimension (LDWORK,N) 102*> \endverbatim 103*> 104*> \param[in] LDWORK 105*> \verbatim 106*> LDWORK is INTEGER 107*> The leading dimension of the array WORK. LDWORK >= max(1,N). 108*> \endverbatim 109*> 110*> \param[out] RWORK 111*> \verbatim 112*> RWORK is REAL array, dimension (N) 113*> \endverbatim 114*> 115*> \param[out] RESID 116*> \verbatim 117*> RESID is REAL 118*> The scaled residual: norm(L*U - A) / (norm(A) * EPS) 119*> \endverbatim 120* 121* Authors: 122* ======== 123* 124*> \author Univ. of Tennessee 125*> \author Univ. of California Berkeley 126*> \author Univ. of Colorado Denver 127*> \author NAG Ltd. 128* 129*> \ingroup single_lin 130* 131* ===================================================================== 132 SUBROUTINE SGTT01( N, DL, D, DU, DLF, DF, DUF, DU2, IPIV, WORK, 133 $ LDWORK, RWORK, RESID ) 134* 135* -- LAPACK test routine -- 136* -- LAPACK is a software package provided by Univ. of Tennessee, -- 137* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 138* 139* .. Scalar Arguments .. 140 INTEGER LDWORK, N 141 REAL RESID 142* .. 143* .. Array Arguments .. 144 INTEGER IPIV( * ) 145 REAL D( * ), DF( * ), DL( * ), DLF( * ), DU( * ), 146 $ DU2( * ), DUF( * ), RWORK( * ), 147 $ WORK( LDWORK, * ) 148* .. 149* 150* ===================================================================== 151* 152* .. Parameters .. 153 REAL ONE, ZERO 154 PARAMETER ( ONE = 1.0E+0, ZERO = 0.0E+0 ) 155* .. 156* .. Local Scalars .. 157 INTEGER I, IP, J, LASTJ 158 REAL ANORM, EPS, LI 159* .. 160* .. External Functions .. 161 REAL SLAMCH, SLANGT, SLANHS 162 EXTERNAL SLAMCH, SLANGT, SLANHS 163* .. 164* .. Intrinsic Functions .. 165 INTRINSIC MIN 166* .. 167* .. External Subroutines .. 168 EXTERNAL SAXPY, SSWAP 169* .. 170* .. Executable Statements .. 171* 172* Quick return if possible 173* 174 IF( N.LE.0 ) THEN 175 RESID = ZERO 176 RETURN 177 END IF 178* 179 EPS = SLAMCH( 'Epsilon' ) 180* 181* Copy the matrix U to WORK. 182* 183 DO 20 J = 1, N 184 DO 10 I = 1, N 185 WORK( I, J ) = ZERO 186 10 CONTINUE 187 20 CONTINUE 188 DO 30 I = 1, N 189 IF( I.EQ.1 ) THEN 190 WORK( I, I ) = DF( I ) 191 IF( N.GE.2 ) 192 $ WORK( I, I+1 ) = DUF( I ) 193 IF( N.GE.3 ) 194 $ WORK( I, I+2 ) = DU2( I ) 195 ELSE IF( I.EQ.N ) THEN 196 WORK( I, I ) = DF( I ) 197 ELSE 198 WORK( I, I ) = DF( I ) 199 WORK( I, I+1 ) = DUF( I ) 200 IF( I.LT.N-1 ) 201 $ WORK( I, I+2 ) = DU2( I ) 202 END IF 203 30 CONTINUE 204* 205* Multiply on the left by L. 206* 207 LASTJ = N 208 DO 40 I = N - 1, 1, -1 209 LI = DLF( I ) 210 CALL SAXPY( LASTJ-I+1, LI, WORK( I, I ), LDWORK, 211 $ WORK( I+1, I ), LDWORK ) 212 IP = IPIV( I ) 213 IF( IP.EQ.I ) THEN 214 LASTJ = MIN( I+2, N ) 215 ELSE 216 CALL SSWAP( LASTJ-I+1, WORK( I, I ), LDWORK, WORK( I+1, I ), 217 $ LDWORK ) 218 END IF 219 40 CONTINUE 220* 221* Subtract the matrix A. 222* 223 WORK( 1, 1 ) = WORK( 1, 1 ) - D( 1 ) 224 IF( N.GT.1 ) THEN 225 WORK( 1, 2 ) = WORK( 1, 2 ) - DU( 1 ) 226 WORK( N, N-1 ) = WORK( N, N-1 ) - DL( N-1 ) 227 WORK( N, N ) = WORK( N, N ) - D( N ) 228 DO 50 I = 2, N - 1 229 WORK( I, I-1 ) = WORK( I, I-1 ) - DL( I-1 ) 230 WORK( I, I ) = WORK( I, I ) - D( I ) 231 WORK( I, I+1 ) = WORK( I, I+1 ) - DU( I ) 232 50 CONTINUE 233 END IF 234* 235* Compute the 1-norm of the tridiagonal matrix A. 236* 237 ANORM = SLANGT( '1', N, DL, D, DU ) 238* 239* Compute the 1-norm of WORK, which is only guaranteed to be 240* upper Hessenberg. 241* 242 RESID = SLANHS( '1', N, WORK, LDWORK, RWORK ) 243* 244* Compute norm(L*U - A) / (norm(A) * EPS) 245* 246 IF( ANORM.LE.ZERO ) THEN 247 IF( RESID.NE.ZERO ) 248 $ RESID = ONE / EPS 249 ELSE 250 RESID = ( RESID / ANORM ) / EPS 251 END IF 252* 253 RETURN 254* 255* End of SGTT01 256* 257 END 258