![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b DLATTR 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 * Definition: 00009 * =========== 00010 * 00011 * SUBROUTINE DLATTR( IMAT, UPLO, TRANS, DIAG, ISEED, N, A, LDA, B, 00012 * WORK, INFO ) 00013 * 00014 * .. Scalar Arguments .. 00015 * CHARACTER DIAG, TRANS, UPLO 00016 * INTEGER IMAT, INFO, LDA, N 00017 * .. 00018 * .. Array Arguments .. 00019 * INTEGER ISEED( 4 ) 00020 * DOUBLE PRECISION A( LDA, * ), B( * ), WORK( * ) 00021 * .. 00022 * 00023 * 00024 *> \par Purpose: 00025 * ============= 00026 *> 00027 *> \verbatim 00028 *> 00029 *> DLATTR generates a triangular test matrix. 00030 *> IMAT and UPLO uniquely specify the properties of the test 00031 *> matrix, which is returned in the array A. 00032 *> \endverbatim 00033 * 00034 * Arguments: 00035 * ========== 00036 * 00037 *> \param[in] IMAT 00038 *> \verbatim 00039 *> IMAT is INTEGER 00040 *> An integer key describing which matrix to generate for this 00041 *> path. 00042 *> \endverbatim 00043 *> 00044 *> \param[in] UPLO 00045 *> \verbatim 00046 *> UPLO is CHARACTER*1 00047 *> Specifies whether the matrix A will be upper or lower 00048 *> triangular. 00049 *> = 'U': Upper triangular 00050 *> = 'L': Lower triangular 00051 *> \endverbatim 00052 *> 00053 *> \param[in] TRANS 00054 *> \verbatim 00055 *> TRANS is CHARACTER*1 00056 *> Specifies whether the matrix or its transpose will be used. 00057 *> = 'N': No transpose 00058 *> = 'T': Transpose 00059 *> = 'C': Conjugate transpose (= Transpose) 00060 *> \endverbatim 00061 *> 00062 *> \param[out] DIAG 00063 *> \verbatim 00064 *> DIAG is CHARACTER*1 00065 *> Specifies whether or not the matrix A is unit triangular. 00066 *> = 'N': Non-unit triangular 00067 *> = 'U': Unit triangular 00068 *> \endverbatim 00069 *> 00070 *> \param[in,out] ISEED 00071 *> \verbatim 00072 *> ISEED is INTEGER array, dimension (4) 00073 *> The seed vector for the random number generator (used in 00074 *> DLATMS). Modified on exit. 00075 *> \endverbatim 00076 *> 00077 *> \param[in] N 00078 *> \verbatim 00079 *> N is INTEGER 00080 *> The order of the matrix to be generated. 00081 *> \endverbatim 00082 *> 00083 *> \param[out] A 00084 *> \verbatim 00085 *> A is DOUBLE PRECISION array, dimension (LDA,N) 00086 *> The triangular matrix A. If UPLO = 'U', the leading n by n 00087 *> upper triangular part of the array A contains the upper 00088 *> triangular matrix, and the strictly lower triangular part of 00089 *> A is not referenced. If UPLO = 'L', the leading n by n lower 00090 *> triangular part of the array A contains the lower triangular 00091 *> matrix, and the strictly upper triangular part of A is not 00092 *> referenced. If DIAG = 'U', the diagonal elements of A are 00093 *> set so that A(k,k) = k for 1 <= k <= n. 00094 *> \endverbatim 00095 *> 00096 *> \param[in] LDA 00097 *> \verbatim 00098 *> LDA is INTEGER 00099 *> The leading dimension of the array A. LDA >= max(1,N). 00100 *> \endverbatim 00101 *> 00102 *> \param[out] B 00103 *> \verbatim 00104 *> B is DOUBLE PRECISION array, dimension (N) 00105 *> The right hand side vector, if IMAT > 10. 00106 *> \endverbatim 00107 *> 00108 *> \param[out] WORK 00109 *> \verbatim 00110 *> WORK is DOUBLE PRECISION array, dimension (3*N) 00111 *> \endverbatim 00112 *> 00113 *> \param[out] INFO 00114 *> \verbatim 00115 *> INFO is INTEGER 00116 *> = 0: successful exit 00117 *> < 0: if INFO = -k, the k-th argument had an illegal value 00118 *> \endverbatim 00119 * 00120 * Authors: 00121 * ======== 00122 * 00123 *> \author Univ. of Tennessee 00124 *> \author Univ. of California Berkeley 00125 *> \author Univ. of Colorado Denver 00126 *> \author NAG Ltd. 00127 * 00128 *> \date November 2011 00129 * 00130 *> \ingroup double_lin 00131 * 00132 * ===================================================================== 00133 SUBROUTINE DLATTR( IMAT, UPLO, TRANS, DIAG, ISEED, N, A, LDA, B, 00134 $ WORK, INFO ) 00135 * 00136 * -- LAPACK test routine (version 3.4.0) -- 00137 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00138 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00139 * November 2011 00140 * 00141 * .. Scalar Arguments .. 00142 CHARACTER DIAG, TRANS, UPLO 00143 INTEGER IMAT, INFO, LDA, N 00144 * .. 00145 * .. Array Arguments .. 00146 INTEGER ISEED( 4 ) 00147 DOUBLE PRECISION A( LDA, * ), B( * ), WORK( * ) 00148 * .. 00149 * 00150 * ===================================================================== 00151 * 00152 * .. Parameters .. 00153 DOUBLE PRECISION ONE, TWO, ZERO 00154 PARAMETER ( ONE = 1.0D+0, TWO = 2.0D+0, ZERO = 0.0D+0 ) 00155 * .. 00156 * .. Local Scalars .. 00157 LOGICAL UPPER 00158 CHARACTER DIST, TYPE 00159 CHARACTER*3 PATH 00160 INTEGER I, IY, J, JCOUNT, KL, KU, MODE 00161 DOUBLE PRECISION ANORM, BIGNUM, BNORM, BSCAL, C, CNDNUM, PLUS1, 00162 $ PLUS2, RA, RB, REXP, S, SFAC, SMLNUM, STAR1, 00163 $ TEXP, TLEFT, TSCAL, ULP, UNFL, X, Y, Z 00164 * .. 00165 * .. External Functions .. 00166 LOGICAL LSAME 00167 INTEGER IDAMAX 00168 DOUBLE PRECISION DLAMCH, DLARND 00169 EXTERNAL LSAME, IDAMAX, DLAMCH, DLARND 00170 * .. 00171 * .. External Subroutines .. 00172 EXTERNAL DCOPY, DLABAD, DLARNV, DLATB4, DLATMS, DROT, 00173 $ DROTG, DSCAL, DSWAP 00174 * .. 00175 * .. Intrinsic Functions .. 00176 INTRINSIC ABS, DBLE, MAX, SIGN, SQRT 00177 * .. 00178 * .. Executable Statements .. 00179 * 00180 PATH( 1: 1 ) = 'Double precision' 00181 PATH( 2: 3 ) = 'TR' 00182 UNFL = DLAMCH( 'Safe minimum' ) 00183 ULP = DLAMCH( 'Epsilon' )*DLAMCH( 'Base' ) 00184 SMLNUM = UNFL 00185 BIGNUM = ( ONE-ULP ) / SMLNUM 00186 CALL DLABAD( SMLNUM, BIGNUM ) 00187 IF( ( IMAT.GE.7 .AND. IMAT.LE.10 ) .OR. IMAT.EQ.18 ) THEN 00188 DIAG = 'U' 00189 ELSE 00190 DIAG = 'N' 00191 END IF 00192 INFO = 0 00193 * 00194 * Quick return if N.LE.0. 00195 * 00196 IF( N.LE.0 ) 00197 $ RETURN 00198 * 00199 * Call DLATB4 to set parameters for SLATMS. 00200 * 00201 UPPER = LSAME( UPLO, 'U' ) 00202 IF( UPPER ) THEN 00203 CALL DLATB4( PATH, IMAT, N, N, TYPE, KL, KU, ANORM, MODE, 00204 $ CNDNUM, DIST ) 00205 ELSE 00206 CALL DLATB4( PATH, -IMAT, N, N, TYPE, KL, KU, ANORM, MODE, 00207 $ CNDNUM, DIST ) 00208 END IF 00209 * 00210 * IMAT <= 6: Non-unit triangular matrix 00211 * 00212 IF( IMAT.LE.6 ) THEN 00213 CALL DLATMS( N, N, DIST, ISEED, TYPE, B, MODE, CNDNUM, ANORM, 00214 $ KL, KU, 'No packing', A, LDA, WORK, INFO ) 00215 * 00216 * IMAT > 6: Unit triangular matrix 00217 * The diagonal is deliberately set to something other than 1. 00218 * 00219 * IMAT = 7: Matrix is the identity 00220 * 00221 ELSE IF( IMAT.EQ.7 ) THEN 00222 IF( UPPER ) THEN 00223 DO 20 J = 1, N 00224 DO 10 I = 1, J - 1 00225 A( I, J ) = ZERO 00226 10 CONTINUE 00227 A( J, J ) = J 00228 20 CONTINUE 00229 ELSE 00230 DO 40 J = 1, N 00231 A( J, J ) = J 00232 DO 30 I = J + 1, N 00233 A( I, J ) = ZERO 00234 30 CONTINUE 00235 40 CONTINUE 00236 END IF 00237 * 00238 * IMAT > 7: Non-trivial unit triangular matrix 00239 * 00240 * Generate a unit triangular matrix T with condition CNDNUM by 00241 * forming a triangular matrix with known singular values and 00242 * filling in the zero entries with Givens rotations. 00243 * 00244 ELSE IF( IMAT.LE.10 ) THEN 00245 IF( UPPER ) THEN 00246 DO 60 J = 1, N 00247 DO 50 I = 1, J - 1 00248 A( I, J ) = ZERO 00249 50 CONTINUE 00250 A( J, J ) = J 00251 60 CONTINUE 00252 ELSE 00253 DO 80 J = 1, N 00254 A( J, J ) = J 00255 DO 70 I = J + 1, N 00256 A( I, J ) = ZERO 00257 70 CONTINUE 00258 80 CONTINUE 00259 END IF 00260 * 00261 * Since the trace of a unit triangular matrix is 1, the product 00262 * of its singular values must be 1. Let s = sqrt(CNDNUM), 00263 * x = sqrt(s) - 1/sqrt(s), y = sqrt(2/(n-2))*x, and z = x**2. 00264 * The following triangular matrix has singular values s, 1, 1, 00265 * ..., 1, 1/s: 00266 * 00267 * 1 y y y ... y y z 00268 * 1 0 0 ... 0 0 y 00269 * 1 0 ... 0 0 y 00270 * . ... . . . 00271 * . . . . 00272 * 1 0 y 00273 * 1 y 00274 * 1 00275 * 00276 * To fill in the zeros, we first multiply by a matrix with small 00277 * condition number of the form 00278 * 00279 * 1 0 0 0 0 ... 00280 * 1 + * 0 0 ... 00281 * 1 + 0 0 0 00282 * 1 + * 0 0 00283 * 1 + 0 0 00284 * ... 00285 * 1 + 0 00286 * 1 0 00287 * 1 00288 * 00289 * Each element marked with a '*' is formed by taking the product 00290 * of the adjacent elements marked with '+'. The '*'s can be 00291 * chosen freely, and the '+'s are chosen so that the inverse of 00292 * T will have elements of the same magnitude as T. If the *'s in 00293 * both T and inv(T) have small magnitude, T is well conditioned. 00294 * The two offdiagonals of T are stored in WORK. 00295 * 00296 * The product of these two matrices has the form 00297 * 00298 * 1 y y y y y . y y z 00299 * 1 + * 0 0 . 0 0 y 00300 * 1 + 0 0 . 0 0 y 00301 * 1 + * . . . . 00302 * 1 + . . . . 00303 * . . . . . 00304 * . . . . 00305 * 1 + y 00306 * 1 y 00307 * 1 00308 * 00309 * Now we multiply by Givens rotations, using the fact that 00310 * 00311 * [ c s ] [ 1 w ] [ -c -s ] = [ 1 -w ] 00312 * [ -s c ] [ 0 1 ] [ s -c ] [ 0 1 ] 00313 * and 00314 * [ -c -s ] [ 1 0 ] [ c s ] = [ 1 0 ] 00315 * [ s -c ] [ w 1 ] [ -s c ] [ -w 1 ] 00316 * 00317 * where c = w / sqrt(w**2+4) and s = 2 / sqrt(w**2+4). 00318 * 00319 STAR1 = 0.25D0 00320 SFAC = 0.5D0 00321 PLUS1 = SFAC 00322 DO 90 J = 1, N, 2 00323 PLUS2 = STAR1 / PLUS1 00324 WORK( J ) = PLUS1 00325 WORK( N+J ) = STAR1 00326 IF( J+1.LE.N ) THEN 00327 WORK( J+1 ) = PLUS2 00328 WORK( N+J+1 ) = ZERO 00329 PLUS1 = STAR1 / PLUS2 00330 REXP = DLARND( 2, ISEED ) 00331 STAR1 = STAR1*( SFAC**REXP ) 00332 IF( REXP.LT.ZERO ) THEN 00333 STAR1 = -SFAC**( ONE-REXP ) 00334 ELSE 00335 STAR1 = SFAC**( ONE+REXP ) 00336 END IF 00337 END IF 00338 90 CONTINUE 00339 * 00340 X = SQRT( CNDNUM ) - 1 / SQRT( CNDNUM ) 00341 IF( N.GT.2 ) THEN 00342 Y = SQRT( 2.D0 / ( N-2 ) )*X 00343 ELSE 00344 Y = ZERO 00345 END IF 00346 Z = X*X 00347 * 00348 IF( UPPER ) THEN 00349 IF( N.GT.3 ) THEN 00350 CALL DCOPY( N-3, WORK, 1, A( 2, 3 ), LDA+1 ) 00351 IF( N.GT.4 ) 00352 $ CALL DCOPY( N-4, WORK( N+1 ), 1, A( 2, 4 ), LDA+1 ) 00353 END IF 00354 DO 100 J = 2, N - 1 00355 A( 1, J ) = Y 00356 A( J, N ) = Y 00357 100 CONTINUE 00358 A( 1, N ) = Z 00359 ELSE 00360 IF( N.GT.3 ) THEN 00361 CALL DCOPY( N-3, WORK, 1, A( 3, 2 ), LDA+1 ) 00362 IF( N.GT.4 ) 00363 $ CALL DCOPY( N-4, WORK( N+1 ), 1, A( 4, 2 ), LDA+1 ) 00364 END IF 00365 DO 110 J = 2, N - 1 00366 A( J, 1 ) = Y 00367 A( N, J ) = Y 00368 110 CONTINUE 00369 A( N, 1 ) = Z 00370 END IF 00371 * 00372 * Fill in the zeros using Givens rotations. 00373 * 00374 IF( UPPER ) THEN 00375 DO 120 J = 1, N - 1 00376 RA = A( J, J+1 ) 00377 RB = 2.0D0 00378 CALL DROTG( RA, RB, C, S ) 00379 * 00380 * Multiply by [ c s; -s c] on the left. 00381 * 00382 IF( N.GT.J+1 ) 00383 $ CALL DROT( N-J-1, A( J, J+2 ), LDA, A( J+1, J+2 ), 00384 $ LDA, C, S ) 00385 * 00386 * Multiply by [-c -s; s -c] on the right. 00387 * 00388 IF( J.GT.1 ) 00389 $ CALL DROT( J-1, A( 1, J+1 ), 1, A( 1, J ), 1, -C, -S ) 00390 * 00391 * Negate A(J,J+1). 00392 * 00393 A( J, J+1 ) = -A( J, J+1 ) 00394 120 CONTINUE 00395 ELSE 00396 DO 130 J = 1, N - 1 00397 RA = A( J+1, J ) 00398 RB = 2.0D0 00399 CALL DROTG( RA, RB, C, S ) 00400 * 00401 * Multiply by [ c -s; s c] on the right. 00402 * 00403 IF( N.GT.J+1 ) 00404 $ CALL DROT( N-J-1, A( J+2, J+1 ), 1, A( J+2, J ), 1, C, 00405 $ -S ) 00406 * 00407 * Multiply by [-c s; -s -c] on the left. 00408 * 00409 IF( J.GT.1 ) 00410 $ CALL DROT( J-1, A( J, 1 ), LDA, A( J+1, 1 ), LDA, -C, 00411 $ S ) 00412 * 00413 * Negate A(J+1,J). 00414 * 00415 A( J+1, J ) = -A( J+1, J ) 00416 130 CONTINUE 00417 END IF 00418 * 00419 * IMAT > 10: Pathological test cases. These triangular matrices 00420 * are badly scaled or badly conditioned, so when used in solving a 00421 * triangular system they may cause overflow in the solution vector. 00422 * 00423 ELSE IF( IMAT.EQ.11 ) THEN 00424 * 00425 * Type 11: Generate a triangular matrix with elements between 00426 * -1 and 1. Give the diagonal norm 2 to make it well-conditioned. 00427 * Make the right hand side large so that it requires scaling. 00428 * 00429 IF( UPPER ) THEN 00430 DO 140 J = 1, N 00431 CALL DLARNV( 2, ISEED, J, A( 1, J ) ) 00432 A( J, J ) = SIGN( TWO, A( J, J ) ) 00433 140 CONTINUE 00434 ELSE 00435 DO 150 J = 1, N 00436 CALL DLARNV( 2, ISEED, N-J+1, A( J, J ) ) 00437 A( J, J ) = SIGN( TWO, A( J, J ) ) 00438 150 CONTINUE 00439 END IF 00440 * 00441 * Set the right hand side so that the largest value is BIGNUM. 00442 * 00443 CALL DLARNV( 2, ISEED, N, B ) 00444 IY = IDAMAX( N, B, 1 ) 00445 BNORM = ABS( B( IY ) ) 00446 BSCAL = BIGNUM / MAX( ONE, BNORM ) 00447 CALL DSCAL( N, BSCAL, B, 1 ) 00448 * 00449 ELSE IF( IMAT.EQ.12 ) THEN 00450 * 00451 * Type 12: Make the first diagonal element in the solve small to 00452 * cause immediate overflow when dividing by T(j,j). 00453 * In type 12, the offdiagonal elements are small (CNORM(j) < 1). 00454 * 00455 CALL DLARNV( 2, ISEED, N, B ) 00456 TSCAL = ONE / MAX( ONE, DBLE( N-1 ) ) 00457 IF( UPPER ) THEN 00458 DO 160 J = 1, N 00459 CALL DLARNV( 2, ISEED, J, A( 1, J ) ) 00460 CALL DSCAL( J-1, TSCAL, A( 1, J ), 1 ) 00461 A( J, J ) = SIGN( ONE, A( J, J ) ) 00462 160 CONTINUE 00463 A( N, N ) = SMLNUM*A( N, N ) 00464 ELSE 00465 DO 170 J = 1, N 00466 CALL DLARNV( 2, ISEED, N-J+1, A( J, J ) ) 00467 IF( N.GT.J ) 00468 $ CALL DSCAL( N-J, TSCAL, A( J+1, J ), 1 ) 00469 A( J, J ) = SIGN( ONE, A( J, J ) ) 00470 170 CONTINUE 00471 A( 1, 1 ) = SMLNUM*A( 1, 1 ) 00472 END IF 00473 * 00474 ELSE IF( IMAT.EQ.13 ) THEN 00475 * 00476 * Type 13: Make the first diagonal element in the solve small to 00477 * cause immediate overflow when dividing by T(j,j). 00478 * In type 13, the offdiagonal elements are O(1) (CNORM(j) > 1). 00479 * 00480 CALL DLARNV( 2, ISEED, N, B ) 00481 IF( UPPER ) THEN 00482 DO 180 J = 1, N 00483 CALL DLARNV( 2, ISEED, J, A( 1, J ) ) 00484 A( J, J ) = SIGN( ONE, A( J, J ) ) 00485 180 CONTINUE 00486 A( N, N ) = SMLNUM*A( N, N ) 00487 ELSE 00488 DO 190 J = 1, N 00489 CALL DLARNV( 2, ISEED, N-J+1, A( J, J ) ) 00490 A( J, J ) = SIGN( ONE, A( J, J ) ) 00491 190 CONTINUE 00492 A( 1, 1 ) = SMLNUM*A( 1, 1 ) 00493 END IF 00494 * 00495 ELSE IF( IMAT.EQ.14 ) THEN 00496 * 00497 * Type 14: T is diagonal with small numbers on the diagonal to 00498 * make the growth factor underflow, but a small right hand side 00499 * chosen so that the solution does not overflow. 00500 * 00501 IF( UPPER ) THEN 00502 JCOUNT = 1 00503 DO 210 J = N, 1, -1 00504 DO 200 I = 1, J - 1 00505 A( I, J ) = ZERO 00506 200 CONTINUE 00507 IF( JCOUNT.LE.2 ) THEN 00508 A( J, J ) = SMLNUM 00509 ELSE 00510 A( J, J ) = ONE 00511 END IF 00512 JCOUNT = JCOUNT + 1 00513 IF( JCOUNT.GT.4 ) 00514 $ JCOUNT = 1 00515 210 CONTINUE 00516 ELSE 00517 JCOUNT = 1 00518 DO 230 J = 1, N 00519 DO 220 I = J + 1, N 00520 A( I, J ) = ZERO 00521 220 CONTINUE 00522 IF( JCOUNT.LE.2 ) THEN 00523 A( J, J ) = SMLNUM 00524 ELSE 00525 A( J, J ) = ONE 00526 END IF 00527 JCOUNT = JCOUNT + 1 00528 IF( JCOUNT.GT.4 ) 00529 $ JCOUNT = 1 00530 230 CONTINUE 00531 END IF 00532 * 00533 * Set the right hand side alternately zero and small. 00534 * 00535 IF( UPPER ) THEN 00536 B( 1 ) = ZERO 00537 DO 240 I = N, 2, -2 00538 B( I ) = ZERO 00539 B( I-1 ) = SMLNUM 00540 240 CONTINUE 00541 ELSE 00542 B( N ) = ZERO 00543 DO 250 I = 1, N - 1, 2 00544 B( I ) = ZERO 00545 B( I+1 ) = SMLNUM 00546 250 CONTINUE 00547 END IF 00548 * 00549 ELSE IF( IMAT.EQ.15 ) THEN 00550 * 00551 * Type 15: Make the diagonal elements small to cause gradual 00552 * overflow when dividing by T(j,j). To control the amount of 00553 * scaling needed, the matrix is bidiagonal. 00554 * 00555 TEXP = ONE / MAX( ONE, DBLE( N-1 ) ) 00556 TSCAL = SMLNUM**TEXP 00557 CALL DLARNV( 2, ISEED, N, B ) 00558 IF( UPPER ) THEN 00559 DO 270 J = 1, N 00560 DO 260 I = 1, J - 2 00561 A( I, J ) = 0.D0 00562 260 CONTINUE 00563 IF( J.GT.1 ) 00564 $ A( J-1, J ) = -ONE 00565 A( J, J ) = TSCAL 00566 270 CONTINUE 00567 B( N ) = ONE 00568 ELSE 00569 DO 290 J = 1, N 00570 DO 280 I = J + 2, N 00571 A( I, J ) = 0.D0 00572 280 CONTINUE 00573 IF( J.LT.N ) 00574 $ A( J+1, J ) = -ONE 00575 A( J, J ) = TSCAL 00576 290 CONTINUE 00577 B( 1 ) = ONE 00578 END IF 00579 * 00580 ELSE IF( IMAT.EQ.16 ) THEN 00581 * 00582 * Type 16: One zero diagonal element. 00583 * 00584 IY = N / 2 + 1 00585 IF( UPPER ) THEN 00586 DO 300 J = 1, N 00587 CALL DLARNV( 2, ISEED, J, A( 1, J ) ) 00588 IF( J.NE.IY ) THEN 00589 A( J, J ) = SIGN( TWO, A( J, J ) ) 00590 ELSE 00591 A( J, J ) = ZERO 00592 END IF 00593 300 CONTINUE 00594 ELSE 00595 DO 310 J = 1, N 00596 CALL DLARNV( 2, ISEED, N-J+1, A( J, J ) ) 00597 IF( J.NE.IY ) THEN 00598 A( J, J ) = SIGN( TWO, A( J, J ) ) 00599 ELSE 00600 A( J, J ) = ZERO 00601 END IF 00602 310 CONTINUE 00603 END IF 00604 CALL DLARNV( 2, ISEED, N, B ) 00605 CALL DSCAL( N, TWO, B, 1 ) 00606 * 00607 ELSE IF( IMAT.EQ.17 ) THEN 00608 * 00609 * Type 17: Make the offdiagonal elements large to cause overflow 00610 * when adding a column of T. In the non-transposed case, the 00611 * matrix is constructed to cause overflow when adding a column in 00612 * every other step. 00613 * 00614 TSCAL = UNFL / ULP 00615 TSCAL = ( ONE-ULP ) / TSCAL 00616 DO 330 J = 1, N 00617 DO 320 I = 1, N 00618 A( I, J ) = 0.D0 00619 320 CONTINUE 00620 330 CONTINUE 00621 TEXP = ONE 00622 IF( UPPER ) THEN 00623 DO 340 J = N, 2, -2 00624 A( 1, J ) = -TSCAL / DBLE( N+1 ) 00625 A( J, J ) = ONE 00626 B( J ) = TEXP*( ONE-ULP ) 00627 A( 1, J-1 ) = -( TSCAL / DBLE( N+1 ) ) / DBLE( N+2 ) 00628 A( J-1, J-1 ) = ONE 00629 B( J-1 ) = TEXP*DBLE( N*N+N-1 ) 00630 TEXP = TEXP*2.D0 00631 340 CONTINUE 00632 B( 1 ) = ( DBLE( N+1 ) / DBLE( N+2 ) )*TSCAL 00633 ELSE 00634 DO 350 J = 1, N - 1, 2 00635 A( N, J ) = -TSCAL / DBLE( N+1 ) 00636 A( J, J ) = ONE 00637 B( J ) = TEXP*( ONE-ULP ) 00638 A( N, J+1 ) = -( TSCAL / DBLE( N+1 ) ) / DBLE( N+2 ) 00639 A( J+1, J+1 ) = ONE 00640 B( J+1 ) = TEXP*DBLE( N*N+N-1 ) 00641 TEXP = TEXP*2.D0 00642 350 CONTINUE 00643 B( N ) = ( DBLE( N+1 ) / DBLE( N+2 ) )*TSCAL 00644 END IF 00645 * 00646 ELSE IF( IMAT.EQ.18 ) THEN 00647 * 00648 * Type 18: Generate a unit triangular matrix with elements 00649 * between -1 and 1, and make the right hand side large so that it 00650 * requires scaling. 00651 * 00652 IF( UPPER ) THEN 00653 DO 360 J = 1, N 00654 CALL DLARNV( 2, ISEED, J-1, A( 1, J ) ) 00655 A( J, J ) = ZERO 00656 360 CONTINUE 00657 ELSE 00658 DO 370 J = 1, N 00659 IF( J.LT.N ) 00660 $ CALL DLARNV( 2, ISEED, N-J, A( J+1, J ) ) 00661 A( J, J ) = ZERO 00662 370 CONTINUE 00663 END IF 00664 * 00665 * Set the right hand side so that the largest value is BIGNUM. 00666 * 00667 CALL DLARNV( 2, ISEED, N, B ) 00668 IY = IDAMAX( N, B, 1 ) 00669 BNORM = ABS( B( IY ) ) 00670 BSCAL = BIGNUM / MAX( ONE, BNORM ) 00671 CALL DSCAL( N, BSCAL, B, 1 ) 00672 * 00673 ELSE IF( IMAT.EQ.19 ) THEN 00674 * 00675 * Type 19: Generate a triangular matrix with elements between 00676 * BIGNUM/(n-1) and BIGNUM so that at least one of the column 00677 * norms will exceed BIGNUM. 00678 * 1/3/91: DLATRS no longer can handle this case 00679 * 00680 TLEFT = BIGNUM / MAX( ONE, DBLE( N-1 ) ) 00681 TSCAL = BIGNUM*( DBLE( N-1 ) / MAX( ONE, DBLE( N ) ) ) 00682 IF( UPPER ) THEN 00683 DO 390 J = 1, N 00684 CALL DLARNV( 2, ISEED, J, A( 1, J ) ) 00685 DO 380 I = 1, J 00686 A( I, J ) = SIGN( TLEFT, A( I, J ) ) + TSCAL*A( I, J ) 00687 380 CONTINUE 00688 390 CONTINUE 00689 ELSE 00690 DO 410 J = 1, N 00691 CALL DLARNV( 2, ISEED, N-J+1, A( J, J ) ) 00692 DO 400 I = J, N 00693 A( I, J ) = SIGN( TLEFT, A( I, J ) ) + TSCAL*A( I, J ) 00694 400 CONTINUE 00695 410 CONTINUE 00696 END IF 00697 CALL DLARNV( 2, ISEED, N, B ) 00698 CALL DSCAL( N, TWO, B, 1 ) 00699 END IF 00700 * 00701 * Flip the matrix if the transpose will be used. 00702 * 00703 IF( .NOT.LSAME( TRANS, 'N' ) ) THEN 00704 IF( UPPER ) THEN 00705 DO 420 J = 1, N / 2 00706 CALL DSWAP( N-2*J+1, A( J, J ), LDA, A( J+1, N-J+1 ), 00707 $ -1 ) 00708 420 CONTINUE 00709 ELSE 00710 DO 430 J = 1, N / 2 00711 CALL DSWAP( N-2*J+1, A( J, J ), 1, A( N-J+1, J+1 ), 00712 $ -LDA ) 00713 430 CONTINUE 00714 END IF 00715 END IF 00716 * 00717 RETURN 00718 * 00719 * End of DLATTR 00720 * 00721 END