199 INTEGER NM, NN, NNB, NNS, NOUT
204 INTEGER MVAL( * ), NBVAL( * ), NSVAL( * ),
205 $ NVAL( * ), NXVAL( * )
206 REAL COPYS( * ), S( * )
207 COMPLEX A( * ), B( * ), C( * ), COPYA( * ), COPYB( * )
214 parameter( ntests = 16 )
216 parameter( smlsiz = 25 )
218 parameter( one = 1.0e+0, zero = 0.0e+0 )
220 parameter( cone = ( 1.0e+0, 0.0e+0 ),
221 $ czero = ( 0.0e+0, 0.0e+0 ) )
226 INTEGER CRANK, I, IM, IMB, IN, INB, INFO, INS, IRANK,
227 $ ISCALE, ITRAN, ITYPE, J, K, LDA, LDB, LDWORK,
228 $ LWLSY, LWORK, M, MNMIN, N, NB, NCOLS, NERRS,
229 $ NFAIL, NRHS, NROWS, NRUN, RANK, MB,
230 $ MMAX, NMAX, NSMAX, LIWORK, LRWORK,
231 $ LWORK_CGELS, LWORK_CGETSLS, LWORK_CGELSS,
232 $ LWORK_CGELSY, LWORK_CGELSD,
233 $ LRWORK_CGELSY, LRWORK_CGELSS, LRWORK_CGELSD
234 REAL EPS, NORMA, NORMB, RCOND
237 INTEGER ISEED( 4 ), ISEEDY( 4 ), IWQ( 1 )
238 REAL RESULT( NTESTS ), RWQ( 1 )
242 COMPLEX,
ALLOCATABLE :: WORK (:)
243 REAL,
ALLOCATABLE :: RWORK (:), WORK2 (:)
244 INTEGER,
ALLOCATABLE :: IWORK (:)
247 REAL CQRT12, CQRT14, CQRT17, SASUM, SLAMCH
257 INTRINSIC max, min, int, real, sqrt
262 INTEGER INFOT, IOUNIT
265 COMMON / infoc / infot, iounit, ok, lerr
266 COMMON / srnamc / srnamt
269 DATA iseedy / 1988, 1989, 1990, 1991 /
275 path( 1: 1 ) =
'Complex precision'
281 iseed( i ) = iseedy( i )
287 rcond = sqrt( eps ) - ( sqrt( eps )-eps ) / 2
293 $
CALL cerrls( path, nout )
297 IF( ( nm.EQ.0 .OR. nn.EQ.0 ) .AND. thresh.EQ.zero )
298 $
CALL alahd( nout, path )
307 IF ( mval( i ).GT.mmax )
THEN
312 IF ( nval( i ).GT.nmax )
THEN
317 IF ( nsval( i ).GT.nsmax )
THEN
324 mnmin = max( min( m, n ), 1 )
329 lwork = max( 1, ( m+n )*nrhs,
330 $ ( n+nrhs )*( m+2 ), ( m+nrhs )*( n+2 ),
331 $ max( m+mnmin, nrhs*mnmin,2*n+m ),
332 $ max( m*n+4*mnmin+max(m,n), m*n+2*mnmin+4*n ) )
344 mnmin = max(min( m, n ),1)
350 itype = ( irank-1 )*3 + iscale
351 IF( dotype( itype ) )
THEN
352 IF( irank.EQ.1 )
THEN
354 IF( itran.EQ.1 )
THEN
361 CALL cgels( trans, m, n, nrhs, a, lda,
362 $ b, ldb, wq, -1, info )
363 lwork_cgels = int( wq( 1 ) )
365 CALL cgetsls( trans, m, n, nrhs, a, lda,
366 $ b, ldb, wq, -1, info )
367 lwork_cgetsls = int( wq( 1 ) )
371 CALL cgelsy( m, n, nrhs, a, lda, b, ldb,
372 $ iwq, rcond, crank, wq, -1, rwq,
374 lwork_cgelsy = int( wq( 1 ) )
377 CALL cgelss( m, n, nrhs, a, lda, b, ldb, s,
378 $ rcond, crank, wq, -1, rwq, info )
379 lwork_cgelss = int( wq( 1 ) )
380 lrwork_cgelss = 5*mnmin
382 CALL cgelsd( m, n, nrhs, a, lda, b, ldb, s,
383 $ rcond, crank, wq, -1, rwq, iwq,
385 lwork_cgelsd = int( wq( 1 ) )
386 lrwork_cgelsd = int( rwq( 1 ) )
388 liwork = max( liwork, n, iwq( 1 ) )
390 lrwork = max( lrwork, lrwork_cgelsy,
391 $ lrwork_cgelss, lrwork_cgelsd )
393 lwork = max( lwork, lwork_cgels, lwork_cgetsls,
394 $ lwork_cgelsy, lwork_cgelss,
405 ALLOCATE( work( lwork ) )
406 ALLOCATE( iwork( liwork ) )
407 ALLOCATE( rwork( lrwork ) )
408 ALLOCATE( work2( 2 * lwork ) )
416 mnmin = max(min( m, n ),1)
425 itype = ( irank-1 )*3 + iscale
426 IF( .NOT.dotype( itype ) )
429 IF( irank.EQ.1 )
THEN
435 CALL cqrt13( iscale, m, n, copya, lda, norma,
440 CALL xlaenv( 3, nxval( inb ) )
443 IF( itran.EQ.1 )
THEN
452 ldwork = max( 1, ncols )
456 IF( ncols.GT.0 )
THEN
457 CALL clarnv( 2, iseed, ncols*nrhs,
460 $ one / real( ncols ), work,
463 CALL cgemm( trans,
'No transpose', nrows,
464 $ nrhs, ncols, cone, copya, lda,
465 $ work, ldwork, czero, b, ldb )
466 CALL clacpy(
'Full', nrows, nrhs, b, ldb,
471 IF( m.GT.0 .AND. n.GT.0 )
THEN
472 CALL clacpy(
'Full', m, n, copya, lda,
474 CALL clacpy(
'Full', nrows, nrhs,
475 $ copyb, ldb, b, ldb )
478 CALL cgels( trans, m, n, nrhs, a, lda, b,
479 $ ldb, work, lwork, info )
482 $
CALL alaerh( path,
'CGELS ', info, 0,
483 $ trans, m, n, nrhs, -1, nb,
484 $ itype, nfail, nerrs,
489 ldwork = max( 1, nrows )
490 IF( nrows.GT.0 .AND. nrhs.GT.0 )
491 $
CALL clacpy(
'Full', nrows, nrhs,
492 $ copyb, ldb, c, ldb )
493 CALL cqrt16( trans, m, n, nrhs, copya,
494 $ lda, b, ldb, c, ldb, rwork,
497 IF( ( itran.EQ.1 .AND. m.GE.n ) .OR.
498 $ ( itran.EQ.2 .AND. m.LT.n ) )
THEN
502 result( 2 ) =
cqrt17( trans, 1, m, n,
503 $ nrhs, copya, lda, b, ldb,
504 $ copyb, ldb, c, work,
510 result( 2 ) =
cqrt14( trans, m, n,
511 $ nrhs, copya, lda, b, ldb,
519 IF( result( k ).GE.thresh )
THEN
520 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
521 $
CALL alahd( nout, path )
522 WRITE( nout, fmt = 9999 )trans, m,
523 $ n, nrhs, nb, itype, k,
537 CALL cqrt13( iscale, m, n, copya, lda, norma,
547 IF( itran.EQ.1 )
THEN
556 ldwork = max( 1, ncols )
560 IF( ncols.GT.0 )
THEN
561 CALL clarnv( 2, iseed, ncols*nrhs,
563 CALL cscal( ncols*nrhs,
564 $ cone / real( ncols ), work,
567 CALL cgemm( trans,
'No transpose', nrows,
568 $ nrhs, ncols, cone, copya, lda,
569 $ work, ldwork, czero, b, ldb )
570 CALL clacpy(
'Full', nrows, nrhs, b, ldb,
575 IF( m.GT.0 .AND. n.GT.0 )
THEN
576 CALL clacpy(
'Full', m, n, copya, lda,
578 CALL clacpy(
'Full', nrows, nrhs,
579 $ copyb, ldb, b, ldb )
582 CALL cgetsls( trans, m, n, nrhs, a,
583 $ lda, b, ldb, work, lwork, info )
585 $
CALL alaerh( path,
'CGETSLS ', info, 0,
586 $ trans, m, n, nrhs, -1, nb,
587 $ itype, nfail, nerrs,
592 ldwork = max( 1, nrows )
593 IF( nrows.GT.0 .AND. nrhs.GT.0 )
594 $
CALL clacpy(
'Full', nrows, nrhs,
595 $ copyb, ldb, c, ldb )
596 CALL cqrt16( trans, m, n, nrhs, copya,
597 $ lda, b, ldb, c, ldb, work2,
600 IF( ( itran.EQ.1 .AND. m.GE.n ) .OR.
601 $ ( itran.EQ.2 .AND. m.LT.n ) )
THEN
605 result( 16 ) =
cqrt17( trans, 1, m, n,
606 $ nrhs, copya, lda, b, ldb,
607 $ copyb, ldb, c, work,
613 result( 16 ) =
cqrt14( trans, m, n,
614 $ nrhs, copya, lda, b, ldb,
622 IF( result( k ).GE.thresh )
THEN
623 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
624 $
CALL alahd( nout, path )
625 WRITE( nout, fmt = 9997 )trans, m,
626 $ n, nrhs, mb, nb, itype, k,
640 CALL cqrt15( iscale, irank, m, n, nrhs, copya, lda,
641 $ copyb, ldb, copys, rank, norma, normb,
642 $ iseed, work, lwork )
653 CALL xlaenv( 3, nxval( inb ) )
662 CALL clacpy(
'Full', m, n, copya, lda, a, lda )
663 CALL clacpy(
'Full', m, nrhs, copyb, ldb, b,
673 CALL cgelsy( m, n, nrhs, a, lda, b, ldb, iwork,
674 $ rcond, crank, work, lwlsy, rwork,
677 $
CALL alaerh( path,
'CGELSY', info, 0,
' ', m,
678 $ n, nrhs, -1, nb, itype, nfail,
686 result( 3 ) =
cqrt12( crank, crank, a, lda,
687 $ copys, work, lwork, rwork )
692 CALL clacpy(
'Full', m, nrhs, copyb, ldb, work,
694 CALL cqrt16(
'No transpose', m, n, nrhs, copya,
695 $ lda, b, ldb, work, ldwork, rwork,
703 $ result( 5 ) =
cqrt17(
'No transpose', 1, m,
704 $ n, nrhs, copya, lda, b, ldb,
705 $ copyb, ldb, c, work, lwork )
713 $ result( 6 ) =
cqrt14(
'No transpose', m, n,
714 $ nrhs, copya, lda, b, ldb,
723 CALL clacpy(
'Full', m, n, copya, lda, a, lda )
724 CALL clacpy(
'Full', m, nrhs, copyb, ldb, b,
727 CALL cgelss( m, n, nrhs, a, lda, b, ldb, s,
728 $ rcond, crank, work, lwork, rwork,
732 $
CALL alaerh( path,
'CGELSS', info, 0,
' ', m,
733 $ n, nrhs, -1, nb, itype, nfail,
742 CALL saxpy( mnmin, -one, copys, 1, s, 1 )
743 result( 7 ) =
sasum( mnmin, s, 1 ) /
744 $
sasum( mnmin, copys, 1 ) /
745 $ ( eps*real( mnmin ) )
752 CALL clacpy(
'Full', m, nrhs, copyb, ldb, work,
754 CALL cqrt16(
'No transpose', m, n, nrhs, copya,
755 $ lda, b, ldb, work, ldwork, rwork,
762 $ result( 9 ) =
cqrt17(
'No transpose', 1, m,
763 $ n, nrhs, copya, lda, b, ldb,
764 $ copyb, ldb, c, work, lwork )
770 $ result( 10 ) =
cqrt14(
'No transpose', m, n,
771 $ nrhs, copya, lda, b, ldb,
782 CALL clacpy(
'Full', m, n, copya, lda, a, lda )
783 CALL clacpy(
'Full', m, nrhs, copyb, ldb, b,
787 CALL cgelsd( m, n, nrhs, a, lda, b, ldb, s,
788 $ rcond, crank, work, lwork, rwork,
791 $
CALL alaerh( path,
'CGELSD', info, 0,
' ', m,
792 $ n, nrhs, -1, nb, itype, nfail,
798 CALL saxpy( mnmin, -one, copys, 1, s, 1 )
799 result( 11 ) =
sasum( mnmin, s, 1 ) /
800 $
sasum( mnmin, copys, 1 ) /
801 $ ( eps*real( mnmin ) )
808 CALL clacpy(
'Full', m, nrhs, copyb, ldb, work,
810 CALL cqrt16(
'No transpose', m, n, nrhs, copya,
811 $ lda, b, ldb, work, ldwork, rwork,
818 $ result( 13 ) =
cqrt17(
'No transpose', 1, m,
819 $ n, nrhs, copya, lda, b, ldb,
820 $ copyb, ldb, c, work, lwork )
826 $ result( 14 ) =
cqrt14(
'No transpose', m, n,
827 $ nrhs, copya, lda, b, ldb,
834 IF( result( k ).GE.thresh )
THEN
835 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
836 $
CALL alahd( nout, path )
837 WRITE( nout, fmt = 9998 )m, n, nrhs, nb,
838 $ itype, k, result( k )
853 CALL alasvm( path, nout, nfail, nrun, nerrs )
855 9999
FORMAT(
' TRANS=''', a1,
''', M=', i5,
', N=', i5,
', NRHS=', i4,
856 $
', NB=', i4,
', type', i2,
', test(', i2,
')=', g12.5 )
857 9998
FORMAT(
' M=', i5,
', N=', i5,
', NRHS=', i4,
', NB=', i4,
858 $
', type', i2,
', test(', i2,
')=', g12.5 )
859 9997
FORMAT(
' TRANS=''', a1,
' M=', i5,
', N=', i5,
', NRHS=', i4,
860 $
', MB=', i4,
', NB=', i4,
', type', i2,
861 $
', test(', i2,
')=', g12.5 )
subroutine alasvm(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASVM
subroutine xlaenv(ISPEC, NVALUE)
XLAENV
subroutine alahd(IOUNIT, PATH)
ALAHD
subroutine alaerh(PATH, SUBNAM, INFO, INFOE, OPTS, M, N, KL, KU, N5, IMAT, NFAIL, NERRS, NOUT)
ALAERH
subroutine csscal(N, SA, CX, INCX)
CSSCAL
subroutine cscal(N, CA, CX, INCX)
CSCAL
subroutine cgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
CGEMM
real function cqrt14(TRANS, M, N, NRHS, A, LDA, X, LDX, WORK, LWORK)
CQRT14
subroutine cqrt16(TRANS, M, N, NRHS, A, LDA, X, LDX, B, LDB, RWORK, RESID)
CQRT16
real function cqrt17(TRANS, IRESID, M, N, NRHS, A, LDA, X, LDX, B, LDB, C, WORK, LWORK)
CQRT17
subroutine cqrt13(SCALE, M, N, A, LDA, NORMA, ISEED)
CQRT13
subroutine cerrls(PATH, NUNIT)
CERRLS
subroutine cqrt15(SCALE, RKSEL, M, N, NRHS, A, LDA, B, LDB, S, RANK, NORMA, NORMB, ISEED, WORK, LWORK)
CQRT15
real function cqrt12(M, N, A, LDA, S, WORK, LWORK, RWORK)
CQRT12
subroutine cgels(TRANS, M, N, NRHS, A, LDA, B, LDB, WORK, LWORK, INFO)
CGELS solves overdetermined or underdetermined systems for GE matrices
subroutine cgelss(M, N, NRHS, A, LDA, B, LDB, S, RCOND, RANK, WORK, LWORK, RWORK, INFO)
CGELSS solves overdetermined or underdetermined systems for GE matrices
subroutine cgelsy(M, N, NRHS, A, LDA, B, LDB, JPVT, RCOND, RANK, WORK, LWORK, RWORK, INFO)
CGELSY solves overdetermined or underdetermined systems for GE matrices
subroutine cgelsd(M, N, NRHS, A, LDA, B, LDB, S, RCOND, RANK, WORK, LWORK, RWORK, IWORK, INFO)
CGELSD computes the minimum-norm solution to a linear least squares problem for GE matrices
subroutine cgetsls(TRANS, M, N, NRHS, A, LDA, B, LDB, WORK, LWORK, INFO)
CGETSLS
subroutine clarnv(IDIST, ISEED, N, X)
CLARNV returns a vector of random numbers from a uniform or normal distribution.
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
subroutine saxpy(N, SA, SX, INCX, SY, INCY)
SAXPY
real function sasum(N, SX, INCX)
SASUM
real function slamch(CMACH)
SLAMCH