184      IMPLICIT NONE
  185
  186
  187
  188
  189
  190
  191      INTEGER            NM, NN, NNB, NNS, NOUT
  192      REAL               THRESH
  193
  194
  195      LOGICAL            DOTYPE( * )
  196      INTEGER            IWORK( * ), NBVAL( * ), MVAL( * ), NVAL( * ),
  197     $                   NSVAL( * ), NXVAL( * )
  198      REAL               A( * ), COPYA( * ), B( * ), COPYB( * ),
  199     $                   S( * ), TAU( * ), WORK( * )
  200
  201
  202
  203
  204
  205      INTEGER            NTYPES
  206      parameter( ntypes = 19 )
  207      INTEGER            NTESTS
  208      parameter( ntests = 5 )
  209      REAL               ONE, ZERO, BIGNUM
  210      parameter( one = 1.0e+0, zero = 0.0e+0,
  211     $                     bignum = 1.0e+38 )
  212
  213
  214      CHARACTER          DIST, TYPE
  215      CHARACTER*3        PATH
  216      INTEGER            I, IHIGH, ILOW, IM, IMAT, IN, INC_ZERO,
  217     $                   INB, IND_OFFSET_GEN,
  218     $                   IND_IN, IND_OUT, INS, INFO,
  219     $                   ISTEP, J, J_INC, J_FIRST_NZ, JB_ZERO,
  220     $                   KFACT, KL, KMAX, KU, LDA, LW, LWORK,
  221     $                   LWORK_MQR, M, MINMN, MINMNB_GEN, MODE, N,
  222     $                   NB, NB_ZERO, NERRS, NFAIL, NB_GEN, NRHS,
  223     $                   NRUN, NX, T
  224      REAL               ANORM, CNDNUM, EPS, ABSTOL, RELTOL,
  225     $                   DTEMP, MAXC2NRMK, RELMAXC2NRMK
  226
  227
  228      INTEGER            ISEED( 4 ), ISEEDY( 4 )
  229      REAL               RESULT( NTESTS ), RDUMMY( 1 )
  230
  231
  232      REAL               SLAMCH, SQPT01, SQRT11, SQRT12, SLANGE
  234
  235
  239
  240
  241      INTRINSIC          abs, max, min, mod, real
  242
  243
  244      LOGICAL            LERR, OK
  245      CHARACTER*32       SRNAMT
  246      INTEGER            INFOT, IOUNIT
  247
  248
  249      COMMON             / infoc / infot, iounit, ok, lerr
  250      COMMON             / srnamc / srnamt
  251
  252
  253      DATA               iseedy / 1988, 1989, 1990, 1991 /
  254
  255
  256
  257
  258
  259      path( 1: 1 ) = 'Single precision'
  260      path( 2: 3 ) = 'QK'
  261      nrun = 0
  262      nfail = 0
  263      nerrs = 0
  264      DO i = 1, 4
  265         iseed( i ) = iseedy( i )
  266      END DO
  268      infot = 0
  269
  270      DO im = 1, nm
  271
  272
  273
  274         m = mval( im )
  275         lda = max( 1, m )
  276
  277         DO in = 1, nn
  278
  279
  280
  281            n = nval( in )
  282            minmn = min( m, n )
  283            lwork = max( 1, m*max( m, n )+4*minmn+max( m, n ),
  284     $                   m*n + 2*minmn + 4*n )
  285
  286            DO ins = 1, nns
  287               nrhs = nsval( ins )
  288
  289
  290
  291
  292
  293
  294
  295                  CALL slatb4( path, 14, m, nrhs, 
TYPE, KL, KU, ANORM,
 
  296     $                         MODE, CNDNUM, DIST )
  297
  298                  srnamt = 'SLATMS'
  299                  CALL slatms( m, nrhs, dist, iseed, 
TYPE, S, MODE,
 
  300     $                         CNDNUM, ANORM, KL, KU, 'No packing',
  301     $                         COPYB, LDA, WORK, INFO )
  302 
  303 
  304
  305
  306
  307                  IF( info.NE.0 ) THEN
  308                     CALL alaerh( path, 
'SLATMS', info, 0, 
' ', m,
 
  309     $                            nrhs, -1, -1, -1, 6, nfail, nerrs,
  310     $                            nout )
  311                     cycle
  312                  END IF
  313
  314               DO imat = 1, ntypes
  315
  316
  317
  318               IF( .NOT.dotype( imat ) )
  319     $            cycle
  320
  321
  322
  323
  324
  325
  326
  327
  328
  329
  330
  331
  332
  333
  334
  335
  336
  337
  338
  339
  340
  341
  342
  343
  344
  345
  346
  347
  348               IF( imat.EQ.1 ) THEN
  349
  350
  351
  352                  CALL slaset( 
'Full', m, n, zero, zero, copya, lda )
 
  353                  DO i = 1, minmn
  354                     s( i ) = zero
  355                  END DO
  356
  357               ELSE IF( (imat.GE.2 .AND. imat.LE.4 )
  358     $                  .OR. (imat.GE.14 .AND. imat.LE.19 ) ) THEN
  359
  360
  361
  362
  363
  364
  365                  CALL slatb4( path, imat, m, n, 
TYPE, KL, KU, ANORM,
 
  366     $                         MODE, CNDNUM, DIST )
  367
  368                  srnamt = 'SLATMS'
  369                  CALL slatms( m, n, dist, iseed, 
TYPE, S, MODE,
 
  370     $                         CNDNUM, ANORM, KL, KU, 'No packing',
  371     $                         COPYA, LDA, WORK, INFO )
  372
  373
  374
  375                  IF( info.NE.0 ) THEN
  376                     CALL alaerh( path, 
'SLATMS', info, 0, 
' ', m, n,
 
  377     $                            -1, -1, -1, imat, nfail, nerrs,
  378     $                            nout )
  379                     cycle
  380                  END IF
  381
  382                  CALL slaord( 
'Decreasing', minmn, s, 1 )
 
  383
  384               ELSE IF( minmn.GE.2
  385     $                  .AND. imat.GE.5 .AND. imat.LE.13 ) THEN
  386
  387
  388
  389
  390
  391
  392
  393
  394
  395
  396
  397
  398
  399                  IF( imat.EQ.5 ) THEN
  400
  401
  402
  403                     jb_zero = 1
  404                     nb_zero = 1
  405                     nb_gen = n - nb_zero
  406
  407                  ELSE IF( imat.EQ.6 ) THEN
  408
  409
  410
  411                     jb_zero = minmn
  412                     nb_zero = 1
  413                     nb_gen = n - nb_zero
  414
  415                  ELSE IF( imat.EQ.7 ) THEN
  416
  417
  418
  419                     jb_zero = n
  420                     nb_zero = 1
  421                     nb_gen = n - nb_zero
  422
  423                  ELSE IF( imat.EQ.8 ) THEN
  424
  425
  426
  427                     jb_zero = minmn / 2 + 1
  428                     nb_zero = 1
  429                     nb_gen = n - nb_zero
  430
  431                  ELSE IF( imat.EQ.9 ) THEN
  432
  433
  434
  435                     jb_zero = 1
  436                     nb_zero = minmn / 2
  437                     nb_gen = n - nb_zero
  438
  439                  ELSE IF( imat.EQ.10 ) THEN
  440
  441
  442
  443
  444                     jb_zero = minmn / 2 + 1
  445                     nb_zero = n - jb_zero + 1
  446                     nb_gen = n - nb_zero
  447
  448                  ELSE IF( imat.EQ.11 ) THEN
  449
  450
  451
  452
  453
  454                     jb_zero = minmn / 2 - (minmn / 2) / 2 + 1
  455                     nb_zero = minmn / 2
  456                     nb_gen = n - nb_zero
  457
  458                  ELSE IF( imat.EQ.12 ) THEN
  459
  460
  461
  462                     nb_gen = n / 2
  463                     nb_zero = n - nb_gen
  464                     j_inc = 2
  465                     j_first_nz = 2
  466
  467                  ELSE IF( imat.EQ.13 ) THEN
  468
  469
  470
  471                     nb_zero = n / 2
  472                     nb_gen = n - nb_zero
  473                     j_inc = 2
  474                     j_first_nz = 1
  475
  476                  END IF
  477
  478
  479
  480
  481
  482                  CALL slaset( 
'Full', m, nb_zero, zero, zero,
 
  483     $                         copya, lda )
  484
  485
  486
  487
  488
  489                  CALL slatb4( path, imat, m, nb_gen, 
TYPE, KL, KU,
 
  490     $                         ANORM, MODE, CNDNUM, DIST )
  491
  492                  srnamt = 'SLATMS'
  493
  494                  ind_offset_gen = nb_zero * lda
  495
  496                  CALL slatms( m, nb_gen, dist, iseed, 
TYPE, S, MODE,
 
  497     $                        CNDNUM, ANORM, KL, KU, 'No packing',
  498     $                        COPYA( IND_OFFSET_GEN + 1 ), LDA,
  499     $                        WORK, INFO )
  500
  501
  502
  503                  IF( info.NE.0 ) THEN
  504                     CALL alaerh( path, 
'SLATMS', info, 0, 
' ', m,
 
  505     $                            nb_gen, -1, -1, -1, imat, nfail,
  506     $                            nerrs, nout )
  507                     cycle
  508                  END IF
  509
  510
  511
  512
  513
  514                  IF( imat.EQ.6
  515     $                    .OR. imat.EQ.7
  516     $                    .OR. imat.EQ.8
  517     $                    .OR. imat.EQ.10
  518     $                    .OR. imat.EQ.11 ) THEN
  519
  520
  521
  522
  523
  524
  525                     DO j = 1, jb_zero-1, 1
  527     $                        copya( ( nb_zero+j-1)*lda+1), 1,
  528     $                        copya( (j-1)*lda + 1 ), 1 )
  529                     END DO
  530
  531                  ELSE IF( imat.EQ.12 .OR. imat.EQ.13 ) THEN
  532
  533
  534
  535
  536
  537
  538
  539
  540
  541
  542
  543                     DO j = 1, nb_gen, 1
  544                        ind_out = ( nb_zero+j-1 )*lda + 1
  545                        ind_in = ( j_inc*(j-1)+(j_first_nz-1) )*lda
  546     $                            + 1
  548     $                              copya( ind_out ), 1,
  549     $                              copya( ind_in), 1 )
  550                        END DO
  551
  552                  END IF
  553
  554
  555
  556
  557
  558
  559                  minmnb_gen = min( m, nb_gen )
  560
  561                  DO i = minmnb_gen+1, minmn
  562                     s( i ) = zero
  563                  END DO
  564
  565               ELSE
  566
  567
  568
  569                     cycle
  570               END IF
  571
  572
  573
  574               DO i = 1, n
  575                  iwork( i ) = 0
  576               END DO
  577
  578               DO inb = 1, nnb
  579
  580
  581
  582                  nb = nbval( inb )
  584                  nx = nxval( inb )
  586
  587
  588
  589
  590
  591                  DO kmax = 0, min(m,n)+1
  592
  593
  594
  595
  596
  597
  598
  599
  600
  601                  CALL slacpy( 
'All', m, n, copya, lda, a, lda )
 
  602                  CALL slacpy( 
'All', m, nrhs, copyb, lda,
 
  603     $                         a( lda*n + 1 ),  lda )
  604                  CALL slacpy( 
'All', m, nrhs, copyb, lda,
 
  605     $                         b,  lda )
  606                  CALL icopy( n, iwork( 1 ), 1, iwork( n+1 ), 1 )
 
  607                  DO i = 1, ntests
  608                    result( i ) = zero
  609                  END DO
  610
  611                  abstol = -1.0
  612                  reltol = -1.0
  613
  614
  615
  616                  lw = max( 1, max( 2*n + nb*( n+nrhs+1 ),
  617     $                              3*n + nrhs - 1 ) )
  618
  619
  620
  621                  srnamt = 'SGEQP3RK'
  622                  CALL sgeqp3rk( m, n, nrhs, kmax, abstol, reltol,
 
  623     $                           a, lda, kfact, maxc2nrmk,
  624     $                           relmaxc2nrmk, iwork( n+1 ), tau,
  625     $                           work, lw, iwork( 2*n+1 ), info )
  626
  627
  628
  629                  IF( info.LT.0 )
  630     $               
