164
165
166
167
168
169
170 LOGICAL TSTERR
171 INTEGER NMAX, NN, NOUT, NRHS
172 DOUBLE PRECISION THRESH
173
174
175 LOGICAL DOTYPE( * )
176 INTEGER IWORK( * ), NVAL( * )
177 DOUBLE PRECISION A( * ), AFAC( * ), ASAV( * ), B( * ),
178 $ BSAV( * ), RWORK( * ), S( * ), WORK( * ),
179 $ X( * ), XACT( * )
180
181
182
183
184
185 DOUBLE PRECISION ONE, ZERO
186 parameter( one = 1.0d+0, zero = 0.0d+0 )
187 INTEGER NTYPES
188 parameter( ntypes = 11 )
189 INTEGER NTESTS
190 parameter( ntests = 7 )
191 INTEGER NTRAN
192 parameter( ntran = 3 )
193
194
195 LOGICAL EQUIL, NOFACT, PREFAC, TRFCON, ZEROT
196 CHARACTER DIST, EQUED, FACT, TRANS, TYPE, XTYPE
197 CHARACTER*3 PATH
198 INTEGER I, IEQUED, IFACT, IMAT, IN, INFO, IOFF, ITRAN,
199 $ IZERO, K, K1, KL, KU, LDA, LWORK, MODE, N, NB,
200 $ NBMIN, NERRS, NFACT, NFAIL, NIMAT, NRUN, NT
201 DOUBLE PRECISION AINVNM, AMAX, ANORM, ANORMI, ANORMO, CNDNUM,
202 $ COLCND, RCOND, RCONDC, RCONDI, RCONDO, ROLDC,
203 $ ROLDI, ROLDO, ROWCND, RPVGRW
204
205
206 CHARACTER EQUEDS( 4 ), FACTS( 3 ), TRANSS( NTRAN )
207 INTEGER ISEED( 4 ), ISEEDY( 4 )
208 DOUBLE PRECISION RESULT( NTESTS )
209
210
211 LOGICAL LSAME
212 DOUBLE PRECISION DGET06, DLAMCH, DLANGE, DLANTR
214
215
220
221
222 INTRINSIC abs, max
223
224
225 LOGICAL LERR, OK
226 CHARACTER*32 SRNAMT
227 INTEGER INFOT, NUNIT
228
229
230 COMMON / infoc / infot, nunit, ok, lerr
231 COMMON / srnamc / srnamt
232
233
234 DATA iseedy / 1988, 1989, 1990, 1991 /
235 DATA transs / 'N', 'T', 'C' /
236 DATA facts / 'F', 'N', 'E' /
237 DATA equeds / 'N', 'R', 'C', 'B' /
238
239
240
241
242
243 path( 1: 1 ) = 'Double precision'
244 path( 2: 3 ) = 'GE'
245 nrun = 0
246 nfail = 0
247 nerrs = 0
248 DO 10 i = 1, 4
249 iseed( i ) = iseedy( i )
250 10 CONTINUE
251
252
253
254 IF( tsterr )
255 $
CALL derrvx( path, nout )
256 infot = 0
257
258
259
260 nb = 1
261 nbmin = 2
264
265
266
267 DO 90 in = 1, nn
268 n = nval( in )
269 lda = max( n, 1 )
270 xtype = 'N'
271 nimat = ntypes
272 IF( n.LE.0 )
273 $ nimat = 1
274
275 DO 80 imat = 1, nimat
276
277
278
279 IF( .NOT.dotype( imat ) )
280 $ GO TO 80
281
282
283
284 zerot = imat.GE.5 .AND. imat.LE.7
285 IF( zerot .AND. n.LT.imat-4 )
286 $ GO TO 80
287
288
289
290
291 CALL dlatb4( path, imat, n, n,
TYPE, KL, KU, ANORM, MODE,
292 $ CNDNUM, DIST )
293 rcondc = one / cndnum
294
295 srnamt = 'DLATMS'
296 CALL dlatms( n, n, dist, iseed,
TYPE, RWORK, MODE, CNDNUM,
297 $ ANORM, KL, KU, 'No packing', A, LDA, WORK,
298 $ INFO )
299
300
301
302 IF( info.NE.0 ) THEN
303 CALL alaerh( path,
'DLATMS', info, 0,
' ', n, n, -1, -1,
304 $ -1, imat, nfail, nerrs, nout )
305 GO TO 80
306 END IF
307
308
309
310
311 IF( zerot ) THEN
312 IF( imat.EQ.5 ) THEN
313 izero = 1
314 ELSE IF( imat.EQ.6 ) THEN
315 izero = n
316 ELSE
317 izero = n / 2 + 1
318 END IF
319 ioff = ( izero-1 )*lda
320 IF( imat.LT.7 ) THEN
321 DO 20 i = 1, n
322 a( ioff+i ) = zero
323 20 CONTINUE
324 ELSE
325 CALL dlaset(
'Full', n, n-izero+1, zero, zero,
326 $ a( ioff+1 ), lda )
327 END IF
328 ELSE
329 izero = 0
330 END IF
331
332
333
334 CALL dlacpy(
'Full', n, n, a, lda, asav, lda )
335
336 DO 70 iequed = 1, 4
337 equed = equeds( iequed )
338 IF( iequed.EQ.1 ) THEN
339 nfact = 3
340 ELSE
341 nfact = 1
342 END IF
343
344 DO 60 ifact = 1, nfact
345 fact = facts( ifact )
346 prefac =
lsame( fact,
'F' )
347 nofact =
lsame( fact,
'N' )
348 equil =
lsame( fact,
'E' )
349
350 IF( zerot ) THEN
351 IF( prefac )
352 $ GO TO 60
353 rcondo = zero
354 rcondi = zero
355
356 ELSE IF( .NOT.nofact ) THEN
357
358
359
360
361
362
363 CALL dlacpy(
'Full', n, n, asav, lda, afac, lda )
364 IF( equil .OR. iequed.GT.1 ) THEN
365
366
367
368
369 CALL dgeequ( n, n, afac, lda, s, s( n+1 ),
370 $ rowcnd, colcnd, amax, info )
371 IF( info.EQ.0 .AND. n.GT.0 ) THEN
372 IF(
lsame( equed,
'R' ) )
THEN
373 rowcnd = zero
374 colcnd = one
375 ELSE IF(
lsame( equed,
'C' ) )
THEN
376 rowcnd = one
377 colcnd = zero
378 ELSE IF(
lsame( equed,
'B' ) )
THEN
379 rowcnd = zero
380 colcnd = zero
381 END IF
382
383
384
385 CALL dlaqge( n, n, afac, lda, s, s( n+1 ),
386 $ rowcnd, colcnd, amax, equed )
387 END IF
388 END IF
389
390
391
392
393 IF( equil ) THEN
394 roldo = rcondo
395 roldi = rcondi
396 END IF
397
398
399
400 anormo =
dlange(
'1', n, n, afac, lda, rwork )
401 anormi =
dlange(
'I', n, n, afac, lda, rwork )
402
403
404
405 srnamt = 'DGETRF'
406 CALL dgetrf( n, n, afac, lda, iwork, info )
407
408
409
410 CALL dlacpy(
'Full', n, n, afac, lda, a, lda )
411 lwork = nmax*max( 3, nrhs )
412 srnamt = 'DGETRI'
413 CALL dgetri( n, a, lda, iwork, work, lwork, info )
414
415
416
417 ainvnm =
dlange(
'1', n, n, a, lda, rwork )
418 IF( anormo.LE.zero .OR. ainvnm.LE.zero ) THEN
419 rcondo = one
420 ELSE
421 rcondo = ( one / anormo ) / ainvnm
422 END IF
423
424
425
426 ainvnm =
dlange(
'I', n, n, a, lda, rwork )
427 IF( anormi.LE.zero .OR. ainvnm.LE.zero ) THEN
428 rcondi = one
429 ELSE
430 rcondi = ( one / anormi ) / ainvnm
431 END IF
432 END IF
433
434 DO 50 itran = 1, ntran
435
436
437
438 trans = transs( itran )
439 IF( itran.EQ.1 ) THEN
440 rcondc = rcondo
441 ELSE
442 rcondc = rcondi
443 END IF
444
445
446
447 CALL dlacpy(
'Full', n, n, asav, lda, a, lda )
448
449
450
451 srnamt = 'DLARHS'
452 CALL dlarhs( path, xtype,
'Full', trans, n, n, kl,
453 $ ku, nrhs, a, lda, xact, lda, b, lda,
454 $ iseed, info )
455 xtype = 'C'
456 CALL dlacpy(
'Full', n, nrhs, b, lda, bsav, lda )
457
458 IF( nofact .AND. itran.EQ.1 ) THEN
459
460
461
462
463
464
465 CALL dlacpy(
'Full', n, n, a, lda, afac, lda )
466 CALL dlacpy(
'Full', n, nrhs, b, lda, x, lda )
467
468 srnamt = 'DGESV '
469 CALL dgesv( n, nrhs, afac, lda, iwork, x, lda,
470 $ info )
471
472
473
474 IF( info.NE.izero )
475 $
CALL alaerh( path,
'DGESV ', info, izero,
476 $ ' ', n, n, -1, -1, nrhs, imat,
477 $ nfail, nerrs, nout )
478
479
480
481
482 CALL dget01( n, n, a, lda, afac, lda, iwork,
483 $ rwork, result( 1 ) )
484 nt = 1
485 IF( izero.EQ.0 ) THEN
486
487
488
489 CALL dlacpy(
'Full', n, nrhs, b, lda, work,
490 $ lda )
491 CALL dget02(
'No transpose', n, n, nrhs, a,
492 $ lda, x, lda, work, lda, rwork,
493 $ result( 2 ) )
494
495
496
497 CALL dget04( n, nrhs, x, lda, xact, lda,
498 $ rcondc, result( 3 ) )
499 nt = 3
500 END IF
501
502
503
504
505 DO 30 k = 1, nt
506 IF( result( k ).GE.thresh ) THEN
507 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
508 $
CALL aladhd( nout, path )
509 WRITE( nout, fmt = 9999 )'DGESV ', n,
510 $ imat, k, result( k )
511 nfail = nfail + 1
512 END IF
513 30 CONTINUE
514 nrun = nrun + nt
515 END IF
516
517
518
519 IF( .NOT.prefac )
520 $
CALL dlaset(
'Full', n, n, zero, zero, afac,
521 $ lda )
522 CALL dlaset(
'Full', n, nrhs, zero, zero, x, lda )
523 IF( iequed.GT.1 .AND. n.GT.0 ) THEN
524
525
526
527
528 CALL dlaqge( n, n, a, lda, s, s( n+1 ), rowcnd,
529 $ colcnd, amax, equed )
530 END IF
531
532
533
534
535 srnamt = 'DGESVX'
536 CALL dgesvx( fact, trans, n, nrhs, a, lda, afac,
537 $ lda, iwork, equed, s, s( n+1 ), b,
538 $ lda, x, lda, rcond, rwork,
539 $ rwork( nrhs+1 ), work, iwork( n+1 ),
540 $ info )
541
542
543
544 IF( info.NE.izero )
545 $
CALL alaerh( path,
'DGESVX', info, izero,
546 $ fact // trans, n, n, -1, -1, nrhs,
547 $ imat, nfail, nerrs, nout )
548
549
550
551
552 IF( info.NE.0 .AND. info.LE.n) THEN
553 rpvgrw =
dlantr(
'M',
'U',
'N', info, info,
554 $ afac, lda, work )
555 IF( rpvgrw.EQ.zero ) THEN
556 rpvgrw = one
557 ELSE
558 rpvgrw =
dlange(
'M', n, info, a, lda,
559 $ work ) / rpvgrw
560 END IF
561 ELSE
562 rpvgrw =
dlantr(
'M',
'U',
'N', n, n, afac, lda,
563 $ work )
564 IF( rpvgrw.EQ.zero ) THEN
565 rpvgrw = one
566 ELSE
567 rpvgrw =
dlange(
'M', n, n, a, lda, work ) /
568 $ rpvgrw
569 END IF
570 END IF
571 result( 7 ) = abs( rpvgrw-work( 1 ) ) /
572 $ max( work( 1 ), rpvgrw ) /
574
575 IF( .NOT.prefac ) THEN
576
577
578
579
580 CALL dget01( n, n, a, lda, afac, lda, iwork,
581 $ rwork( 2*nrhs+1 ), result( 1 ) )
582 k1 = 1
583 ELSE
584 k1 = 2
585 END IF
586
587 IF( info.EQ.0 ) THEN
588 trfcon = .false.
589
590
591
592 CALL dlacpy(
'Full', n, nrhs, bsav, lda, work,
593 $ lda )
594 CALL dget02( trans, n, n, nrhs, asav, lda, x,
595 $ lda, work, lda, rwork( 2*nrhs+1 ),
596 $ result( 2 ) )
597
598
599
600 IF( nofact .OR. ( prefac .AND.
lsame( equed,
601 $ 'N' ) ) ) THEN
602 CALL dget04( n, nrhs, x, lda, xact, lda,
603 $ rcondc, result( 3 ) )
604 ELSE
605 IF( itran.EQ.1 ) THEN
606 roldc = roldo
607 ELSE
608 roldc = roldi
609 END IF
610 CALL dget04( n, nrhs, x, lda, xact, lda,
611 $ roldc, result( 3 ) )
612 END IF
613
614
615
616
617 CALL dget07( trans, n, nrhs, asav, lda, b, lda,
618 $ x, lda, xact, lda, rwork, .true.,
619 $ rwork( nrhs+1 ), result( 4 ) )
620 ELSE
621 trfcon = .true.
622 END IF
623
624
625
626
627 result( 6 ) =
dget06( rcond, rcondc )
628
629
630
631
632 IF( .NOT.trfcon ) THEN
633 DO 40 k = k1, ntests
634 IF( result( k ).GE.thresh ) THEN
635 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
636 $
CALL aladhd( nout, path )
637 IF( prefac ) THEN
638 WRITE( nout, fmt = 9997 )'DGESVX',
639 $ fact, trans, n, equed, imat, k,
640 $ result( k )
641 ELSE
642 WRITE( nout, fmt = 9998 )'DGESVX',
643 $ fact, trans, n, imat, k, result( k )
644 END IF
645 nfail = nfail + 1
646 END IF
647 40 CONTINUE
648 nrun = nrun + ntests - k1 + 1
649 ELSE
650 IF( result( 1 ).GE.thresh .AND. .NOT.prefac )
651 $ THEN
652 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
653 $
CALL aladhd( nout, path )
654 IF( prefac ) THEN
655 WRITE( nout, fmt = 9997 )'DGESVX', fact,
656 $ trans, n, equed, imat, 1, result( 1 )
657 ELSE
658 WRITE( nout, fmt = 9998 )'DGESVX', fact,
659 $ trans, n, imat, 1, result( 1 )
660 END IF
661 nfail = nfail + 1
662 nrun = nrun + 1
663 END IF
664 IF( result( 6 ).GE.thresh ) THEN
665 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
666 $
CALL aladhd( nout, path )
667 IF( prefac ) THEN
668 WRITE( nout, fmt = 9997 )'DGESVX', fact,
669 $ trans, n, equed, imat, 6, result( 6 )
670 ELSE
671 WRITE( nout, fmt = 9998 )'DGESVX', fact,
672 $ trans, n, imat, 6, result( 6 )
673 END IF
674 nfail = nfail + 1
675 nrun = nrun + 1
676 END IF
677 IF( result( 7 ).GE.thresh ) THEN
678 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
679 $
CALL aladhd( nout, path )
680 IF( prefac ) THEN
681 WRITE( nout, fmt = 9997 )'DGESVX', fact,
682 $ trans, n, equed, imat, 7, result( 7 )
683 ELSE
684 WRITE( nout, fmt = 9998 )'DGESVX', fact,
685 $ trans, n, imat, 7, result( 7 )
686 END IF
687 nfail = nfail + 1
688 nrun = nrun + 1
689 END IF
690
691 END IF
692
693 50 CONTINUE
694 60 CONTINUE
695 70 CONTINUE
696 80 CONTINUE
697 90 CONTINUE
698
699
700
701 CALL alasvm( path, nout, nfail, nrun, nerrs )
702
703 9999 FORMAT( 1x, a, ', N =', i5, ', type ', i2, ', test(', i2, ') =',
704 $ g12.5 )
705 9998 FORMAT( 1x, a, ', FACT=''', a1, ''', TRANS=''', a1, ''', N=', i5,
706 $ ', type ', i2, ', test(', i1, ')=', g12.5 )
707 9997 FORMAT( 1x, a, ', FACT=''', a1, ''', TRANS=''', a1, ''', N=', i5,
708 $ ', EQUED=''', a1, ''', type ', i2, ', test(', i1, ')=',
709 $ g12.5 )
710 RETURN
711
712
713
double precision function dlamch(CMACH)
DLAMCH
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
subroutine dlaset(UPLO, M, N, ALPHA, BETA, A, LDA)
DLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
logical function lsame(CA, CB)
LSAME
subroutine alasvm(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASVM
subroutine xlaenv(ISPEC, NVALUE)
XLAENV
subroutine aladhd(IOUNIT, PATH)
ALADHD
subroutine alaerh(PATH, SUBNAM, INFO, INFOE, OPTS, M, N, KL, KU, N5, IMAT, NFAIL, NERRS, NOUT)
ALAERH
subroutine dlarhs(PATH, XTYPE, UPLO, TRANS, M, N, KL, KU, NRHS, A, LDA, X, LDX, B, LDB, ISEED, INFO)
DLARHS
subroutine dget02(TRANS, M, N, NRHS, A, LDA, X, LDX, B, LDB, RWORK, RESID)
DGET02
subroutine dget07(TRANS, N, NRHS, A, LDA, B, LDB, X, LDX, XACT, LDXACT, FERR, CHKFERR, BERR, RESLTS)
DGET07
subroutine dget04(N, NRHS, X, LDX, XACT, LDXACT, RCOND, RESID)
DGET04
subroutine dget01(M, N, A, LDA, AFAC, LDAFAC, IPIV, RWORK, RESID)
DGET01
subroutine dlatb4(PATH, IMAT, M, N, TYPE, KL, KU, ANORM, MODE, CNDNUM, DIST)
DLATB4
subroutine derrvx(PATH, NUNIT)
DERRVX
double precision function dget06(RCOND, RCONDC)
DGET06
subroutine dlatms(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, KL, KU, PACK, A, LDA, WORK, INFO)
DLATMS
subroutine dlaqge(M, N, A, LDA, R, C, ROWCND, COLCND, AMAX, EQUED)
DLAQGE scales a general rectangular matrix, using row and column scaling factors computed by sgeequ.
double precision function dlange(NORM, M, N, A, LDA, WORK)
DLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
subroutine dgetrf(M, N, A, LDA, IPIV, INFO)
DGETRF
subroutine dgeequ(M, N, A, LDA, R, C, ROWCND, COLCND, AMAX, INFO)
DGEEQU
subroutine dgetri(N, A, LDA, IPIV, WORK, LWORK, INFO)
DGETRI
subroutine dgesv(N, NRHS, A, LDA, IPIV, B, LDB, INFO)
DGESV computes the solution to system of linear equations A * X = B for GE matrices
subroutine dgesvx(FACT, TRANS, N, NRHS, A, LDA, AF, LDAF, IPIV, EQUED, R, C, B, LDB, X, LDX, RCOND, FERR, BERR, WORK, IWORK, INFO)
DGESVX computes the solution to system of linear equations A * X = B for GE matrices
double precision function dlantr(NORM, UPLO, DIAG, M, N, A, LDA, WORK)
DLANTR returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...