189 SUBROUTINE ddrvls( DOTYPE, NM, MVAL, NN, NVAL, NNS, NSVAL, NNB,
190 $ NBVAL, NXVAL, THRESH, TSTERR, A, COPYA, B,
191 $ COPYB, C, S, COPYS, NOUT )
199 INTEGER NM, NN, NNB, NNS, NOUT
200 DOUBLE PRECISION THRESH
204 INTEGER MVAL( * ), NBVAL( * ), NSVAL( * ),
205 $ nval( * ), nxval( * )
206 DOUBLE PRECISION A( * ), B( * ), C( * ), COPYA( * ), COPYB( * ),
214 PARAMETER ( NTESTS = 18 )
216 parameter( smlsiz = 25 )
217 DOUBLE PRECISION ONE, TWO, ZERO
218 parameter( one = 1.0d0, two = 2.0d0, zero = 0.0d0 )
223 INTEGER CRANK, I, IM, IMB, IN, INB, INFO, INS, IRANK,
224 $ iscale, itran, itype, j, k, lda, ldb, ldwork,
225 $ lwlsy, lwork, m, mnmin, n, nb, ncols, nerrs,
226 $ nfail, nrhs, nrows, nrun, rank, mb,
227 $ mmax, nmax, nsmax, liwork,
228 $ lwork_dgels, lwork_dgelst, lwork_dgetsls,
229 $ lwork_dgelss, lwork_dgelsy, lwork_dgelsd
230 DOUBLE PRECISION EPS, NORMA, NORMB, RCOND
233 INTEGER ISEED( 4 ), ISEEDY( 4 ), IWQ( 1 )
234 DOUBLE PRECISION RESULT( NTESTS ), WQ( 1 )
237 DOUBLE PRECISION,
ALLOCATABLE :: WORK (:)
238 INTEGER,
ALLOCATABLE :: IWORK (:)
241 DOUBLE PRECISION DASUM, DLAMCH, DQRT12, DQRT14, DQRT17
242 EXTERNAL DASUM, DLAMCH, DQRT12, DQRT14, DQRT17
251 INTRINSIC dble, int, max, min, sqrt
256 INTEGER INFOT, IOUNIT
259 COMMON / infoc / infot, iounit, ok, lerr
260 COMMON / srnamc / srnamt
263 DATA iseedy / 1988, 1989, 1990, 1991 /
269 path( 1: 1 ) =
'Double precision'
275 iseed( i ) = iseedy( i )
277 eps = dlamch(
'Epsilon' )
281 rcond = sqrt( eps ) - ( sqrt( eps )-eps ) / 2
288 $
CALL derrls( path, nout )
292 IF( ( nm.EQ.0 .OR. nn.EQ.0 ) .AND. thresh.EQ.zero )
293 $
CALL alahd( nout, path )
304 IF ( mval( i ).GT.mmax )
THEN
309 IF ( nval( i ).GT.nmax )
THEN
314 IF ( nsval( i ).GT.nsmax )
THEN
321 mnmin = max( min( m, n ), 1 )
326 lwork = max( 1, ( m+n )*nrhs,
327 $ ( n+nrhs )*( m+2 ), ( m+nrhs )*( n+2 ),
328 $ max( m+mnmin, nrhs*mnmin,2*n+m ),
329 $ max( m*n+4*mnmin+max(m,n), m*n+2*mnmin+4*n ) )
341 mnmin = max(min( m, n ),1)
347 itype = ( irank-1 )*3 + iscale
348 IF( dotype( itype ) )
THEN
349 IF( irank.EQ.1 )
THEN
351 IF( itran.EQ.1 )
THEN
358 CALL dgels( trans, m, n, nrhs, a, lda,
359 $ b, ldb, wq, -1, info )
360 lwork_dgels = int( wq( 1 ) )
362 CALL dgelst( trans, m, n, nrhs, a, lda,
363 $ b, ldb, wq, -1, info )
364 lwork_dgelst = int( wq( 1 ) )
366 CALL dgetsls( trans, m, n, nrhs, a, lda,
367 $ b, ldb, wq, -1, info )
368 lwork_dgetsls = int( wq( 1 ) )
372 CALL dgelsy( m, n, nrhs, a, lda, b, ldb, iwq,
373 $ rcond, crank, wq, -1, info )
374 lwork_dgelsy = int( wq( 1 ) )
376 CALL dgelss( m, n, nrhs, a, lda, b, ldb, s,
377 $ rcond, crank, wq, -1 , info )
378 lwork_dgelss = int( wq( 1 ) )
380 CALL dgelsd( m, n, nrhs, a, lda, b, ldb, s,
381 $ rcond, crank, wq, -1, iwq, info )
382 lwork_dgelsd = int( wq( 1 ) )
384 liwork = max( liwork, n, iwq( 1 ) )
386 lwork = max( lwork, lwork_dgels, lwork_dgelst,
387 $ lwork_dgetsls, lwork_dgelsy,
388 $ lwork_dgelss, lwork_dgelsd )
398 ALLOCATE( work( lwork ) )
399 ALLOCATE( iwork( liwork ) )
407 mnmin = max(min( m, n ),1)
416 itype = ( irank-1 )*3 + iscale
417 IF( .NOT.dotype( itype ) )
422 IF( irank.EQ.1 )
THEN
426 CALL dqrt13( iscale, m, n, copya, lda, norma,
434 CALL xlaenv( 3, nxval( inb ) )
439 IF( itran.EQ.1 )
THEN
448 ldwork = max( 1, ncols )
452 IF( ncols.GT.0 )
THEN
453 CALL dlarnv( 2, iseed, ncols*nrhs,
455 CALL dscal( ncols*nrhs,
456 $ one / dble( ncols ), work,
459 CALL dgemm( trans,
'No transpose', nrows,
460 $ nrhs, ncols, one, copya, lda,
461 $ work, ldwork, zero, b, ldb )
462 CALL dlacpy(
'Full', nrows, nrhs, b, ldb,
467 IF( m.GT.0 .AND. n.GT.0 )
THEN
468 CALL dlacpy(
'Full', m, n, copya, lda,
470 CALL dlacpy(
'Full', nrows, nrhs,
471 $ copyb, ldb, b, ldb )
474 CALL dgels( trans, m, n, nrhs, a, lda, b,
475 $ ldb, work, lwork, info )
477 $
CALL alaerh( path,
'DGELS ', info, 0,
478 $ trans, m, n, nrhs, -1, nb,
479 $ itype, nfail, nerrs,
487 IF( nrows.GT.0 .AND. nrhs.GT.0 )
488 $
CALL dlacpy(
'Full', nrows, nrhs,
489 $ copyb, ldb, c, ldb )
490 CALL dqrt16( trans, m, n, nrhs, copya,
491 $ lda, b, ldb, c, ldb, work,
497 IF( ( itran.EQ.1 .AND. m.GE.n ) .OR.
498 $ ( itran.EQ.2 .AND. m.LT.n ) )
THEN
504 result( 2 ) = dqrt17( trans, 1, m, n,
505 $ nrhs, copya, lda, b, ldb,
506 $ copyb, ldb, c, work,
512 result( 2 ) = dqrt14( trans, m, n,
513 $ nrhs, copya, lda, b, ldb,
521 IF( result( k ).GE.thresh )
THEN
522 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
523 $
CALL alahd( nout, path )
524 WRITE( nout, fmt = 9999 ) trans, m,
525 $ n, nrhs, nb, itype, k,
540 IF( irank.EQ.1 )
THEN
544 CALL dqrt13( iscale, m, n, copya, lda, norma,
556 IF( itran.EQ.1 )
THEN
565 ldwork = max( 1, ncols )
569 IF( ncols.GT.0 )
THEN
570 CALL dlarnv( 2, iseed, ncols*nrhs,
572 CALL dscal( ncols*nrhs,
573 $ one / dble( ncols ), work,
576 CALL dgemm( trans,
'No transpose', nrows,
577 $ nrhs, ncols, one, copya, lda,
578 $ work, ldwork, zero, b, ldb )
579 CALL dlacpy(
'Full', nrows, nrhs, b, ldb,
584 IF( m.GT.0 .AND. n.GT.0 )
THEN
585 CALL dlacpy(
'Full', m, n, copya, lda,
587 CALL dlacpy(
'Full', nrows, nrhs,
588 $ copyb, ldb, b, ldb )
591 CALL dgelst( trans, m, n, nrhs, a, lda, b,
592 $ ldb, work, lwork, info )
594 $
CALL alaerh( path,
'DGELST', info, 0,
595 $ trans, m, n, nrhs, -1, nb,
596 $ itype, nfail, nerrs,
604 IF( nrows.GT.0 .AND. nrhs.GT.0 )
605 $
CALL dlacpy(
'Full', nrows, nrhs,
606 $ copyb, ldb, c, ldb )
607 CALL dqrt16( trans, m, n, nrhs, copya,
608 $ lda, b, ldb, c, ldb, work,
614 IF( ( itran.EQ.1 .AND. m.GE.n ) .OR.
615 $ ( itran.EQ.2 .AND. m.LT.n ) )
THEN
621 result( 4 ) = dqrt17( trans, 1, m, n,
622 $ nrhs, copya, lda, b, ldb,
623 $ copyb, ldb, c, work,
629 result( 4 ) = dqrt14( trans, m, n,
630 $ nrhs, copya, lda, b, ldb,
638 IF( result( k ).GE.thresh )
THEN
639 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
640 $
CALL alahd( nout, path )
641 WRITE( nout, fmt = 9999 ) trans, m,
642 $ n, nrhs, nb, itype, k,
657 IF( irank.EQ.1 )
THEN
661 CALL dqrt13( iscale, m, n, copya, lda, norma,
680 IF( itran.EQ.1 )
THEN
689 ldwork = max( 1, ncols )
693 IF( ncols.GT.0 )
THEN
694 CALL dlarnv( 2, iseed, ncols*nrhs,
696 CALL dscal( ncols*nrhs,
697 $ one / dble( ncols ),
700 CALL dgemm( trans,
'No transpose',
701 $ nrows, nrhs, ncols, one,
702 $ copya, lda, work, ldwork,
704 CALL dlacpy(
'Full', nrows, nrhs,
705 $ b, ldb, copyb, ldb )
709 IF( m.GT.0 .AND. n.GT.0 )
THEN
710 CALL dlacpy(
'Full', m, n,
711 $ copya, lda, a, lda )
712 CALL dlacpy(
'Full', nrows, nrhs,
713 $ copyb, ldb, b, ldb )
716 CALL dgetsls( trans, m, n, nrhs,
717 $ a, lda, b, ldb, work, lwork,
720 $
CALL alaerh( path,
'DGETSLS', info,
721 $ 0, trans, m, n, nrhs,
722 $ -1, nb, itype, nfail,
730 IF( nrows.GT.0 .AND. nrhs.GT.0 )
731 $
CALL dlacpy(
'Full', nrows, nrhs,
732 $ copyb, ldb, c, ldb )
733 CALL dqrt16( trans, m, n, nrhs,
734 $ copya, lda, b, ldb,
741 IF( ( itran.EQ.1 .AND. m.GE.n ) .OR.
742 $ ( itran.EQ.2 .AND. m.LT.n ) )
THEN
748 result( 6 ) = dqrt17( trans, 1, m,
749 $ n, nrhs, copya, lda,
750 $ b, ldb, copyb, ldb,
756 result( 6 ) = dqrt14( trans, m, n,
758 $ b, ldb, work, lwork )
765 IF( result( k ).GE.thresh )
THEN
766 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
767 $
CALL alahd( nout, path )
768 WRITE( nout, fmt = 9997 ) trans,
769 $ m, n, nrhs, mb, nb, itype,
786 CALL dqrt15( iscale, irank, m, n, nrhs, copya, lda,
787 $ copyb, ldb, copys, rank, norma, normb,
788 $ iseed, work, lwork )
799 CALL xlaenv( 3, nxval( inb ) )
814 CALL dlacpy(
'Full', m, n, copya, lda, a, lda )
815 CALL dlacpy(
'Full', m, nrhs, copyb, ldb, b,
819 CALL dgelsy( m, n, nrhs, a, lda, b, ldb, iwork,
820 $ rcond, crank, work, lwlsy, info )
822 $
CALL alaerh( path,
'DGELSY', info, 0,
' ', m,
823 $ n, nrhs, -1, nb, itype, nfail,
829 result( 7 ) = dqrt12( crank, crank, a, lda,
830 $ copys, work, lwork )
835 CALL dlacpy(
'Full', m, nrhs, copyb, ldb, work,
837 CALL dqrt16(
'No transpose', m, n, nrhs, copya,
838 $ lda, b, ldb, work, ldwork,
839 $ work( m*nrhs+1 ), result( 8 ) )
846 $ result( 9 ) = dqrt17(
'No transpose', 1, m,
847 $ n, nrhs, copya, lda, b, ldb,
848 $ copyb, ldb, c, work, lwork )
856 $ result( 10 ) = dqrt14(
'No transpose', m, n,
857 $ nrhs, copya, lda, b, ldb,
866 CALL dlacpy(
'Full', m, n, copya, lda, a, lda )
867 CALL dlacpy(
'Full', m, nrhs, copyb, ldb, b,
870 CALL dgelss( m, n, nrhs, a, lda, b, ldb, s,
871 $ rcond, crank, work, lwork, info )
873 $
CALL alaerh( path,
'DGELSS', info, 0,
' ', m,
874 $ n, nrhs, -1, nb, itype, nfail,
883 CALL daxpy( mnmin, -one, copys, 1, s, 1 )
884 result( 11 ) = dasum( mnmin, s, 1 ) /
885 $ dasum( mnmin, copys, 1 ) /
886 $ ( eps*dble( mnmin ) )
893 CALL dlacpy(
'Full', m, nrhs, copyb, ldb, work,
895 CALL dqrt16(
'No transpose', m, n, nrhs, copya,
896 $ lda, b, ldb, work, ldwork,
897 $ work( m*nrhs+1 ), result( 12 ) )
903 $ result( 13 ) = dqrt17(
'No transpose', 1, m,
904 $ n, nrhs, copya, lda, b, ldb,
905 $ copyb, ldb, c, work, lwork )
911 $ result( 14 ) = dqrt14(
'No transpose', m, n,
912 $ nrhs, copya, lda, b, ldb,
927 CALL dlacpy(
'Full', m, n, copya, lda, a, lda )
928 CALL dlacpy(
'Full', m, nrhs, copyb, ldb, b,
932 CALL dgelsd( m, n, nrhs, a, lda, b, ldb, s,
933 $ rcond, crank, work, lwork, iwork,
936 $
CALL alaerh( path,
'DGELSD', info, 0,
' ', m,
937 $ n, nrhs, -1, nb, itype, nfail,
943 CALL daxpy( mnmin, -one, copys, 1, s, 1 )
944 result( 15 ) = dasum( mnmin, s, 1 ) /
945 $ dasum( mnmin, copys, 1 ) /
946 $ ( eps*dble( mnmin ) )
953 CALL dlacpy(
'Full', m, nrhs, copyb, ldb, work,
955 CALL dqrt16(
'No transpose', m, n, nrhs, copya,
956 $ lda, b, ldb, work, ldwork,
957 $ work( m*nrhs+1 ), result( 16 ) )
963 $ result( 17 ) = dqrt17(
'No transpose', 1, m,
964 $ n, nrhs, copya, lda, b, ldb,
965 $ copyb, ldb, c, work, lwork )
971 $ result( 18 ) = dqrt14(
'No transpose', m, n,
972 $ nrhs, copya, lda, b, ldb,
979 IF( result( k ).GE.thresh )
THEN
980 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
981 $
CALL alahd( nout, path )
982 WRITE( nout, fmt = 9998 )m, n, nrhs, nb,
983 $ itype, k, result( k )
1004 CALL alasvm( path, nout, nfail, nrun, nerrs )
1006 9999
FORMAT(
' TRANS=''', a1,
''', M=', i5,
', N=', i5,
', NRHS=', i4,
1007 $
', NB=', i4,
', type', i2,
', test(', i2,
')=', g12.5 )
1008 9998
FORMAT(
' M=', i5,
', N=', i5,
', NRHS=', i4,
', NB=', i4,
1009 $
', type', i2,
', test(', i2,
')=', g12.5 )
1010 9997
FORMAT(
' TRANS=''', a1,
' M=', i5,
', N=', i5,
', NRHS=', i4,
1011 $
', MB=', i4,
', NB=', i4,
', type', i2,
1012 $
', test(', i2,
')=', g12.5 )