CALL alaerh( path, 
'SGEQP3RK', info, 0, 
' ',
 
  631     $                            m, n, nx, -1, nb, imat,
  632     $                            nfail, nerrs, nout )
  633
  634
  635
  636
  637
  638
  639
  640
  641
  642
  643
  644
  645
  646
  647
  648                  IF( kfact.EQ.minmn ) THEN
  649
  650                     result( 1 ) = 
sqrt12( m, n, a, lda, s, work,
 
  651     $                                     lwork )
  652
  653                     nrun = nrun + 1
  654
  655
  656
  657                  END IF
  658
  659
  660
  661
  662
  663
  664
  665                  result( 2 ) = 
sqpt01( m, n, kfact, copya, a, lda, tau,
 
  666     $                                  iwork( n+1 ), work, lwork )
  667
  668
  669
  670
  671
  672
  673
  674                  result( 3 ) = 
sqrt11( m, kfact, a, lda, tau, work,
 
  675     $                                  lwork )
  676
  677                  nrun = nrun + 2
  678
  679
  680
  681
  682
  683
  684
  685
  686
  687
  688
  689
  690
  691                  IF( min(kfact, minmn).GE.2 ) THEN
  692
  693                     DO j = 1, kfact-1, 1
  694 
  695                        dtemp = (( abs( a( (j-1)*lda+j ) ) -
  696     $                          abs( a( (j)*lda+j+1 ) ) ) /
  697     $                          abs( a(1) ) )
  698
  699                        IF( dtemp.LT.zero ) THEN
  700                           result( 4 ) = bignum
  701                        END IF
  702
  703                     END DO
  704
  705                     nrun = nrun + 1
  706
  707
  708
  709                  END IF
  710
  711
  712
  713
  714
  715
  716
  717
  718
  719
  720
  721
  722                  IF( minmn.GT.0 ) THEN
  723
  724                     lwork_mqr = max(1, nrhs)
  725                     CALL sormqr( 
'Left', 
'Transpose',
 
  726     $                            m, nrhs, kfact, a, lda, tau, b, lda,
  727     $                            work, lwork_mqr, info )
  728
  729                     DO i = 1, nrhs
  730
  731
  732
  733                        CALL saxpy( m, -one, a( ( n+i-1 )*lda+1 ), 1,
 
  734     $                              b( ( i-1 )*lda+1 ), 1 )
  735                     END DO
  736
  737                     result( 5 ) = abs(
  738     $                  
slange( 
'One-norm', m, nrhs, b, lda, rdummy ) /
 
  739     $                  ( real( m )*
slamch( 
'Epsilon' ) ) )
 
  740
  741                     nrun = nrun + 1
  742
  743
  744
  745                  END IF
  746
  747
  748
  749
  750                  DO t = 1, ntests
  751                     IF( result( t ).GE.thresh ) THEN
  752                        IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
  753     $                     
