135
136
137
138
139
140
141 CHARACTER DIAG, TRANS, UPLO
142 INTEGER IMAT, INFO, KD, LDAB, N
143
144
145 INTEGER ISEED( 4 )
146 DOUBLE PRECISION AB( LDAB, * ), B( * ), WORK( * )
147
148
149
150
151
152 DOUBLE PRECISION ONE, TWO, ZERO
153 parameter( one = 1.0d+0, two = 2.0d+0, zero = 0.0d+0 )
154
155
156 LOGICAL UPPER
157 CHARACTER DIST, PACKIT, TYPE
158 CHARACTER*3 PATH
159 INTEGER I, IOFF, IY, J, JCOUNT, KL, KU, LENJ, MODE
160 DOUBLE PRECISION ANORM, BIGNUM, BNORM, BSCAL, CNDNUM, PLUS1,
161 $ PLUS2, REXP, SFAC, SMLNUM, STAR1, TEXP, TLEFT,
162 $ TNORM, TSCAL, ULP, UNFL
163
164
165 LOGICAL LSAME
166 INTEGER IDAMAX
167 DOUBLE PRECISION DLAMCH, DLARND
169
170
172
173
174 INTRINSIC abs, dble, max, min, sign, sqrt
175
176
177
178 path( 1: 1 ) = 'Double precision'
179 path( 2: 3 ) = 'TB'
180 unfl =
dlamch(
'Safe minimum' )
182 smlnum = unfl
183 bignum = ( one-ulp ) / smlnum
184 IF( ( imat.GE.6 .AND. imat.LE.9 ) .OR. imat.EQ.17 ) THEN
185 diag = 'U'
186 ELSE
187 diag = 'N'
188 END IF
189 info = 0
190
191
192
193 IF( n.LE.0 )
194 $ RETURN
195
196
197
198 upper =
lsame( uplo,
'U' )
199 IF( upper ) THEN
200 CALL dlatb4( path, imat, n, n,
TYPE, KL, KU, ANORM, MODE,
201 $ CNDNUM, DIST )
202 ku = kd
203 ioff = 1 + max( 0, kd-n+1 )
204 kl = 0
205 packit = 'Q'
206 ELSE
207 CALL dlatb4( path, -imat, n, n,
TYPE, KL, KU, ANORM, MODE,
208 $ CNDNUM, DIST )
209 kl = kd
210 ioff = 1
211 ku = 0
212 packit = 'B'
213 END IF
214
215
216
217 IF( imat.LE.5 ) THEN
218 CALL dlatms( n, n, dist, iseed,
TYPE, B, MODE, CNDNUM, ANORM,
219 $ KL, KU, PACKIT, AB( IOFF, 1 ), LDAB, WORK, INFO )
220
221
222
223
224
225
226 ELSE IF( imat.EQ.6 ) THEN
227 IF( upper ) THEN
228 DO 20 j = 1, n
229 DO 10 i = max( 1, kd+2-j ), kd
230 ab( i, j ) = zero
231 10 CONTINUE
232 ab( kd+1, j ) = j
233 20 CONTINUE
234 ELSE
235 DO 40 j = 1, n
236 ab( 1, j ) = j
237 DO 30 i = 2, min( kd+1, n-j+1 )
238 ab( i, j ) = zero
239 30 CONTINUE
240 40 CONTINUE
241 END IF
242
243
244
245
246
247
248 ELSE IF( imat.LE.9 ) THEN
249 tnorm = sqrt( cndnum )
250
251
252
253 IF( upper ) THEN
254 DO 60 j = 1, n
255 DO 50 i = max( 1, kd+2-j ), kd
256 ab( i, j ) = zero
257 50 CONTINUE
258 ab( kd+1, j ) = dble( j )
259 60 CONTINUE
260 ELSE
261 DO 80 j = 1, n
262 DO 70 i = 2, min( kd+1, n-j+1 )
263 ab( i, j ) = zero
264 70 CONTINUE
265 ab( 1, j ) = dble( j )
266 80 CONTINUE
267 END IF
268
269
270
271
272 IF( kd.EQ.1 ) THEN
273 IF( upper ) THEN
274 ab( 1, 2 ) = sign( tnorm,
dlarnd( 2, iseed ) )
275 lenj = ( n-3 ) / 2
276 CALL dlarnv( 2, iseed, lenj, work )
277 DO 90 j = 1, lenj
278 ab( 1, 2*( j+1 ) ) = tnorm*work( j )
279 90 CONTINUE
280 ELSE
281 ab( 2, 1 ) = sign( tnorm,
dlarnd( 2, iseed ) )
282 lenj = ( n-3 ) / 2
283 CALL dlarnv( 2, iseed, lenj, work )
284 DO 100 j = 1, lenj
285 ab( 2, 2*j+1 ) = tnorm*work( j )
286 100 CONTINUE
287 END IF
288 ELSE IF( kd.GT.1 ) THEN
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306 star1 = sign( tnorm,
dlarnd( 2, iseed ) )
307 sfac = sqrt( tnorm )
308 plus1 = sign( sfac,
dlarnd( 2, iseed ) )
309 DO 110 j = 1, n, 2
310 plus2 = star1 / plus1
311 work( j ) = plus1
312 work( n+j ) = star1
313 IF( j+1.LE.n ) THEN
314 work( j+1 ) = plus2
315 work( n+j+1 ) = zero
316 plus1 = star1 / plus2
317
318
319
320
322 IF( rexp.LT.zero ) THEN
323 star1 = -sfac**( one-rexp )
324 ELSE
325 star1 = sfac**( one+rexp )
326 END IF
327 END IF
328 110 CONTINUE
329
330
331
332 IF( upper ) THEN
333 CALL dcopy( n-1, work, 1, ab( kd, 2 ), ldab )
334 CALL dcopy( n-2, work( n+1 ), 1, ab( kd-1, 3 ), ldab )
335 ELSE
336 CALL dcopy( n-1, work, 1, ab( 2, 1 ), ldab )
337 CALL dcopy( n-2, work( n+1 ), 1, ab( 3, 1 ), ldab )
338 END IF
339 END IF
340
341
342
343
344
345 ELSE IF( imat.EQ.10 ) THEN
346
347
348
349
350
351 IF( upper ) THEN
352 DO 120 j = 1, n
353 lenj = min( j, kd+1 )
354 CALL dlarnv( 2, iseed, lenj, ab( kd+2-lenj, j ) )
355 ab( kd+1, j ) = sign( two, ab( kd+1, j ) )
356 120 CONTINUE
357 ELSE
358 DO 130 j = 1, n
359 lenj = min( n-j+1, kd+1 )
360 IF( lenj.GT.0 )
361 $
CALL dlarnv( 2, iseed, lenj, ab( 1, j ) )
362 ab( 1, j ) = sign( two, ab( 1, j ) )
363 130 CONTINUE
364 END IF
365
366
367
368 CALL dlarnv( 2, iseed, n, b )
370 bnorm = abs( b( iy ) )
371 bscal = bignum / max( one, bnorm )
372 CALL dscal( n, bscal, b, 1 )
373
374 ELSE IF( imat.EQ.11 ) THEN
375
376
377
378
379
380 CALL dlarnv( 2, iseed, n, b )
381 tscal = one / dble( kd+1 )
382 IF( upper ) THEN
383 DO 140 j = 1, n
384 lenj = min( j, kd+1 )
385 CALL dlarnv( 2, iseed, lenj, ab( kd+2-lenj, j ) )
386 CALL dscal( lenj-1, tscal, ab( kd+2-lenj, j ), 1 )
387 ab( kd+1, j ) = sign( one, ab( kd+1, j ) )
388 140 CONTINUE
389 ab( kd+1, n ) = smlnum*ab( kd+1, n )
390 ELSE
391 DO 150 j = 1, n
392 lenj = min( n-j+1, kd+1 )
393 CALL dlarnv( 2, iseed, lenj, ab( 1, j ) )
394 IF( lenj.GT.1 )
395 $
CALL dscal( lenj-1, tscal, ab( 2, j ), 1 )
396 ab( 1, j ) = sign( one, ab( 1, j ) )
397 150 CONTINUE
398 ab( 1, 1 ) = smlnum*ab( 1, 1 )
399 END IF
400
401 ELSE IF( imat.EQ.12 ) THEN
402
403
404
405
406
407 CALL dlarnv( 2, iseed, n, b )
408 IF( upper ) THEN
409 DO 160 j = 1, n
410 lenj = min( j, kd+1 )
411 CALL dlarnv( 2, iseed, lenj, ab( kd+2-lenj, j ) )
412 ab( kd+1, j ) = sign( one, ab( kd+1, j ) )
413 160 CONTINUE
414 ab( kd+1, n ) = smlnum*ab( kd+1, n )
415 ELSE
416 DO 170 j = 1, n
417 lenj = min( n-j+1, kd+1 )
418 CALL dlarnv( 2, iseed, lenj, ab( 1, j ) )
419 ab( 1, j ) = sign( one, ab( 1, j ) )
420 170 CONTINUE
421 ab( 1, 1 ) = smlnum*ab( 1, 1 )
422 END IF
423
424 ELSE IF( imat.EQ.13 ) THEN
425
426
427
428
429
430 IF( upper ) THEN
431 jcount = 1
432 DO 190 j = n, 1, -1
433 DO 180 i = max( 1, kd+1-( j-1 ) ), kd
434 ab( i, j ) = zero
435 180 CONTINUE
436 IF( jcount.LE.2 ) THEN
437 ab( kd+1, j ) = smlnum
438 ELSE
439 ab( kd+1, j ) = one
440 END IF
441 jcount = jcount + 1
442 IF( jcount.GT.4 )
443 $ jcount = 1
444 190 CONTINUE
445 ELSE
446 jcount = 1
447 DO 210 j = 1, n
448 DO 200 i = 2, min( n-j+1, kd+1 )
449 ab( i, j ) = zero
450 200 CONTINUE
451 IF( jcount.LE.2 ) THEN
452 ab( 1, j ) = smlnum
453 ELSE
454 ab( 1, j ) = one
455 END IF
456 jcount = jcount + 1
457 IF( jcount.GT.4 )
458 $ jcount = 1
459 210 CONTINUE
460 END IF
461
462
463
464 IF( upper ) THEN
465 b( 1 ) = zero
466 DO 220 i = n, 2, -2
467 b( i ) = zero
468 b( i-1 ) = smlnum
469 220 CONTINUE
470 ELSE
471 b( n ) = zero
472 DO 230 i = 1, n - 1, 2
473 b( i ) = zero
474 b( i+1 ) = smlnum
475 230 CONTINUE
476 END IF
477
478 ELSE IF( imat.EQ.14 ) THEN
479
480
481
482
483
484 texp = one / dble( kd+1 )
485 tscal = smlnum**texp
486 CALL dlarnv( 2, iseed, n, b )
487 IF( upper ) THEN
488 DO 250 j = 1, n
489 DO 240 i = max( 1, kd+2-j ), kd
490 ab( i, j ) = zero
491 240 CONTINUE
492 IF( j.GT.1 .AND. kd.GT.0 )
493 $ ab( kd, j ) = -one
494 ab( kd+1, j ) = tscal
495 250 CONTINUE
496 b( n ) = one
497 ELSE
498 DO 270 j = 1, n
499 DO 260 i = 3, min( n-j+1, kd+1 )
500 ab( i, j ) = zero
501 260 CONTINUE
502 IF( j.LT.n .AND. kd.GT.0 )
503 $ ab( 2, j ) = -one
504 ab( 1, j ) = tscal
505 270 CONTINUE
506 b( 1 ) = one
507 END IF
508
509 ELSE IF( imat.EQ.15 ) THEN
510
511
512
513 iy = n / 2 + 1
514 IF( upper ) THEN
515 DO 280 j = 1, n
516 lenj = min( j, kd+1 )
517 CALL dlarnv( 2, iseed, lenj, ab( kd+2-lenj, j ) )
518 IF( j.NE.iy ) THEN
519 ab( kd+1, j ) = sign( two, ab( kd+1, j ) )
520 ELSE
521 ab( kd+1, j ) = zero
522 END IF
523 280 CONTINUE
524 ELSE
525 DO 290 j = 1, n
526 lenj = min( n-j+1, kd+1 )
527 CALL dlarnv( 2, iseed, lenj, ab( 1, j ) )
528 IF( j.NE.iy ) THEN
529 ab( 1, j ) = sign( two, ab( 1, j ) )
530 ELSE
531 ab( 1, j ) = zero
532 END IF
533 290 CONTINUE
534 END IF
535 CALL dlarnv( 2, iseed, n, b )
536 CALL dscal( n, two, b, 1 )
537
538 ELSE IF( imat.EQ.16 ) THEN
539
540
541
542
543
544
545 tscal = unfl / ulp
546 tscal = ( one-ulp ) / tscal
547 DO 310 j = 1, n
548 DO 300 i = 1, kd + 1
549 ab( i, j ) = zero
550 300 CONTINUE
551 310 CONTINUE
552 texp = one
553 IF( kd.GT.0 ) THEN
554 IF( upper ) THEN
555 DO 330 j = n, 1, -kd
556 DO 320 i = j, max( 1, j-kd+1 ), -2
557 ab( 1+( j-i ), i ) = -tscal / dble( kd+2 )
558 ab( kd+1, i ) = one
559 b( i ) = texp*( one-ulp )
560 IF( i.GT.max( 1, j-kd+1 ) ) THEN
561 ab( 2+( j-i ), i-1 ) = -( tscal / dble( kd+2 ) )
562 $ / dble( kd+3 )
563 ab( kd+1, i-1 ) = one
564 b( i-1 ) = texp*dble( ( kd+1 )*( kd+1 )+kd )
565 END IF
566 texp = texp*two
567 320 CONTINUE
568 b( max( 1, j-kd+1 ) ) = ( dble( kd+2 ) /
569 $ dble( kd+3 ) )*tscal
570 330 CONTINUE
571 ELSE
572 DO 350 j = 1, n, kd
573 texp = one
574 lenj = min( kd+1, n-j+1 )
575 DO 340 i = j, min( n, j+kd-1 ), 2
576 ab( lenj-( i-j ), j ) = -tscal / dble( kd+2 )
577 ab( 1, j ) = one
578 b( j ) = texp*( one-ulp )
579 IF( i.LT.min( n, j+kd-1 ) ) THEN
580 ab( lenj-( i-j+1 ), i+1 ) = -( tscal /
581 $ dble( kd+2 ) ) / dble( kd+3 )
582 ab( 1, i+1 ) = one
583 b( i+1 ) = texp*dble( ( kd+1 )*( kd+1 )+kd )
584 END IF
585 texp = texp*two
586 340 CONTINUE
587 b( min( n, j+kd-1 ) ) = ( dble( kd+2 ) /
588 $ dble( kd+3 ) )*tscal
589 350 CONTINUE
590 END IF
591 ELSE
592 DO 360 j = 1, n
593 ab( 1, j ) = one
594 b( j ) = dble( j )
595 360 CONTINUE
596 END IF
597
598 ELSE IF( imat.EQ.17 ) THEN
599
600
601
602
603
604 IF( upper ) THEN
605 DO 370 j = 1, n
606 lenj = min( j-1, kd )
607 CALL dlarnv( 2, iseed, lenj, ab( kd+1-lenj, j ) )
608 ab( kd+1, j ) = dble( j )
609 370 CONTINUE
610 ELSE
611 DO 380 j = 1, n
612 lenj = min( n-j, kd )
613 IF( lenj.GT.0 )
614 $
CALL dlarnv( 2, iseed, lenj, ab( 2, j ) )
615 ab( 1, j ) = dble( j )
616 380 CONTINUE
617 END IF
618
619
620
621 CALL dlarnv( 2, iseed, n, b )
623 bnorm = abs( b( iy ) )
624 bscal = bignum / max( one, bnorm )
625 CALL dscal( n, bscal, b, 1 )
626
627 ELSE IF( imat.EQ.18 ) THEN
628
629
630
631
632
633 tleft = bignum / max( one, dble( kd ) )
634 tscal = bignum*( dble( kd ) / dble( kd+1 ) )
635 IF( upper ) THEN
636 DO 400 j = 1, n
637 lenj = min( j, kd+1 )
638 CALL dlarnv( 2, iseed, lenj, ab( kd+2-lenj, j ) )
639 DO 390 i = kd + 2 - lenj, kd + 1
640 ab( i, j ) = sign( tleft, ab( i, j ) ) +
641 $ tscal*ab( i, j )
642 390 CONTINUE
643 400 CONTINUE
644 ELSE
645 DO 420 j = 1, n
646 lenj = min( n-j+1, kd+1 )
647 CALL dlarnv( 2, iseed, lenj, ab( 1, j ) )
648 DO 410 i = 1, lenj
649 ab( i, j ) = sign( tleft, ab( i, j ) ) +
650 $ tscal*ab( i, j )
651 410 CONTINUE
652 420 CONTINUE
653 END IF
654 CALL dlarnv( 2, iseed, n, b )
655 CALL dscal( n, two, b, 1 )
656 END IF
657
658
659
660 IF( .NOT.
lsame( trans,
'N' ) )
THEN
661 IF( upper ) THEN
662 DO 430 j = 1, n / 2
663 lenj = min( n-2*j+1, kd+1 )
664 CALL dswap( lenj, ab( kd+1, j ), ldab-1,
665 $ ab( kd+2-lenj, n-j+1 ), -1 )
666 430 CONTINUE
667 ELSE
668 DO 440 j = 1, n / 2
669 lenj = min( n-2*j+1, kd+1 )
670 CALL dswap( lenj, ab( 1, j ), 1, ab( lenj, n-j+2-lenj ),
671 $ -ldab+1 )
672 440 CONTINUE
673 END IF
674 END IF
675
676 RETURN
677
678
679
double precision function dlarnd(idist, iseed)
DLARND
subroutine dlatb4(path, imat, m, n, type, kl, ku, anorm, mode, cndnum, dist)
DLATB4
subroutine dlatms(m, n, dist, iseed, sym, d, mode, cond, dmax, kl, ku, pack, a, lda, work, info)
DLATMS
subroutine dcopy(n, dx, incx, dy, incy)
DCOPY
integer function idamax(n, dx, incx)
IDAMAX
double precision function dlamch(cmach)
DLAMCH
subroutine dlarnv(idist, iseed, n, x)
DLARNV returns a vector of random numbers from a uniform or normal distribution.
logical function lsame(ca, cb)
LSAME
subroutine dscal(n, da, dx, incx)
DSCAL
subroutine dswap(n, dx, incx, dy, incy)
DSWAP