LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine cdrvls ( logical, dimension( * )  DOTYPE,
integer  NM,
integer, dimension( * )  MVAL,
integer  NN,
integer, dimension( * )  NVAL,
integer  NNS,
integer, dimension( * )  NSVAL,
integer  NNB,
integer, dimension( * )  NBVAL,
integer, dimension( * )  NXVAL,
real  THRESH,
logical  TSTERR,
complex, dimension( * )  A,
complex, dimension( * )  COPYA,
complex, dimension( * )  B,
complex, dimension( * )  COPYB,
complex, dimension( * )  C,
real, dimension( * )  S,
real, dimension( * )  COPYS,
complex, dimension( * )  WORK,
real, dimension( * )  RWORK,
integer, dimension( * )  IWORK,
integer  NOUT 
)

CDRVLS

Purpose:
 CDRVLS tests the least squares driver routines CGELS, CGELSS, CGELSY
 and CGELSD.
Parameters
[in]DOTYPE
          DOTYPE is LOGICAL array, dimension (NTYPES)
          The matrix types to be used for testing.  Matrices of type j
          (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
          .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
          The matrix of type j is generated as follows:
          j=1: A = U*D*V where U and V are random unitary matrices
               and D has random entries (> 0.1) taken from a uniform
               distribution (0,1). A is full rank.
          j=2: The same of 1, but A is scaled up.
          j=3: The same of 1, but A is scaled down.
          j=4: A = U*D*V where U and V are random unitary matrices
               and D has 3*min(M,N)/4 random entries (> 0.1) taken
               from a uniform distribution (0,1) and the remaining
               entries set to 0. A is rank-deficient.
          j=5: The same of 4, but A is scaled up.
          j=6: The same of 5, but A is scaled down.
[in]NM
          NM is INTEGER
          The number of values of M contained in the vector MVAL.
[in]MVAL
          MVAL is INTEGER array, dimension (NM)
          The values of the matrix row dimension M.
[in]NN
          NN is INTEGER
          The number of values of N contained in the vector NVAL.
[in]NVAL
          NVAL is INTEGER array, dimension (NN)
          The values of the matrix column dimension N.
[in]NNB
          NNB is INTEGER
          The number of values of NB and NX contained in the
          vectors NBVAL and NXVAL.  The blocking parameters are used
          in pairs (NB,NX).
[in]NBVAL
          NBVAL is INTEGER array, dimension (NNB)
          The values of the blocksize NB.
[in]NXVAL
          NXVAL is INTEGER array, dimension (NNB)
          The values of the crossover point NX.
[in]NNS
          NNS is INTEGER
          The number of values of NRHS contained in the vector NSVAL.
[in]NSVAL
          NSVAL is INTEGER array, dimension (NNS)
          The values of the number of right hand sides NRHS.
[in]THRESH
          THRESH is REAL
          The threshold value for the test ratios.  A result is
          included in the output file if RESULT >= THRESH.  To have
          every test ratio printed, use THRESH = 0.
[in]TSTERR
          TSTERR is LOGICAL
          Flag that indicates whether error exits are to be tested.
[out]A
          A is COMPLEX array, dimension (MMAX*NMAX)
          where MMAX is the maximum value of M in MVAL and NMAX is the
          maximum value of N in NVAL.
[out]COPYA
          COPYA is COMPLEX array, dimension (MMAX*NMAX)
[out]B
          B is COMPLEX array, dimension (MMAX*NSMAX)
          where MMAX is the maximum value of M in MVAL and NSMAX is the
          maximum value of NRHS in NSVAL.
[out]COPYB
          COPYB is COMPLEX array, dimension (MMAX*NSMAX)
[out]C
          C is COMPLEX array, dimension (MMAX*NSMAX)
[out]S
          S is REAL array, dimension
                      (min(MMAX,NMAX))
[out]COPYS
          COPYS is REAL array, dimension
                      (min(MMAX,NMAX))
[out]WORK
          WORK is COMPLEX array, dimension
                      (MMAX*NMAX + 4*NMAX + MMAX).
[out]RWORK
          RWORK is REAL array, dimension (5*NMAX-1)
[out]IWORK
          IWORK is INTEGER array, dimension (15*NMAX)
[in]NOUT
          NOUT is INTEGER
          The unit number for output.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2015

Definition at line 213 of file cdrvls.f.

213 *
214 * -- LAPACK test routine (version 3.6.0) --
215 * -- LAPACK is a software package provided by Univ. of Tennessee, --
216 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
217 * November 2015
218 *
219 * .. Scalar Arguments ..
220  LOGICAL tsterr
221  INTEGER nm, nn, nnb, nns, nout
222  REAL thresh
223 * ..
224 * .. Array Arguments ..
225  LOGICAL dotype( * )
226  INTEGER iwork( * ), mval( * ), nbval( * ), nsval( * ),
227  $ nval( * ), nxval( * )
228  REAL copys( * ), rwork( * ), s( * )
229  COMPLEX a( * ), b( * ), c( * ), copya( * ), copyb( * ),
230  $ work( * )
231 * ..
232 *
233 * =====================================================================
234 *
235 * .. Parameters ..
236  INTEGER ntests
237  parameter ( ntests = 14 )
238  INTEGER smlsiz
239  parameter ( smlsiz = 25 )
240  REAL one, zero
241  parameter ( one = 1.0e+0, zero = 0.0e+0 )
242  COMPLEX cone, czero
243  parameter ( cone = ( 1.0e+0, 0.0e+0 ),
244  $ czero = ( 0.0e+0, 0.0e+0 ) )
245 * ..
246 * .. Local Scalars ..
247  CHARACTER trans
248  CHARACTER*3 path
249  INTEGER crank, i, im, in, inb, info, ins, irank,
250  $ iscale, itran, itype, j, k, lda, ldb, ldwork,
251  $ lwlsy, lwork, m, mnmin, n, nb, ncols, nerrs,
252  $ nfail, nrhs, nrows, nrun, rank
253  REAL eps, norma, normb, rcond
254 * ..
255 * .. Local Arrays ..
256  INTEGER iseed( 4 ), iseedy( 4 )
257  REAL result( ntests )
258 * ..
259 * .. External Functions ..
260  REAL cqrt12, cqrt14, cqrt17, sasum, slamch
261  EXTERNAL cqrt12, cqrt14, cqrt17, sasum, slamch
262 * ..
263 * .. External Subroutines ..
264  EXTERNAL alaerh, alahd, alasvm, cerrls, cgels, cgelsd,
267  $ xlaenv
268 * ..
269 * .. Intrinsic Functions ..
270  INTRINSIC max, min, REAL, sqrt
271 * ..
272 * .. Scalars in Common ..
273  LOGICAL lerr, ok
274  CHARACTER*32 srnamt
275  INTEGER infot, iounit
276 * ..
277 * .. Common blocks ..
278  COMMON / infoc / infot, iounit, ok, lerr
279  COMMON / srnamc / srnamt
280 * ..
281 * .. Data statements ..
282  DATA iseedy / 1988, 1989, 1990, 1991 /
283 * ..
284 * .. Executable Statements ..
285 *
286 * Initialize constants and the random number seed.
287 *
288  path( 1: 1 ) = 'Complex precision'
289  path( 2: 3 ) = 'LS'
290  nrun = 0
291  nfail = 0
292  nerrs = 0
293  DO 10 i = 1, 4
294  iseed( i ) = iseedy( i )
295  10 CONTINUE
296  eps = slamch( 'Epsilon' )
297 *
298 * Threshold for rank estimation
299 *
300  rcond = sqrt( eps ) - ( sqrt( eps )-eps ) / 2
301 *
302 * Test the error exits
303 *
304  CALL xlaenv( 9, smlsiz )
305  IF( tsterr )
306  $ CALL cerrls( path, nout )
307 *
308 * Print the header if NM = 0 or NN = 0 and THRESH = 0.
309 *
310  IF( ( nm.EQ.0 .OR. nn.EQ.0 ) .AND. thresh.EQ.zero )
311  $ CALL alahd( nout, path )
312  infot = 0
313 *
314  DO 140 im = 1, nm
315  m = mval( im )
316  lda = max( 1, m )
317 *
318  DO 130 in = 1, nn
319  n = nval( in )
320  mnmin = min( m, n )
321  ldb = max( 1, m, n )
322 *
323  DO 120 ins = 1, nns
324  nrhs = nsval( ins )
325  lwork = max( 1, ( m+nrhs )*( n+2 ), ( n+nrhs )*( m+2 ),
326  $ m*n+4*mnmin+max( m, n ), 2*n+m )
327 *
328  DO 110 irank = 1, 2
329  DO 100 iscale = 1, 3
330  itype = ( irank-1 )*3 + iscale
331  IF( .NOT.dotype( itype ) )
332  $ GO TO 100
333 *
334  IF( irank.EQ.1 ) THEN
335 *
336 * Test CGELS
337 *
338 * Generate a matrix of scaling type ISCALE
339 *
340  CALL cqrt13( iscale, m, n, copya, lda, norma,
341  $ iseed )
342  DO 40 inb = 1, nnb
343  nb = nbval( inb )
344  CALL xlaenv( 1, nb )
345  CALL xlaenv( 3, nxval( inb ) )
346 *
347  DO 30 itran = 1, 2
348  IF( itran.EQ.1 ) THEN
349  trans = 'N'
350  nrows = m
351  ncols = n
352  ELSE
353  trans = 'C'
354  nrows = n
355  ncols = m
356  END IF
357  ldwork = max( 1, ncols )
358 *
359 * Set up a consistent rhs
360 *
361  IF( ncols.GT.0 ) THEN
362  CALL clarnv( 2, iseed, ncols*nrhs,
363  $ work )
364  CALL csscal( ncols*nrhs,
365  $ one / REAL( NCOLS ), work,
366  $ 1 )
367  END IF
368  CALL cgemm( trans, 'No transpose', nrows,
369  $ nrhs, ncols, cone, copya, lda,
370  $ work, ldwork, czero, b, ldb )
371  CALL clacpy( 'Full', nrows, nrhs, b, ldb,
372  $ copyb, ldb )
373 *
374 * Solve LS or overdetermined system
375 *
376  IF( m.GT.0 .AND. n.GT.0 ) THEN
377  CALL clacpy( 'Full', m, n, copya, lda,
378  $ a, lda )
379  CALL clacpy( 'Full', nrows, nrhs,
380  $ copyb, ldb, b, ldb )
381  END IF
382  srnamt = 'CGELS '
383  CALL cgels( trans, m, n, nrhs, a, lda, b,
384  $ ldb, work, lwork, info )
385 *
386  IF( info.NE.0 )
387  $ CALL alaerh( path, 'CGELS ', info, 0,
388  $ trans, m, n, nrhs, -1, nb,
389  $ itype, nfail, nerrs,
390  $ nout )
391 *
392 * Check correctness of results
393 *
394  ldwork = max( 1, nrows )
395  IF( nrows.GT.0 .AND. nrhs.GT.0 )
396  $ CALL clacpy( 'Full', nrows, nrhs,
397  $ copyb, ldb, c, ldb )
398  CALL cqrt16( trans, m, n, nrhs, copya,
399  $ lda, b, ldb, c, ldb, rwork,
400  $ result( 1 ) )
401 *
402  IF( ( itran.EQ.1 .AND. m.GE.n ) .OR.
403  $ ( itran.EQ.2 .AND. m.LT.n ) ) THEN
404 *
405 * Solving LS system
406 *
407  result( 2 ) = cqrt17( trans, 1, m, n,
408  $ nrhs, copya, lda, b, ldb,
409  $ copyb, ldb, c, work,
410  $ lwork )
411  ELSE
412 *
413 * Solving overdetermined system
414 *
415  result( 2 ) = cqrt14( trans, m, n,
416  $ nrhs, copya, lda, b, ldb,
417  $ work, lwork )
418  END IF
419 *
420 * Print information about the tests that
421 * did not pass the threshold.
422 *
423  DO 20 k = 1, 2
424  IF( result( k ).GE.thresh ) THEN
425  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
426  $ CALL alahd( nout, path )
427  WRITE( nout, fmt = 9999 )trans, m,
428  $ n, nrhs, nb, itype, k,
429  $ result( k )
430  nfail = nfail + 1
431  END IF
432  20 CONTINUE
433  nrun = nrun + 2
434  30 CONTINUE
435  40 CONTINUE
436  END IF
437 *
438 * Generate a matrix of scaling type ISCALE and rank
439 * type IRANK.
440 *
441  CALL cqrt15( iscale, irank, m, n, nrhs, copya, lda,
442  $ copyb, ldb, copys, rank, norma, normb,
443  $ iseed, work, lwork )
444 *
445 * workspace used: MAX(M+MIN(M,N),NRHS*MIN(M,N),2*N+M)
446 *
447  ldwork = max( 1, m )
448 *
449 * Loop for testing different block sizes.
450 *
451  DO 90 inb = 1, nnb
452  nb = nbval( inb )
453  CALL xlaenv( 1, nb )
454  CALL xlaenv( 3, nxval( inb ) )
455 *
456 * Test CGELSY
457 *
458 * CGELSY: Compute the minimum-norm solution
459 * X to min( norm( A * X - B ) )
460 * using the rank-revealing orthogonal
461 * factorization.
462 *
463  CALL clacpy( 'Full', m, n, copya, lda, a, lda )
464  CALL clacpy( 'Full', m, nrhs, copyb, ldb, b,
465  $ ldb )
466 *
467 * Initialize vector IWORK.
468 *
469  DO 70 j = 1, n
470  iwork( j ) = 0
471  70 CONTINUE
472 *
473 * Set LWLSY to the adequate value.
474 *
475  lwlsy = mnmin + max( 2*mnmin, nb*( n+1 ),
476  $ mnmin+nb*nrhs )
477  lwlsy = max( 1, lwlsy )
478 *
479  srnamt = 'CGELSY'
480  CALL cgelsy( m, n, nrhs, a, lda, b, ldb, iwork,
481  $ rcond, crank, work, lwlsy, rwork,
482  $ info )
483  IF( info.NE.0 )
484  $ CALL alaerh( path, 'CGELSY', info, 0, ' ', m,
485  $ n, nrhs, -1, nb, itype, nfail,
486  $ nerrs, nout )
487 *
488 * workspace used: 2*MNMIN+NB*NB+NB*MAX(N,NRHS)
489 *
490 * Test 3: Compute relative error in svd
491 * workspace: M*N + 4*MIN(M,N) + MAX(M,N)
492 *
493  result( 3 ) = cqrt12( crank, crank, a, lda,
494  $ copys, work, lwork, rwork )
495 *
496 * Test 4: Compute error in solution
497 * workspace: M*NRHS + M
498 *
499  CALL clacpy( 'Full', m, nrhs, copyb, ldb, work,
500  $ ldwork )
501  CALL cqrt16( 'No transpose', m, n, nrhs, copya,
502  $ lda, b, ldb, work, ldwork, rwork,
503  $ result( 4 ) )
504 *
505 * Test 5: Check norm of r'*A
506 * workspace: NRHS*(M+N)
507 *
508  result( 5 ) = zero
509  IF( m.GT.crank )
510  $ result( 5 ) = cqrt17( 'No transpose', 1, m,
511  $ n, nrhs, copya, lda, b, ldb,
512  $ copyb, ldb, c, work, lwork )
513 *
514 * Test 6: Check if x is in the rowspace of A
515 * workspace: (M+NRHS)*(N+2)
516 *
517  result( 6 ) = zero
518 *
519  IF( n.GT.crank )
520  $ result( 6 ) = cqrt14( 'No transpose', m, n,
521  $ nrhs, copya, lda, b, ldb,
522  $ work, lwork )
523 *
524 * Test CGELSS
525 *
526 * CGELSS: Compute the minimum-norm solution
527 * X to min( norm( A * X - B ) )
528 * using the SVD.
529 *
530  CALL clacpy( 'Full', m, n, copya, lda, a, lda )
531  CALL clacpy( 'Full', m, nrhs, copyb, ldb, b,
532  $ ldb )
533  srnamt = 'CGELSS'
534  CALL cgelss( m, n, nrhs, a, lda, b, ldb, s,
535  $ rcond, crank, work, lwork, rwork,
536  $ info )
537 *
538  IF( info.NE.0 )
539  $ CALL alaerh( path, 'CGELSS', info, 0, ' ', m,
540  $ n, nrhs, -1, nb, itype, nfail,
541  $ nerrs, nout )
542 *
543 * workspace used: 3*min(m,n) +
544 * max(2*min(m,n),nrhs,max(m,n))
545 *
546 * Test 7: Compute relative error in svd
547 *
548  IF( rank.GT.0 ) THEN
549  CALL saxpy( mnmin, -one, copys, 1, s, 1 )
550  result( 7 ) = sasum( mnmin, s, 1 ) /
551  $ sasum( mnmin, copys, 1 ) /
552  $ ( eps*REAL( MNMIN ) )
553  ELSE
554  result( 7 ) = zero
555  END IF
556 *
557 * Test 8: Compute error in solution
558 *
559  CALL clacpy( 'Full', m, nrhs, copyb, ldb, work,
560  $ ldwork )
561  CALL cqrt16( 'No transpose', m, n, nrhs, copya,
562  $ lda, b, ldb, work, ldwork, rwork,
563  $ result( 8 ) )
564 *
565 * Test 9: Check norm of r'*A
566 *
567  result( 9 ) = zero
568  IF( m.GT.crank )
569  $ result( 9 ) = cqrt17( 'No transpose', 1, m,
570  $ n, nrhs, copya, lda, b, ldb,
571  $ copyb, ldb, c, work, lwork )
572 *
573 * Test 10: Check if x is in the rowspace of A
574 *
575  result( 10 ) = zero
576  IF( n.GT.crank )
577  $ result( 10 ) = cqrt14( 'No transpose', m, n,
578  $ nrhs, copya, lda, b, ldb,
579  $ work, lwork )
580 *
581 * Test CGELSD
582 *
583 * CGELSD: Compute the minimum-norm solution X
584 * to min( norm( A * X - B ) ) using a
585 * divide and conquer SVD.
586 *
587  CALL xlaenv( 9, 25 )
588 *
589  CALL clacpy( 'Full', m, n, copya, lda, a, lda )
590  CALL clacpy( 'Full', m, nrhs, copyb, ldb, b,
591  $ ldb )
592 *
593  srnamt = 'CGELSD'
594  CALL cgelsd( m, n, nrhs, a, lda, b, ldb, s,
595  $ rcond, crank, work, lwork, rwork,
596  $ iwork, info )
597  IF( info.NE.0 )
598  $ CALL alaerh( path, 'CGELSD', info, 0, ' ', m,
599  $ n, nrhs, -1, nb, itype, nfail,
600  $ nerrs, nout )
601 *
602 * Test 11: Compute relative error in svd
603 *
604  IF( rank.GT.0 ) THEN
605  CALL saxpy( mnmin, -one, copys, 1, s, 1 )
606  result( 11 ) = sasum( mnmin, s, 1 ) /
607  $ sasum( mnmin, copys, 1 ) /
608  $ ( eps*REAL( MNMIN ) )
609  ELSE
610  result( 11 ) = zero
611  END IF
612 *
613 * Test 12: Compute error in solution
614 *
615  CALL clacpy( 'Full', m, nrhs, copyb, ldb, work,
616  $ ldwork )
617  CALL cqrt16( 'No transpose', m, n, nrhs, copya,
618  $ lda, b, ldb, work, ldwork, rwork,
619  $ result( 12 ) )
620 *
621 * Test 13: Check norm of r'*A
622 *
623  result( 13 ) = zero
624  IF( m.GT.crank )
625  $ result( 13 ) = cqrt17( 'No transpose', 1, m,
626  $ n, nrhs, copya, lda, b, ldb,
627  $ copyb, ldb, c, work, lwork )
628 *
629 * Test 14: Check if x is in the rowspace of A
630 *
631  result( 14 ) = zero
632  IF( n.GT.crank )
633  $ result( 14 ) = cqrt14( 'No transpose', m, n,
634  $ nrhs, copya, lda, b, ldb,
635  $ work, lwork )
636 *
637 * Print information about the tests that did not
638 * pass the threshold.
639 *
640  DO 80 k = 3, ntests
641  IF( result( k ).GE.thresh ) THEN
642  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
643  $ CALL alahd( nout, path )
644  WRITE( nout, fmt = 9998 )m, n, nrhs, nb,
645  $ itype, k, result( k )
646  nfail = nfail + 1
647  END IF
648  80 CONTINUE
649  nrun = nrun + 12
650 *
651  90 CONTINUE
652  100 CONTINUE
653  110 CONTINUE
654  120 CONTINUE
655  130 CONTINUE
656  140 CONTINUE
657 *
658 * Print a summary of the results.
659 *
660  CALL alasvm( path, nout, nfail, nrun, nerrs )
661 *
662  9999 FORMAT( ' TRANS=''', a1, ''', M=', i5, ', N=', i5, ', NRHS=', i4,
663  $ ', NB=', i4, ', type', i2, ', test(', i2, ')=', g12.5 )
664  9998 FORMAT( ' M=', i5, ', N=', i5, ', NRHS=', i4, ', NB=', i4,
665  $ ', type', i2, ', test(', i2, ')=', g12.5 )
666  RETURN
667 *
668 * End of CDRVLS
669 *
subroutine cqrt16(TRANS, M, N, NRHS, A, LDA, X, LDX, B, LDB, RWORK, RESID)
CQRT16
Definition: cqrt16.f:135
subroutine alasvm(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASVM
Definition: alasvm.f:75
subroutine alahd(IOUNIT, PATH)
ALAHD
Definition: alahd.f:95
subroutine alaerh(PATH, SUBNAM, INFO, INFOE, OPTS, M, N, KL, KU, N5, IMAT, NFAIL, NERRS, NOUT)
ALAERH
Definition: alaerh.f:149
real function cqrt14(TRANS, M, N, NRHS, A, LDA, X, LDX, WORK, LWORK)
CQRT14
Definition: cqrt14.f:118
subroutine cqrt15(SCALE, RKSEL, M, N, NRHS, A, LDA, B, LDB, S, RANK, NORMA, NORMB, ISEED, WORK, LWORK)
CQRT15
Definition: cqrt15.f:151
subroutine cqrt13(SCALE, M, N, A, LDA, NORMA, ISEED)
CQRT13
Definition: cqrt13.f:93
subroutine xlaenv(ISPEC, NVALUE)
XLAENV
Definition: xlaenv.f:83
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 ...
Definition: cgelsd.f:227
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
Definition: cgelss.f:180
subroutine cgels(TRANS, M, N, NRHS, A, LDA, B, LDB, WORK, LWORK, INFO)
CGELS solves overdetermined or underdetermined systems for GE matrices
Definition: cgels.f:184
real function cqrt17(TRANS, IRESID, M, N, NRHS, A, LDA, X, LDX, B, LDB, C, WORK, LWORK)
CQRT17
Definition: cqrt17.f:152
subroutine clarnv(IDIST, ISEED, N, X)
CLARNV returns a vector of random numbers from a uniform or normal distribution.
Definition: clarnv.f:101
real function sasum(N, SX, INCX)
SASUM
Definition: sasum.f:54
subroutine saxpy(N, SA, SX, INCX, SY, INCY)
SAXPY
Definition: saxpy.f:54
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:105
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
Definition: cgelsy.f:212
real function cqrt12(M, N, A, LDA, S, WORK, LWORK, RWORK)
CQRT12
Definition: cqrt12.f:99
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
subroutine cerrls(PATH, NUNIT)
CERRLS
Definition: cerrls.f:57
subroutine cgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
CGEMM
Definition: cgemm.f:189
subroutine csscal(N, SA, CX, INCX)
CSSCAL
Definition: csscal.f:54

Here is the call graph for this function:

Here is the caller graph for this function: