slatbs.f.html

来自「famous linear algebra library (LAPACK) p」· HTML 代码 · 共 748 行 · 第 1/4 页

HTML
748 行
字号
1234
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  =====================================================================
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     .. Parameters ..
</span>      REAL               ZERO, HALF, ONE
      PARAMETER          ( ZERO = 0.0E+0, HALF = 0.5E+0, ONE = 1.0E+0 )
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Local Scalars ..
</span>      LOGICAL            NOTRAN, NOUNIT, UPPER
      INTEGER            I, IMAX, J, JFIRST, JINC, JLAST, JLEN, MAIND
      REAL               BIGNUM, GROW, REC, SMLNUM, SUMJ, TJJ, TJJS,
     $                   TMAX, TSCAL, USCAL, XBND, XJ, XMAX
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. External Functions ..
</span>      LOGICAL            <a name="LSAME.183"></a><a href="lsame.f.html#LSAME.1">LSAME</a>
      INTEGER            ISAMAX
      REAL               SASUM, SDOT, <a name="SLAMCH.185"></a><a href="slamch.f.html#SLAMCH.1">SLAMCH</a>
      EXTERNAL           <a name="LSAME.186"></a><a href="lsame.f.html#LSAME.1">LSAME</a>, ISAMAX, SASUM, SDOT, <a name="SLAMCH.186"></a><a href="slamch.f.html#SLAMCH.1">SLAMCH</a>
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. External Subroutines ..
</span>      EXTERNAL           SAXPY, SSCAL, STBSV, <a name="XERBLA.189"></a><a href="xerbla.f.html#XERBLA.1">XERBLA</a>
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Intrinsic Functions ..
</span>      INTRINSIC          ABS, MAX, MIN
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Executable Statements ..
</span><span class="comment">*</span><span class="comment">
</span>      INFO = 0
      UPPER = <a name="LSAME.197"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( UPLO, <span class="string">'U'</span> )
      NOTRAN = <a name="LSAME.198"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( TRANS, <span class="string">'N'</span> )
      NOUNIT = <a name="LSAME.199"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( DIAG, <span class="string">'N'</span> )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Test the input parameters.
</span><span class="comment">*</span><span class="comment">
</span>      IF( .NOT.UPPER .AND. .NOT.<a name="LSAME.203"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( UPLO, <span class="string">'L'</span> ) ) THEN
         INFO = -1
      ELSE IF( .NOT.NOTRAN .AND. .NOT.<a name="LSAME.205"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( TRANS, <span class="string">'T'</span> ) .AND. .NOT.
     $         <a name="LSAME.206"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( TRANS, <span class="string">'C'</span> ) ) THEN
         INFO = -2
      ELSE IF( .NOT.NOUNIT .AND. .NOT.<a name="LSAME.208"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( DIAG, <span class="string">'U'</span> ) ) THEN
         INFO = -3
      ELSE IF( .NOT.<a name="LSAME.210"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( NORMIN, <span class="string">'Y'</span> ) .AND. .NOT.
     $         <a name="LSAME.211"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( NORMIN, <span class="string">'N'</span> ) ) THEN
         INFO = -4
      ELSE IF( N.LT.0 ) THEN
         INFO = -5
      ELSE IF( KD.LT.0 ) THEN
         INFO = -6
      ELSE IF( LDAB.LT.KD+1 ) THEN
         INFO = -8
      END IF
      IF( INFO.NE.0 ) THEN
         CALL <a name="XERBLA.221"></a><a href="xerbla.f.html#XERBLA.1">XERBLA</a>( <span class="string">'<a name="SLATBS.221"></a><a href="slatbs.f.html#SLATBS.1">SLATBS</a>'</span>, -INFO )
         RETURN
      END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Quick return if possible
</span><span class="comment">*</span><span class="comment">
</span>      IF( N.EQ.0 )
     $   RETURN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Determine machine dependent parameters to control overflow.
</span><span class="comment">*</span><span class="comment">
</span>      SMLNUM = <a name="SLAMCH.232"></a><a href="slamch.f.html#SLAMCH.1">SLAMCH</a>( <span class="string">'Safe minimum'</span> ) / <a name="SLAMCH.232"></a><a href="slamch.f.html#SLAMCH.1">SLAMCH</a>( <span class="string">'Precision'</span> )
      BIGNUM = ONE / SMLNUM
      SCALE = ONE
