LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ slatrs()

subroutine slatrs ( character  UPLO,
character  TRANS,
character  DIAG,
character  NORMIN,
integer  N,
real, dimension( lda, * )  A,
integer  LDA,
real, dimension( * )  X,
real  SCALE,
real, dimension( * )  CNORM,
integer  INFO 
)

SLATRS solves a triangular system of equations with the scale factor set to prevent overflow.

Download SLATRS + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 SLATRS solves one of the triangular systems

    A *x = s*b  or  A**T*x = s*b

 with scaling to prevent overflow.  Here A is an upper or lower
 triangular matrix, A**T denotes the transpose of A, x and b are
 n-element vectors, and s is a scaling factor, usually less than
 or equal to 1, chosen so that the components of x will be less than
 the overflow threshold.  If the unscaled problem will not cause
 overflow, the Level 2 BLAS routine STRSV is called.  If the matrix A
 is singular (A(j,j) = 0 for some j), then s is set to 0 and a
 non-trivial solution to A*x = 0 is returned.
Parameters
[in]UPLO
          UPLO is CHARACTER*1
          Specifies whether the matrix A is upper or lower triangular.
          = 'U':  Upper triangular
          = 'L':  Lower triangular
[in]TRANS
          TRANS is CHARACTER*1
          Specifies the operation applied to A.
          = 'N':  Solve A * x = s*b  (No transpose)
          = 'T':  Solve A**T* x = s*b  (Transpose)
          = 'C':  Solve A**T* x = s*b  (Conjugate transpose = Transpose)
[in]DIAG
          DIAG is CHARACTER*1
          Specifies whether or not the matrix A is unit triangular.
          = 'N':  Non-unit triangular
          = 'U':  Unit triangular
[in]NORMIN
          NORMIN is CHARACTER*1
          Specifies whether CNORM has been set or not.
          = 'Y':  CNORM contains the column norms on entry
          = 'N':  CNORM is not set on entry.  On exit, the norms will
                  be computed and stored in CNORM.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in]A
          A is REAL array, dimension (LDA,N)
          The triangular matrix A.  If UPLO = 'U', the leading n by n
          upper triangular part of the array A contains the upper
          triangular matrix, and the strictly lower triangular part of
          A is not referenced.  If UPLO = 'L', the leading n by n lower
          triangular part of the array A contains the lower triangular
          matrix, and the strictly upper triangular part of A is not
          referenced.  If DIAG = 'U', the diagonal elements of A are
          also not referenced and are assumed to be 1.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max (1,N).
[in,out]X
          X is REAL array, dimension (N)
          On entry, the right hand side b of the triangular system.
          On exit, X is overwritten by the solution vector x.
[out]SCALE
          SCALE is REAL
          The scaling factor s for the triangular system
             A * x = s*b  or  A**T* x = s*b.
          If SCALE = 0, the matrix A is singular or badly scaled, and
          the vector x is an exact or approximate solution to A*x = 0.