CALL alahd( nout, path )
 
  754                        WRITE( nout, fmt = 9999 ) 'SGEQP3RK', m, n,
  755     $                      nrhs, kmax, abstol, reltol,
  756     $                      nb, nx, imat, t, result( t )
  757                        nfail = nfail + 1
  758                     END IF
  759                  END DO
  760
  761
  762
  763                  END DO
  764
  765
  766
  767               END DO
  768
  769
  770
  771               END DO
  772
  773
  774
  775            END DO
  776
  777
  778
  779         END DO
  780
  781
  782
  783      END DO
  784
  785
  786
  787      CALL alasum( path, nout, nfail, nrun, nerrs )
 
  788
  789 9999 FORMAT( 1x, a, ' M =', i5, ', N =', i5, ', NRHS =', i5,
  790     $        ', KMAX =', i5, ', ABSTOL =', g12.5,
  791     $        ', RELTOL =', g12.5, ', NB =', i4, ', NX =', i4,
  792     $        ', type ', i2, ', test ', i2, ', ratio =', g12.5 )
  793
  794
  795
subroutine alasum(type, nout, nfail, nrun, nerrs)
ALASUM
subroutine xlaenv(ispec, nvalue)
XLAENV
subroutine alaerh(path, subnam, info, infoe, opts, m, n, kl, ku, n5, imat, nfail, nerrs, nout)
ALAERH
subroutine alahd(iounit, path)
ALAHD
subroutine saxpy(n, sa, sx, incx, sy, incy)
SAXPY
subroutine slacpy(uplo, m, n, a, lda, b, ldb)
SLACPY copies all or part of one two-dimensional array to another.
real function slamch(cmach)
SLAMCH
real function slange(norm, m, n, a, lda, work)
SLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
subroutine slaset(uplo, m, n, alpha, beta, a, lda)
SLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
subroutine sswap(n, sx, incx, sy, incy)
SSWAP
subroutine sormqr(side, trans, m, n, k, a, lda, tau, c, ldc, work, lwork, info)
SORMQR
subroutine icopy(n, sx, incx, sy, incy)
ICOPY
subroutine sgeqp3rk(m, n, nrhs, kmax, abstol, reltol, a, lda, k, maxc2nrmk, relmaxc2nrmk, jpiv, tau, work, lwork, iwork, info)
SGEQP3RK computes a truncated Householder QR factorization with column pivoting of a real m-by-n matr...
subroutine slaord(job, n, x, incx)
SLAORD
subroutine slatb4(path, imat, m, n, type, kl, ku, anorm, mode, cndnum, dist)
SLATB4
subroutine slatms(m, n, dist, iseed, sym, d, mode, cond, dmax, kl, ku, pack, a, lda, work, info)
SLATMS
real function sqpt01(m, n, k, a, af, lda, tau, jpvt, work, lwork)
SQPT01
real function sqrt11(m, k, a, lda, tau, work, lwork)
SQRT11
real function sqrt12(m, n, a, lda, s, work, lwork)
SQRT12