<span class="comment">*</span><span class="comment">
</span>      IF( <a name="LSAME.236"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( NORMIN, <span class="string">'N'</span> ) ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Compute the 1-norm of each column, not including the diagonal.
</span><span class="comment">*</span><span class="comment">
</span>         IF( UPPER ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           A is upper triangular.
</span><span class="comment">*</span><span class="comment">
</span>            DO 10 J = 1, N
               JLEN = MIN( KD, J-1 )
               CNORM( J ) = SASUM( JLEN, AB( KD+1-JLEN, J ), 1 )
   10       CONTINUE
         ELSE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           A is lower triangular.
</span><span class="comment">*</span><span class="comment">
</span>            DO 20 J = 1, N
               JLEN = MIN( KD, N-J )
               IF( JLEN.GT.0 ) THEN
                  CNORM( J ) = SASUM( JLEN, AB( 2, J ), 1 )
               ELSE
                  CNORM( J ) = ZERO
               END IF
   20       CONTINUE
         END IF
      END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Scale the column norms by TSCAL if the maximum element in CNORM is
</span><span class="comment">*</span><span class="comment">     greater than BIGNUM.
</span><span class="comment">*</span><span class="comment">
</span>      IMAX = ISAMAX( N, CNORM, 1 )
      TMAX = CNORM( IMAX )
      IF( TMAX.LE.BIGNUM ) THEN
         TSCAL = ONE
      ELSE
         TSCAL = ONE / ( SMLNUM*TMAX )
         CALL SSCAL( N, TSCAL, CNORM, 1 )
      END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Compute a bound on the computed solution vector to see if the
</span><span class="comment">*</span><span class="comment">     Level 2 BLAS routine STBSV can be used.
</span><span class="comment">*</span><span class="comment">
</span>      J = ISAMAX( N, X, 1 )
      XMAX = ABS( X( J ) )
      XBND = XMAX
      IF( NOTRAN ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Compute the growth in A * x = b.
</span><span class="comment">*</span><span class="comment">
</span>         IF( UPPER ) THEN
            JFIRST = N
            JLAST = 1
            JINC = -1
            MAIND = KD + 1
         ELSE
            JFIRST = 1
            JLAST = N
            JINC = 1
            MAIND = 1
         END IF
<span class="comment">*</span><span class="comment">
</span>         IF( TSCAL.NE.ONE ) THEN
            GROW = ZERO
            GO TO 50
         END IF
<span class="comment">*</span><span class="comment">
</span>         IF( NOUNIT ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           A is non-unit triangular.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Compute GROW = 1/G(j) and XBND = 1/M(j).
</span><span class="comment">*</span><span class="comment">           Initially, G(0) = max{x(i), i=1,...,n}.
</span><span class="comment">*</span><span class="comment">
</span>            GROW = ONE / MAX( XBND, SMLNUM )
            XBND = GROW
            DO 30 J = JFIRST, JLAST, JINC
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">              Exit the loop if the growth factor is too small.
</span><span class="comment">*</span><span class="comment">
</span>               IF( GROW.LE.SMLNUM )
     $            GO TO 50
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">              M(j) = G(j-1) / abs(A(j,j))
</span><span class="comment">*</span><span class="comment">
</span>               TJJ = ABS( AB( MAIND, J ) )
               XBND = MIN( XBND, MIN( ONE, TJJ )*GROW )
               IF( TJJ+CNORM( J ).GE.SMLNUM ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                 G(j) = G(j-1)*( 1 + CNORM(j) / abs(A(j,j)) )
</span><span class="comment">*</span><span class="comment">
</span>                  GROW = GROW*( TJJ / ( TJJ+CNORM( J ) ) )
               ELSE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                 G(j) could overflow, set GROW to 0.
</span><span class="comment">*</span><span class="comment">
</span>                  GROW = ZERO
               END IF
   30       CONTINUE
            GROW = XBND
         ELSE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           A is unit triangular.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Compute GROW = 1/G(j), where G(0) = max{x(i), i=1,...,n}.
</span><span class="comment">*</span><span class="comment">
</span>            GROW = MIN( ONE, ONE / MAX( XBND, SMLNUM ) )
            DO 40 J = JFIRST, JLAST, JINC
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">              Exit the loop if the growth factor is too small.
</span><span class="comment">*</span><span class="comment">
</span>               IF( GROW.LE.SMLNUM )
     $            GO TO 50
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">              G(j) = G(j-1)*( 1 + CNORM(j) )
</span><span class="comment">*</span><span class="comment">
</span>               GROW = GROW*( ONE / ( ONE+CNORM( J ) ) )
   40       CONTINUE
         END IF
   50    CONTINUE
<span class="comment">*</span><span class="comment">
1234

slatbs.f.html - 源码说明

本页面展示了「famous linear algebra library (LAPACK) ports to windows」中的 slatbs.f.html 源码文件，采用 HTML 编程语言编写，共 748 行代码。您可以在线阅读完整代码内容，也可以返回资源详情页下载完整源码包进行本地学习和开发。

虫虫下载站收录了大量与algebra相关的技术资源，包括源代码、技术文档、电路图等，是电子工程师和嵌入式开发者的专业学习平台。

⌨️ 快捷键说明

复制代码Ctrl + C
搜索代码Ctrl + F
全屏模式F11
增大字号Ctrl + =
减小字号Ctrl + -
显示快捷键?