[in,out]CNORM
          CNORM is REAL array, dimension (N)

          If NORMIN = 'Y', CNORM is an input argument and CNORM(j)
          contains the norm of the off-diagonal part of the j-th column
          of A.  If TRANS = 'N', CNORM(j) must be greater than or equal
          to the infinity-norm, and if TRANS = 'T' or 'C', CNORM(j)
          must be greater than or equal to the 1-norm.

          If NORMIN = 'N', CNORM is an output argument and CNORM(j)
          returns the 1-norm of the offdiagonal part of the j-th column
          of A.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -k, the k-th argument had an illegal value
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
  A rough bound on x is computed; if that is less than overflow, STRSV
  is called, otherwise, specific code is used which checks for possible
  overflow or divide-by-zero at every operation.

  A columnwise scheme is used for solving A*x = b.  The basic algorithm
  if A is lower triangular is

       x[1:n] := b[1:n]
       for j = 1, ..., n
            x(j) := x(j) / A(j,j)
            x[j+1:n] := x[j+1:n] - x(j) * A[j+1:n,j]
       end

  Define bounds on the components of x after j iterations of the loop:
     M(j) = bound on x[1:j]
     G(j) = bound on x[j+1:n]
  Initially, let M(0) = 0 and G(0) = max{x(i), i=1,...,n}.

  Then for iteration j+1 we have
     M(j+1) <= G(j) / | A(j+1,j+1) |
     G(j+1) <= G(j) + M(j+1) * | A[j+2:n,j+1] |
            <= G(j) ( 1 + CNORM(j+1) / | A(j+1,j+1) | )

  where CNORM(j+1) is greater than or equal to the infinity-norm of
  column j+1 of A, not counting the diagonal.  Hence

     G(j) <= G(0) product ( 1 + CNORM(i) / | A(i,i) | )
                  1<=i<=j
  and

     |x(j)| <= ( G(0) / |A(j,j)| ) product ( 1 + CNORM(i) / |A(i,i)| )
                                   1<=i< j

  Since |x(j)| <= M(j), we use the Level 2 BLAS routine STRSV if the
  reciprocal of the largest M(j), j=1,..,n, is larger than
  max(underflow, 1/overflow).

  The bound on x(j) is also used to determine when a step in the
  columnwise method can be performed without fear of overflow.  If
  the computed bound is greater than a large constant, x is scaled to
  prevent overflow, but if the bound overflows, x is set to 0, x(j) to
  1, and scale to 0, and a non-trivial solution to A*x = 0 is found.

  Similarly, a row-wise scheme is used to solve A**T*x = b.  The basic
  algorithm for A upper triangular is

       for j = 1, ..., n
            x(j) := ( b(j) - A[1:j-1,j]**T * x[1:j-1] ) / A(j,j)
       end

  We simultaneously compute two bounds
       G(j) = bound on ( b(i) - A[1:i-1,i]**T * x[1:i-1] ), 1<=i<=j
       M(j) = bound on x(i), 1<=i<=j

  The initial values are G(0) = 0, M(0) = max{b(i), i=1,..,n}, and we
  add the constraint G(j) >= G(j-1) and M(j) >= M(j-1) for j >= 1.
  Then the bound on x(j) is

       M(j) <= M(j-1) * ( 1 + CNORM(j) ) / | A(j,j) |

            <= M(0) * product ( ( 1 + CNORM(i) ) / |A(i,i)| )
                      1<=i<=j

  and we can safely call STRSV if 1/M(n) and 1/G(n) are both greater
  than max(underflow, 1/overflow).

Definition at line 236 of file slatrs.f.

238 *
239 * -- LAPACK auxiliary routine --
240 * -- LAPACK is a software package provided by Univ. of Tennessee, --
241 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
242 *
243 * .. Scalar Arguments ..
244  CHARACTER DIAG, NORMIN, TRANS, UPLO
245  INTEGER INFO, LDA, N
246  REAL SCALE
247 * ..
248 * .. Array Arguments ..
249  REAL A( LDA, * ), CNORM( * ), X( * )
250 * ..
251 *
252 * =====================================================================
253 *
254 * .. Parameters ..
255  REAL ZERO, HALF, ONE
256  parameter( zero = 0.0e+0, half = 0.5e+0, one = 1.0e+0 )
257 * ..
258 * .. Local Scalars ..
259  LOGICAL NOTRAN, NOUNIT, UPPER
260  INTEGER I, IMAX, J, JFIRST, JINC, JLAST
261  REAL BIGNUM, GROW, REC, SMLNUM, SUMJ, TJJ, TJJS,
262  $ TMAX, TSCAL, USCAL, XBND, XJ, XMAX
263 * ..
264 * .. External Functions ..
265  LOGICAL LSAME
266  INTEGER ISAMAX
267  REAL SASUM, SDOT, SLAMCH
268  EXTERNAL lsame, isamax, sasum, sdot, slamch
269 * ..
270 * .. External Subroutines ..
271  EXTERNAL saxpy, sscal, strsv, xerbla
272 * ..
273 * .. Intrinsic Functions ..
274  INTRINSIC abs, max, min
275 * ..
276 * .. Executable Statements ..
277 *
278  info = 0
279  upper = lsame( uplo, 'U' )
280  notran = lsame( trans, 'N' )
281  nounit = lsame( diag, 'N' )
282 *
283 * Test the input parameters.
284 *
285  IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
286  info = -1
287  ELSE IF( .NOT.notran .AND. .NOT.lsame( trans, 'T' ) .AND. .NOT.
288  $ lsame( trans, 'C' ) ) THEN
289  info = -2
290  ELSE IF( .NOT.nounit .AND. .NOT.lsame( diag, 'U' ) ) THEN
291  info = -3
292  ELSE IF( .NOT.lsame( normin, 'Y' ) .AND. .NOT.
293  $ lsame( normin, 'N' ) ) THEN
294  info = -4
295  ELSE IF( n.LT.0 ) THEN
296  info = -5
297  ELSE IF( lda.LT.max( 1, n ) ) THEN
298  info = -7
299  END IF
300  IF( info.NE.0 ) THEN
301  CALL xerbla( 'SLATRS', -info )
302  RETURN
303  END IF
304 *
305 * Quick return if possible
306 *
307  IF( n.EQ.0 )
308  $ RETURN
309 *
310 * Determine machine dependent parameters to control overflow.
311 *
312  smlnum = slamch( 'Safe minimum' ) / slamch( 'Precision' )
313  bignum = one / smlnum
314  scale = one
315 *
316  IF( lsame( normin, 'N' ) ) THEN
317 *
318 * Compute the 1-norm of each column, not including the diagonal.
319 *
320  IF( upper ) THEN
321 *
322 * A is upper triangular.
323 *
324  DO 10 j = 1, n
325  cnorm( j ) = sasum( j-1, a( 1, j ), 1 )
326  10 CONTINUE
327  ELSE
328 *
329 * A is lower triangular.
330 *
331  DO 20 j = 1, n - 1
332  cnorm( j ) = sasum( n-j, a( j+1, j ), 1 )
333  20 CONTINUE
334  cnorm( n ) = zero
335  END IF
336  END IF
337 *
338 * Scale the column norms by TSCAL if the maximum element in CNORM is
339 * greater than BIGNUM.
340 *
341  imax = isamax( n, cnorm, 1 )
342  tmax = cnorm( imax )
343  IF( tmax.LE.bignum ) THEN
344  tscal = one
345  ELSE
346  tscal = one / ( smlnum*tmax )
347  CALL sscal( n, tscal, cnorm, 1 )
348  END IF
349 *
350 * Compute a bound on the computed solution vector to see if the
351 * Level 2 BLAS routine STRSV can be used.
352 *
353  j = isamax( n, x, 1 )
354  xmax = abs( x( j ) )
355  xbnd = xmax
356  IF( notran ) THEN
357 *
358 * Compute the growth in A * x = b.
359 *
360  IF( upper ) THEN
361  jfirst = n
362  jlast = 1
363  jinc = -1
364  ELSE
365  jfirst = 1
366  jlast = n
367  jinc = 1
368  END IF
369 *
370  IF( tscal.NE.one ) THEN
371  grow = zero
372  GO TO 50
373  END IF
374 *
375  IF( nounit ) THEN
376 *
377 * A is non-unit triangular.
378 *
379 * Compute GROW = 1/G(j) and XBND = 1/M(j).
380 * Initially, G(0) = max{x(i), i=1,...,n}.
381 *
382  grow = one / max( xbnd, smlnum )
383  xbnd = grow
384  DO 30 j = jfirst, jlast, jinc
385 *
386 * Exit the loop if the growth factor is too small.
387 *
388  IF( grow.LE.smlnum )
389  $ GO TO 50
390 *
391 * M(j) = G(j-1) / abs(A(j,j))
392 *
393  tjj = abs( a( j, j ) )
394  xbnd = min( xbnd, min( one, tjj )*grow )
395  IF( tjj+cnorm( j ).GE.smlnum ) THEN
396 *
397 * G(j) = G(j-1)*( 1 + CNORM(j) / abs(A(j,j)) )
398 *
399  grow = grow*( tjj / ( tjj+cnorm( j ) ) )
400  ELSE
401 *
402 * G(j) could overflow, set GROW to 0.
403 *
404  grow = zero
405  END IF
406  30 CONTINUE
407  grow = xbnd
408  ELSE
409 *
410 * A is unit triangular.
411 *
412 * Compute GROW = 1/G(j), where G(0) = max{x(i), i=1,...,n}.
413 *
414  grow = min( one, one / max( xbnd, smlnum ) )
415  DO 40 j = jfirst, jlast, jinc
416 *
417 * Exit the loop if the growth factor is too small.
418 *
419  IF( grow.LE.smlnum )
420  $ GO TO 50
421 *
422 * G(j) = G(j-1)*( 1 + CNORM(j) )
423 *
424  grow = grow*( one / ( one+cnorm( j ) ) )
425  40 CONTINUE
426  END IF
427  50 CONTINUE
428 *
429  ELSE
430 *
431 * Compute the growth in A**T * x = b.
432 *
433  IF( upper ) THEN
434  jfirst = 1
435  jlast = n
436  jinc = 1
437  ELSE
438  jfirst = n
439  jlast = 1
440  jinc = -1
441  END IF
442 *
443  IF( tscal.NE.one ) THEN
444  grow = zero
445  GO TO 80
446  END IF
447 *
448  IF( nounit ) THEN
449 *
450 * A is non-unit triangular.
451 *
452 * Compute GROW = 1/G(j) and XBND = 1/M(j).
453 * Initially, M(0) = max{x(i), i=1,...,n}.
454 *
455  grow = one / max( xbnd, smlnum )
456  xbnd = grow
457  DO 60 j = jfirst, jlast, jinc
458 *
459 * Exit the loop if the growth factor is too small.
460 *
461  IF( grow.LE.smlnum )
462  $ GO TO 80
463 *
464 * G(j) = max( G(j-1), M(j-1)*( 1 + CNORM(j) ) )
465 *
466  xj = one + cnorm( j )
467  grow = min( grow, xbnd / xj )
468 *
469 * M(j) = M(j-1)*( 1 + CNORM(j) ) / abs(A(j,j))
470 *
471  tjj = abs( a( j, j ) )
472  IF( xj.GT.tjj )
473  $ xbnd = xbnd*( tjj / xj )
474  60 CONTINUE
475  grow = min( grow, xbnd )
476  ELSE
477 *
478 * A is unit triangular.
479 *
480 * Compute GROW = 1/G(j), where G(0) = max{x(i), i=1,...,n}.
481 *
482  grow = min( one, one / max( xbnd, smlnum ) )
483  DO 70 j = jfirst, jlast, jinc
484 *
485 * Exit the loop if the growth factor is too small.
486 *
487  IF( grow.LE.smlnum )
488  $ GO TO 80
489 *
490 * G(j) = ( 1 + CNORM(j) )*G(j-1)
491 *
492  xj = one + cnorm( j )
493  grow = grow / xj
494  70 CONTINUE
495  END IF
496  80 CONTINUE
497  END IF
498 *
499  IF( ( grow*tscal ).GT.smlnum ) THEN
500 *
501 * Use the Level 2 BLAS solve if the reciprocal of the bound on
502 * elements of X is not too small.
503 *
504  CALL strsv( uplo, trans, diag, n, a, lda, x, 1 )
505  ELSE
506 *
507 * Use a Level 1 BLAS solve, scaling intermediate results.
508 *
509  IF( xmax.GT.bignum ) THEN
510 *
511 * Scale X so that its components are less than or equal to
512 * BIGNUM in absolute value.
513 *
514  scale = bignum / xmax
515  CALL sscal( n, scale, x, 1 )
516  xmax = bignum
517  END IF
518 *
519  IF( notran ) THEN
520 *
521 * Solve A * x = b
522 *
523  DO 100 j = jfirst, jlast, jinc
524 *
525 * Compute x(j) = b(j) / A(j,j), scaling x if necessary.
526 *
527  xj = abs( x( j ) )
528  IF( nounit ) THEN
529  tjjs = a( j, j )*tscal
530  ELSE
531  tjjs = tscal
532  IF( tscal.EQ.one )
533  $ GO TO 95
534  END IF
535  tjj = abs( tjjs )
536  IF( tjj.GT.smlnum ) THEN
537 *
538 * abs(A(j,j)) > SMLNUM:
539 *
540  IF( tjj.LT.one ) THEN
541  IF( xj.GT.tjj*bignum ) THEN
542 *
543 * Scale x by 1/b(j).
544 *
545  rec = one / xj
546  CALL sscal( n, rec, x, 1 )
547  scale = scale*rec
548  xmax = xmax*rec
549  END IF
550  END IF
551  x( j ) = x( j ) / tjjs
552  xj = abs( x( j ) )
553  ELSE IF( tjj.GT.zero ) THEN
554 *
555 * 0 < abs(A(j,j)) <= SMLNUM:
556 *
557  IF( xj.GT.tjj*bignum ) THEN
558 *
559 * Scale x by (1/abs(x(j)))*abs(A(j,j))*BIGNUM
560 * to avoid overflow when dividing by A(j,j).
561 *
562  rec = ( tjj*bignum ) / xj
563  IF( cnorm( j ).GT.one ) THEN
564 *
565 * Scale by 1/CNORM(j) to avoid overflow when
566 * multiplying x(j) times column j.
567 *
568  rec = rec / cnorm( j )
569  END IF
570  CALL sscal( n, rec, x, 1 )
571  scale = scale*rec
572  xmax = xmax*rec
573  END IF
574  x( j ) = x( j ) / tjjs
575  xj = abs( x( j ) )
576  ELSE
577 *
578 * A(j,j) = 0: Set x(1:n) = 0, x(j) = 1, and
579 * scale = 0, and compute a solution to A*x = 0.
580 *
581  DO 90 i = 1, n
582  x( i ) = zero
583  90 CONTINUE
584  x( j ) = one
585  xj = one
586  scale = zero
587  xmax = zero
588  END IF
589  95 CONTINUE
590 *
591 * Scale x if necessary to avoid overflow when adding a
592 * multiple of column j of A.
593 *
594  IF( xj.GT.one ) THEN
595  rec = one / xj
596  IF( cnorm( j ).GT.( bignum-xmax )*rec ) THEN
597 *
598 * Scale x by 1/(2*abs(x(j))).
599 *
600  rec = rec*half
601  CALL sscal( n, rec, x, 1 )
602  scale = scale*rec
603  END IF
604  ELSE IF( xj*cnorm( j ).GT.( bignum-xmax ) ) THEN
605 *
606 * Scale x by 1/2.
607 *
608  CALL sscal( n, half, x, 1 )
609  scale = scale*half
610  END IF
611 *
612  IF( upper ) THEN
613  IF( j.GT.1 ) THEN
614 *
615 * Compute the update
616 * x(1:j-1) := x(1:j-1) - x(j) * A(1:j-1,j)
617 *
618  CALL saxpy( j-1, -x( j )*tscal, a( 1, j ), 1, x,
619  $ 1 )
620  i = isamax( j-1, x, 1 )
621  xmax = abs( x( i ) )
622  END IF
623  ELSE
624  IF( j.LT.n ) THEN
625 *
626 * Compute the update
627 * x(j+1:n) := x(j+1:n) - x(j) * A(j+1:n,j)
628 *
629  CALL saxpy( n-j, -x( j )*tscal, a( j+1, j ), 1,
630  $ x( j+1 ), 1 )
631  i = j + isamax( n-j, x( j+1 ), 1 )
632  xmax = abs( x( i ) )
633  END IF
634  END IF
635  100 CONTINUE
636 *
637  ELSE
638 *
639 * Solve A**T * x = b
640 *
641  DO 140 j = jfirst, jlast, jinc
642 *
643 * Compute x(j) = b(j) - sum A(k,j)*x(k).
644 * k<>j
645 *
646  xj = abs( x( j ) )
647  uscal = tscal
648  rec = one / max( xmax, one )
649  IF( cnorm( j ).GT.( bignum-xj )*rec ) THEN
650 *
651 * If x(j) could overflow, scale x by 1/(2*XMAX).
652 *
653  rec = rec*half
654  IF( nounit ) THEN
655  tjjs = a( j, j )*tscal
656  ELSE
657  tjjs = tscal
658  END IF
659  tjj = abs( tjjs )
660  IF( tjj.GT.one ) THEN
661 *
662 * Divide by A(j,j) when scaling x if A(j,j) > 1.
663 *
664  rec = min( one, rec*tjj )
665  uscal = uscal / tjjs
666  END IF
667  IF( rec.LT.one ) THEN
668  CALL sscal( n, rec, x, 1 )
669  scale = scale*rec
670  xmax = xmax*rec
671  END IF
672  END IF
673 *
674  sumj = zero
675  IF( uscal.EQ.one ) THEN
676 *
677 * If the scaling needed for A in the dot product is 1,
678 * call SDOT to perform the dot product.
679 *
680  IF( upper ) THEN
681  sumj = sdot( j-1, a( 1, j ), 1, x, 1 )
682  ELSE IF( j.LT.n ) THEN
683  sumj = sdot( n-j, a( j+1, j ), 1, x( j+1 ), 1 )
684  END IF
685  ELSE
686 *
687 * Otherwise, use in-line code for the dot product.
688 *
689  IF( upper ) THEN
690  DO 110 i = 1, j - 1
691  sumj = sumj + ( a( i, j )*uscal )*x( i )
692  110 CONTINUE
693  ELSE IF( j.LT.n ) THEN
694  DO 120 i = j + 1, n
695  sumj = sumj + ( a( i, j )*uscal )*x( i )
696  120 CONTINUE
697  END IF
698  END IF
699 *
700  IF( uscal.EQ.tscal ) THEN
701 *
702 * Compute x(j) := ( x(j) - sumj ) / A(j,j) if 1/A(j,j)
703 * was not used to scale the dotproduct.
704 *
705  x( j ) = x( j ) - sumj
706  xj = abs( x( j ) )
707  IF( nounit ) THEN
708  tjjs = a( j, j )*tscal
709  ELSE
710  tjjs = tscal
711  IF( tscal.EQ.one )
712  $ GO TO 135
713  END IF
714 *
715 * Compute x(j) = x(j) / A(j,j), scaling if necessary.
716 *
717  tjj = abs( tjjs )
718  IF( tjj.GT.smlnum ) THEN
719 *
720 * abs(A(j,j)) > SMLNUM:
721 *
722  IF( tjj.LT.one ) THEN
723  IF( xj.GT.tjj*bignum ) THEN
724 *
725 * Scale X by 1/abs(x(j)).
726 *
727  rec = one / xj
728  CALL sscal( n, rec, x, 1 )
729  scale = scale*rec
730  xmax = xmax*rec
731  END IF
732  END IF
733  x( j ) = x( j ) / tjjs
734  ELSE IF( tjj.GT.zero ) THEN
735 *
736 * 0 < abs(A(j,j)) <= SMLNUM:
737 *
738  IF( xj.GT.tjj*bignum ) THEN
739 *
740 * Scale x by (1/abs(x(j)))*abs(A(j,j))*BIGNUM.
741 *
742  rec = ( tjj*bignum ) / xj
743  CALL sscal( n, rec, x, 1 )
744  scale = scale*rec
745  xmax = xmax*rec
746  END IF
747  x( j ) = x( j ) / tjjs
748  ELSE
749 *
750 * A(j,j) = 0: Set x(1:n) = 0, x(j) = 1, and
751 * scale = 0, and compute a solution to A**T*x = 0.
752 *
753  DO 130 i = 1, n
754  x( i ) = zero
755  130 CONTINUE
756  x( j ) = one
757  scale = zero
758  xmax = zero
759  END IF
760  135 CONTINUE
761  ELSE
762 *
763 * Compute x(j) := x(j) / A(j,j) - sumj if the dot
764 * product has already been divided by 1/A(j,j).
765 *
766  x( j ) = x( j ) / tjjs - sumj
767  END IF
768  xmax = max( xmax, abs( x( j ) ) )
769  140 CONTINUE
770  END IF
771  scale = scale / tscal
772  END IF
773 *
774 * Scale the column norms by 1/TSCAL for return.
775 *
776  IF( tscal.NE.one ) THEN
777  CALL sscal( n, one / tscal, cnorm, 1 )
778  END IF
779 *
780  RETURN
781 *
782 * End of SLATRS
783 *
integer function isamax(N, SX, INCX)
ISAMAX
Definition: isamax.f:71
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine sscal(N, SA, SX, INCX)
SSCAL
Definition: sscal.f:79
real function sdot(N, SX, INCX, SY, INCY)
SDOT
Definition: sdot.f:82
subroutine saxpy(N, SA, SX, INCX, SY, INCY)
SAXPY
Definition: saxpy.f:89
real function sasum(N, SX, INCX)
SASUM
Definition: sasum.f:72
subroutine strsv(UPLO, TRANS, DIAG, N, A, LDA, X, INCX)
STRSV
Definition: strsv.f:149
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:68
Here is the call graph for this function:
Here is the caller graph for this function: