LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine schkgb ( logical, dimension( * )  DOTYPE,
integer  NM,
integer, dimension( * )  MVAL,
integer  NN,
integer, dimension( * )  NVAL,
integer  NNB,
integer, dimension( * )  NBVAL,
integer  NNS,
integer, dimension( * )  NSVAL,
real  THRESH,
logical  TSTERR,
real, dimension( * )  A,
integer  LA,
real, dimension( * )  AFAC,
integer  LAFAC,
real, dimension( * )  B,
real, dimension( * )  X,
real, dimension( * )  XACT,
real, dimension( * )  WORK,
real, dimension( * )  RWORK,
integer, dimension( * )  IWORK,
integer  NOUT 
)

SCHKGB

Purpose:
 SCHKGB tests SGBTRF, -TRS, -RFS, and -CON
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.
[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 contained in the vector NBVAL.
[in]NBVAL
          NBVAL is INTEGER array, dimension (NNB)
          The values of the blocksize NB.
[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 REAL array, dimension (LA)
[in]LA
          LA is INTEGER
          The length of the array A.  LA >= (KLMAX+KUMAX+1)*NMAX
          where KLMAX is the largest entry in the local array KLVAL,
                KUMAX is the largest entry in the local array KUVAL and
                NMAX is the largest entry in the input array NVAL.
[out]AFAC
          AFAC is REAL array, dimension (LAFAC)
[in]LAFAC
          LAFAC is INTEGER
          The length of the array AFAC. LAFAC >= (2*KLMAX+KUMAX+1)*NMAX
          where KLMAX is the largest entry in the local array KLVAL,
                KUMAX is the largest entry in the local array KUVAL and
                NMAX is the largest entry in the input array NVAL.
[out]B
          B is REAL array, dimension (NMAX*NSMAX)
          where NSMAX is the largest entry in NSVAL.
[out]X
          X is REAL array, dimension (NMAX*NSMAX)
[out]XACT
          XACT is REAL array, dimension (NMAX*NSMAX)
[out]WORK
          WORK is REAL array, dimension
                      (NMAX*max(3,NSMAX,NMAX))
[out]RWORK
          RWORK is REAL array, dimension
                      (max(NMAX,2*NSMAX))
[out]IWORK
          IWORK is INTEGER array, dimension (2*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 2011

Definition at line 193 of file schkgb.f.

193 *
194 * -- LAPACK test routine (version 3.4.0) --
195 * -- LAPACK is a software package provided by Univ. of Tennessee, --
196 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
197 * November 2011
198 *
199 * .. Scalar Arguments ..
200  LOGICAL tsterr
201  INTEGER la, lafac, nm, nn, nnb, nns, nout
202  REAL thresh
203 * ..
204 * .. Array Arguments ..
205  LOGICAL dotype( * )
206  INTEGER iwork( * ), mval( * ), nbval( * ), nsval( * ),
207  $ nval( * )
208  REAL a( * ), afac( * ), b( * ), rwork( * ),
209  $ work( * ), x( * ), xact( * )
210 * ..
211 *
212 * =====================================================================
213 *
214 * .. Parameters ..
215  REAL one, zero
216  parameter ( one = 1.0e+0, zero = 0.0e+0 )
217  INTEGER ntypes, ntests
218  parameter ( ntypes = 8, ntests = 7 )
219  INTEGER nbw, ntran
220  parameter ( nbw = 4, ntran = 3 )
221 * ..
222 * .. Local Scalars ..
223  LOGICAL trfcon, zerot
224  CHARACTER dist, norm, trans, TYPE, xtype
225  CHARACTER*3 path
226  INTEGER i, i1, i2, ikl, iku, im, imat, in, inb, info,
227  $ ioff, irhs, itran, izero, j, k, kl, koff, ku,
228  $ lda, ldafac, ldb, m, mode, n, nb, nerrs, nfail,
229  $ nimat, nkl, nku, nrhs, nrun
230  REAL ainvnm, anorm, anormi, anormo, cndnum, rcond,
231  $ rcondc, rcondi, rcondo
232 * ..
233 * .. Local Arrays ..
234  CHARACTER transs( ntran )
235  INTEGER iseed( 4 ), iseedy( 4 ), klval( nbw ),
236  $ kuval( nbw )
237  REAL result( ntests )
238 * ..
239 * .. External Functions ..
240  REAL sget06, slangb, slange
241  EXTERNAL sget06, slangb, slange
242 * ..
243 * .. External Subroutines ..
244  EXTERNAL alaerh, alahd, alasum, scopy, serrge, sgbcon,
247  $ xlaenv
248 * ..
249 * .. Intrinsic Functions ..
250  INTRINSIC max, min
251 * ..
252 * .. Scalars in Common ..
253  LOGICAL lerr, ok
254  CHARACTER*32 srnamt
255  INTEGER infot, nunit
256 * ..
257 * .. Common blocks ..
258  COMMON / infoc / infot, nunit, ok, lerr
259  COMMON / srnamc / srnamt
260 * ..
261 * .. Data statements ..
262  DATA iseedy / 1988, 1989, 1990, 1991 / ,
263  $ transs / 'N', 'T', 'C' /
264 * ..
265 * .. Executable Statements ..
266 *
267 * Initialize constants and the random number seed.
268 *
269  path( 1: 1 ) = 'Single precision'
270  path( 2: 3 ) = 'GB'
271  nrun = 0
272  nfail = 0
273  nerrs = 0
274  DO 10 i = 1, 4
275  iseed( i ) = iseedy( i )
276  10 CONTINUE
277 *
278 * Test the error exits
279 *
280  IF( tsterr )
281  $ CALL serrge( path, nout )
282  infot = 0
283  CALL xlaenv( 2, 2 )
284 *
285 * Initialize the first value for the lower and upper bandwidths.
286 *
287  klval( 1 ) = 0
288  kuval( 1 ) = 0
289 *
290 * Do for each value of M in MVAL
291 *
292  DO 160 im = 1, nm
293  m = mval( im )
294 *
295 * Set values to use for the lower bandwidth.
296 *
297  klval( 2 ) = m + ( m+1 ) / 4
298 *
299 * KLVAL( 2 ) = MAX( M-1, 0 )
300 *
301  klval( 3 ) = ( 3*m-1 ) / 4
302  klval( 4 ) = ( m+1 ) / 4
303 *
304 * Do for each value of N in NVAL
305 *
306  DO 150 in = 1, nn
307  n = nval( in )
308  xtype = 'N'
309 *
310 * Set values to use for the upper bandwidth.
311 *
312  kuval( 2 ) = n + ( n+1 ) / 4
313 *
314 * KUVAL( 2 ) = MAX( N-1, 0 )
315 *
316  kuval( 3 ) = ( 3*n-1 ) / 4
317  kuval( 4 ) = ( n+1 ) / 4
318 *
319 * Set limits on the number of loop iterations.
320 *
321  nkl = min( m+1, 4 )
322  IF( n.EQ.0 )
323  $ nkl = 2
324  nku = min( n+1, 4 )
325  IF( m.EQ.0 )
326  $ nku = 2
327  nimat = ntypes
328  IF( m.LE.0 .OR. n.LE.0 )
329  $ nimat = 1
330 *
331  DO 140 ikl = 1, nkl
332 *
333 * Do for KL = 0, (5*M+1)/4, (3M-1)/4, and (M+1)/4. This
334 * order makes it easier to skip redundant values for small
335 * values of M.
336 *
337  kl = klval( ikl )
338  DO 130 iku = 1, nku
339 *
340 * Do for KU = 0, (5*N+1)/4, (3N-1)/4, and (N+1)/4. This
341 * order makes it easier to skip redundant values for
342 * small values of N.
343 *
344  ku = kuval( iku )
345 *
346 * Check that A and AFAC are big enough to generate this
347 * matrix.
348 *
349  lda = kl + ku + 1
350  ldafac = 2*kl + ku + 1
351  IF( ( lda*n ).GT.la .OR. ( ldafac*n ).GT.lafac ) THEN
352  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
353  $ CALL alahd( nout, path )
354  IF( n*( kl+ku+1 ).GT.la ) THEN
355  WRITE( nout, fmt = 9999 )la, m, n, kl, ku,
356  $ n*( kl+ku+1 )
357  nerrs = nerrs + 1
358  END IF
359  IF( n*( 2*kl+ku+1 ).GT.lafac ) THEN
360  WRITE( nout, fmt = 9998 )lafac, m, n, kl, ku,
361  $ n*( 2*kl+ku+1 )
362  nerrs = nerrs + 1
363  END IF
364  GO TO 130
365  END IF
366 *
367  DO 120 imat = 1, nimat
368 *
369 * Do the tests only if DOTYPE( IMAT ) is true.
370 *
371  IF( .NOT.dotype( imat ) )
372  $ GO TO 120
373 *
374 * Skip types 2, 3, or 4 if the matrix size is too
375 * small.
376 *
377  zerot = imat.GE.2 .AND. imat.LE.4
378  IF( zerot .AND. n.LT.imat-1 )
379  $ GO TO 120
380 *
381  IF( .NOT.zerot .OR. .NOT.dotype( 1 ) ) THEN
382 *
383 * Set up parameters with SLATB4 and generate a
384 * test matrix with SLATMS.
385 *
386  CALL slatb4( path, imat, m, n, TYPE, kl, ku,
387  $ anorm, mode, cndnum, dist )
388 *
389  koff = max( 1, ku+2-n )
390  DO 20 i = 1, koff - 1
391  a( i ) = zero
392  20 CONTINUE
393  srnamt = 'SLATMS'
394  CALL slatms( m, n, dist, iseed, TYPE, rwork,
395  $ mode, cndnum, anorm, kl, ku, 'Z',
396  $ a( koff ), lda, work, info )
397 *
398 * Check the error code from SLATMS.
399 *
400  IF( info.NE.0 ) THEN
401  CALL alaerh( path, 'SLATMS', info, 0, ' ', m,
402  $ n, kl, ku, -1, imat, nfail,
403  $ nerrs, nout )
404  GO TO 120
405  END IF
406  ELSE IF( izero.GT.0 ) THEN
407 *
408 * Use the same matrix for types 3 and 4 as for
409 * type 2 by copying back the zeroed out column.
410 *
411  CALL scopy( i2-i1+1, b, 1, a( ioff+i1 ), 1 )
412  END IF
413 *
414 * For types 2, 3, and 4, zero one or more columns of
415 * the matrix to test that INFO is returned correctly.
416 *
417  izero = 0
418  IF( zerot ) THEN
419  IF( imat.EQ.2 ) THEN
420  izero = 1
421  ELSE IF( imat.EQ.3 ) THEN
422  izero = min( m, n )
423  ELSE
424  izero = min( m, n ) / 2 + 1
425  END IF
426  ioff = ( izero-1 )*lda
427  IF( imat.LT.4 ) THEN
428 *
429 * Store the column to be zeroed out in B.
430 *
431  i1 = max( 1, ku+2-izero )
432  i2 = min( kl+ku+1, ku+1+( m-izero ) )
433  CALL scopy( i2-i1+1, a( ioff+i1 ), 1, b, 1 )
434 *
435  DO 30 i = i1, i2
436  a( ioff+i ) = zero
437  30 CONTINUE
438  ELSE
439  DO 50 j = izero, n
440  DO 40 i = max( 1, ku+2-j ),
441  $ min( kl+ku+1, ku+1+( m-j ) )
442  a( ioff+i ) = zero
443  40 CONTINUE
444  ioff = ioff + lda
445  50 CONTINUE
446  END IF
447  END IF
448 *
449 * These lines, if used in place of the calls in the
450 * loop over INB, cause the code to bomb on a Sun
451 * SPARCstation.
452 *
453 * ANORMO = SLANGB( 'O', N, KL, KU, A, LDA, RWORK )
454 * ANORMI = SLANGB( 'I', N, KL, KU, A, LDA, RWORK )
455 *
456 * Do for each blocksize in NBVAL
457 *
458  DO 110 inb = 1, nnb
459  nb = nbval( inb )
460  CALL xlaenv( 1, nb )
461 *
462 * Compute the LU factorization of the band matrix.
463 *
464  IF( m.GT.0 .AND. n.GT.0 )
465  $ CALL slacpy( 'Full', kl+ku+1, n, a, lda,
466  $ afac( kl+1 ), ldafac )
467  srnamt = 'SGBTRF'
468  CALL sgbtrf( m, n, kl, ku, afac, ldafac, iwork,
469  $ info )
470 *
471 * Check error code from SGBTRF.
472 *
473  IF( info.NE.izero )
474  $ CALL alaerh( path, 'SGBTRF', info, izero,
475  $ ' ', m, n, kl, ku, nb, imat,
476  $ nfail, nerrs, nout )
477  trfcon = .false.
478 *
479 *+ TEST 1
480 * Reconstruct matrix from factors and compute
481 * residual.
482 *
483  CALL sgbt01( m, n, kl, ku, a, lda, afac, ldafac,
484  $ iwork, work, result( 1 ) )
485 *
486 * Print information about the tests so far that
487 * did not pass the threshold.
488 *
489  IF( result( 1 ).GE.thresh ) THEN
490  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
491  $ CALL alahd( nout, path )
492  WRITE( nout, fmt = 9997 )m, n, kl, ku, nb,
493  $ imat, 1, result( 1 )
494  nfail = nfail + 1
495  END IF
496  nrun = nrun + 1
497 *
498 * Skip the remaining tests if this is not the
499 * first block size or if M .ne. N.
500 *
501  IF( inb.GT.1 .OR. m.NE.n )
502  $ GO TO 110
503 *
504  anormo = slangb( 'O', n, kl, ku, a, lda, rwork )
505  anormi = slangb( 'I', n, kl, ku, a, lda, rwork )
506 *
507  IF( info.EQ.0 ) THEN
508 *
509 * Form the inverse of A so we can get a good
510 * estimate of CNDNUM = norm(A) * norm(inv(A)).
511 *
512  ldb = max( 1, n )
513  CALL slaset( 'Full', n, n, zero, one, work,
514  $ ldb )
515  srnamt = 'SGBTRS'
516  CALL sgbtrs( 'No transpose', n, kl, ku, n,
517  $ afac, ldafac, iwork, work, ldb,
518  $ info )
519 *
520 * Compute the 1-norm condition number of A.
521 *
522  ainvnm = slange( 'O', n, n, work, ldb,
523  $ rwork )
524  IF( anormo.LE.zero .OR. ainvnm.LE.zero ) THEN
525  rcondo = one
526  ELSE
527  rcondo = ( one / anormo ) / ainvnm
528  END IF
529 *
530 * Compute the infinity-norm condition number of
531 * A.
532 *
533  ainvnm = slange( 'I', n, n, work, ldb,
534  $ rwork )
535  IF( anormi.LE.zero .OR. ainvnm.LE.zero ) THEN
536  rcondi = one
537  ELSE
538  rcondi = ( one / anormi ) / ainvnm
539  END IF
540  ELSE
541 *
542 * Do only the condition estimate if INFO.NE.0.
543 *
544  trfcon = .true.
545  rcondo = zero
546  rcondi = zero
547  END IF
548 *
549 * Skip the solve tests if the matrix is singular.
550 *
551  IF( trfcon )
552  $ GO TO 90
553 *
554  DO 80 irhs = 1, nns
555  nrhs = nsval( irhs )
556  xtype = 'N'
557 *
558  DO 70 itran = 1, ntran
559  trans = transs( itran )
560  IF( itran.EQ.1 ) THEN
561  rcondc = rcondo
562  norm = 'O'
563  ELSE
564  rcondc = rcondi
565  norm = 'I'
566  END IF
567 *
568 *+ TEST 2:
569 * Solve and compute residual for A * X = B.
570 *
571  srnamt = 'SLARHS'
572  CALL slarhs( path, xtype, ' ', trans, n,
573  $ n, kl, ku, nrhs, a, lda,
574  $ xact, ldb, b, ldb, iseed,
575  $ info )
576  xtype = 'C'
577  CALL slacpy( 'Full', n, nrhs, b, ldb, x,
578  $ ldb )
579 *
580  srnamt = 'SGBTRS'
581  CALL sgbtrs( trans, n, kl, ku, nrhs, afac,
582  $ ldafac, iwork, x, ldb, info )
583 *
584 * Check error code from SGBTRS.
585 *
586  IF( info.NE.0 )
587  $ CALL alaerh( path, 'SGBTRS', info, 0,
588  $ trans, n, n, kl, ku, -1,
589  $ imat, nfail, nerrs, nout )
590 *
591  CALL slacpy( 'Full', n, nrhs, b, ldb,
592  $ work, ldb )
593  CALL sgbt02( trans, m, n, kl, ku, nrhs, a,
594  $ lda, x, ldb, work, ldb,
595  $ result( 2 ) )
596 *
597 *+ TEST 3:
598 * Check solution from generated exact
599 * solution.
600 *
601  CALL sget04( n, nrhs, x, ldb, xact, ldb,
602  $ rcondc, result( 3 ) )
603 *
604 *+ TESTS 4, 5, 6:
605 * Use iterative refinement to improve the
606 * solution.
607 *
608  srnamt = 'SGBRFS'
609  CALL sgbrfs( trans, n, kl, ku, nrhs, a,
610  $ lda, afac, ldafac, iwork, b,
611  $ ldb, x, ldb, rwork,
612  $ rwork( nrhs+1 ), work,
613  $ iwork( n+1 ), info )
614 *
615 * Check error code from SGBRFS.
616 *
617  IF( info.NE.0 )
618  $ CALL alaerh( path, 'SGBRFS', info, 0,
619  $ trans, n, n, kl, ku, nrhs,
620  $ imat, nfail, nerrs, nout )
621 *
622  CALL sget04( n, nrhs, x, ldb, xact, ldb,
623  $ rcondc, result( 4 ) )
624  CALL sgbt05( trans, n, kl, ku, nrhs, a,
625  $ lda, b, ldb, x, ldb, xact,
626  $ ldb, rwork, rwork( nrhs+1 ),
627  $ result( 5 ) )
628  DO 60 k = 2, 6
629  IF( result( k ).GE.thresh ) THEN
630  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
631  $ CALL alahd( nout, path )
632  WRITE( nout, fmt = 9996 )trans, n,
633  $ kl, ku, nrhs, imat, k,
634  $ result( k )
635  nfail = nfail + 1
636  END IF
637  60 CONTINUE
638  nrun = nrun + 5
639  70 CONTINUE
640  80 CONTINUE
641 *
642 *+ TEST 7:
643 * Get an estimate of RCOND = 1/CNDNUM.
644 *
645  90 CONTINUE
646  DO 100 itran = 1, 2
647  IF( itran.EQ.1 ) THEN
648  anorm = anormo
649  rcondc = rcondo
650  norm = 'O'
651  ELSE
652  anorm = anormi
653  rcondc = rcondi
654  norm = 'I'
655  END IF
656  srnamt = 'SGBCON'
657  CALL sgbcon( norm, n, kl, ku, afac, ldafac,
658  $ iwork, anorm, rcond, work,
659  $ iwork( n+1 ), info )
660 *
661 * Check error code from SGBCON.
662 *
663  IF( info.NE.0 )
664  $ CALL alaerh( path, 'SGBCON', info, 0,
665  $ norm, n, n, kl, ku, -1, imat,
666  $ nfail, nerrs, nout )
667 *
668  result( 7 ) = sget06( rcond, rcondc )
669 *
670 * Print information about the tests that did
671 * not pass the threshold.
672 *
673  IF( result( 7 ).GE.thresh ) THEN
674  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
675  $ CALL alahd( nout, path )
676  WRITE( nout, fmt = 9995 )norm, n, kl, ku,
677  $ imat, 7, result( 7 )
678  nfail = nfail + 1
679  END IF
680  nrun = nrun + 1
681  100 CONTINUE
682 *
683  110 CONTINUE
684  120 CONTINUE
685  130 CONTINUE
686  140 CONTINUE
687  150 CONTINUE
688  160 CONTINUE
689 *
690 * Print a summary of the results.
691 *
692  CALL alasum( path, nout, nfail, nrun, nerrs )
693 *
694  9999 FORMAT( ' *** In SCHKGB, LA=', i5, ' is too small for M=', i5,
695  $ ', N=', i5, ', KL=', i4, ', KU=', i4,
696  $ / ' ==> Increase LA to at least ', i5 )
697  9998 FORMAT( ' *** In SCHKGB, LAFAC=', i5, ' is too small for M=', i5,
698  $ ', N=', i5, ', KL=', i4, ', KU=', i4,
699  $ / ' ==> Increase LAFAC to at least ', i5 )
700  9997 FORMAT( ' M =', i5, ', N =', i5, ', KL=', i5, ', KU=', i5,
701  $ ', NB =', i4, ', type ', i1, ', test(', i1, ')=', g12.5 )
702  9996 FORMAT( ' TRANS=''', a1, ''', N=', i5, ', KL=', i5, ', KU=', i5,
703  $ ', NRHS=', i3, ', type ', i1, ', test(', i1, ')=', g12.5 )
704  9995 FORMAT( ' NORM =''', a1, ''', N=', i5, ', KL=', i5, ', KU=', i5,
705  $ ',', 10x, ' type ', i1, ', test(', i1, ')=', g12.5 )
706 *
707  RETURN
708 *
709 * End of SCHKGB
710 *
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
subroutine sgbt05(TRANS, N, KL, KU, NRHS, AB, LDAB, B, LDB, X, LDX, XACT, LDXACT, FERR, BERR, RESLTS)
SGBT05
Definition: sgbt05.f:178
subroutine slatb4(PATH, IMAT, M, N, TYPE, KL, KU, ANORM, MODE, CNDNUM, DIST)
SLATB4
Definition: slatb4.f:122
subroutine slarhs(PATH, XTYPE, UPLO, TRANS, M, N, KL, KU, NRHS, A, LDA, X, LDX, B, LDB, ISEED, INFO)
SLARHS
Definition: slarhs.f:206
subroutine serrge(PATH, NUNIT)
SERRGE
Definition: serrge.f:57
subroutine sgbt02(TRANS, M, N, KL, KU, NRHS, A, LDA, X, LDX, B, LDB, RESID)
SGBT02
Definition: sgbt02.f:141
real function sget06(RCOND, RCONDC)
SGET06
Definition: sget06.f:57
subroutine xlaenv(ISPEC, NVALUE)
XLAENV
Definition: xlaenv.f:83
subroutine slacpy(UPLO, M, N, A, LDA, B, LDB)
SLACPY copies all or part of one two-dimensional array to another.
Definition: slacpy.f:105
subroutine slatms(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, KL, KU, PACK, A, LDA, WORK, INFO)
SLATMS
Definition: slatms.f:323
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...
Definition: slaset.f:112
subroutine sgbcon(NORM, N, KL, KU, AB, LDAB, IPIV, ANORM, RCOND, WORK, IWORK, INFO)
SGBCON
Definition: sgbcon.f:148
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 ...
Definition: slange.f:116
subroutine sget04(N, NRHS, X, LDX, XACT, LDXACT, RCOND, RESID)
SGET04
Definition: sget04.f:104
subroutine sgbtrs(TRANS, N, KL, KU, NRHS, AB, LDAB, IPIV, B, LDB, INFO)
SGBTRS
Definition: sgbtrs.f:140
real function slangb(NORM, N, KL, KU, AB, LDAB, WORK)
SLANGB returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: slangb.f:126
subroutine sgbt01(M, N, KL, KU, A, LDA, AFAC, LDAFAC, IPIV, WORK, RESID)
SGBT01
Definition: sgbt01.f:128
subroutine sgbrfs(TRANS, N, KL, KU, NRHS, AB, LDAB, AFB, LDAFB, IPIV, B, LDB, X, LDX, FERR, BERR, WORK, IWORK, INFO)
SGBRFS
Definition: sgbrfs.f:207
subroutine scopy(N, SX, INCX, SY, INCY)
SCOPY
Definition: scopy.f:53
subroutine alasum(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASUM
Definition: alasum.f:75
subroutine sgbtrf(M, N, KL, KU, AB, LDAB, IPIV, INFO)
SGBTRF
Definition: sgbtrf.f:146

Here is the call graph for this function:

Here is the caller graph for this function: