⭐ 欢迎来到虫虫下载站！ | 📦 资源下载 📁 资源专辑 ℹ️ 关于我们
⭐ 虫虫下载站
📤 上传资源

📄 slaqtr.f.html

📁 famous linear algebra library (LAPACK) ports to windows
💻 HTML
📖 第 1 页 / 共 3 页
字号:
🔍
123 
<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN">
<html>
 <head>
  <title>slaqtr.f</title>
 <meta name="generator" content="emacs 21.3.1; htmlfontify 0.20">
<style type="text/css"><!-- 
body { background: rgb(255, 255, 255);  color: rgb(0, 0, 0);  font-style: normal;  font-weight: 500;  font-stretch: normal;  font-family: adobe-courier;  font-size: 11pt;  text-decoration: none; }
span.default   { background: rgb(255, 255, 255);  color: rgb(0, 0, 0);  font-style: normal;  font-weight: 500;  font-stretch: normal;  font-family: adobe-courier;  font-size: 11pt;  text-decoration: none; }
span.default a { background: rgb(255, 255, 255);  color: rgb(0, 0, 0);  font-style: normal;  font-weight: 500;  font-stretch: normal;  font-family: adobe-courier;  font-size: 11pt;  text-decoration: underline; }
span.string   { color: rgb(188, 143, 143);  background: rgb(255, 255, 255);  font-style: normal;  font-weight: 500;  font-stretch: normal;  font-family: adobe-courier;  font-size: 11pt;  text-decoration: none; }
span.string a { color: rgb(188, 143, 143);  background: rgb(255, 255, 255);  font-style: normal;  font-weight: 500;  font-stretch: normal;  font-family: adobe-courier;  font-size: 11pt;  text-decoration: underline; }
span.comment   { color: rgb(178, 34, 34);  background: rgb(255, 255, 255);  font-style: normal;  font-weight: 500;  font-stretch: normal;  font-family: adobe-courier;  font-size: 11pt;  text-decoration: none; }
span.comment a { color: rgb(178, 34, 34);  background: rgb(255, 255, 255);  font-style: normal;  font-weight: 500;  font-stretch: normal;  font-family: adobe-courier;  font-size: 11pt;  text-decoration: underline; }
 --></style>

 </head>
  <body>

<pre>
      SUBROUTINE <a name="SLAQTR.1"></a><a href="slaqtr.f.html#SLAQTR.1">SLAQTR</a>( LTRAN, LREAL, N, T, LDT, B, W, SCALE, X, WORK,
     $                   INFO )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  -- LAPACK auxiliary routine (version 3.1) --
</span><span class="comment">*</span><span class="comment">     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
</span><span class="comment">*</span><span class="comment">     November 2006
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     .. Scalar Arguments ..
</span>      LOGICAL            LREAL, LTRAN
      INTEGER            INFO, LDT, N
      REAL               SCALE, W
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Array Arguments ..
</span>      REAL               B( * ), T( LDT, * ), WORK( * ), X( * )
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  Purpose
</span><span class="comment">*</span><span class="comment">  =======
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  <a name="SLAQTR.20"></a><a href="slaqtr.f.html#SLAQTR.1">SLAQTR</a> solves the real quasi-triangular system
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">               op(T)*p = scale*c,               if LREAL = .TRUE.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  or the complex quasi-triangular systems
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">             op(T + iB)*(p+iq) = scale*(c+id),  if LREAL = .FALSE.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  in real arithmetic, where T is upper quasi-triangular.
</span><span class="comment">*</span><span class="comment">  If LREAL = .FALSE., then the first diagonal block of T must be
</span><span class="comment">*</span><span class="comment">  1 by 1, B is the specially structured matrix
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                 B = [ b(1) b(2) ... b(n) ]
</span><span class="comment">*</span><span class="comment">                     [       w            ]
</span><span class="comment">*</span><span class="comment">                     [           w        ]
</span><span class="comment">*</span><span class="comment">                     [              .     ]
</span><span class="comment">*</span><span class="comment">                     [                 w  ]
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  op(A) = A or A', A' denotes the conjugate transpose of
</span><span class="comment">*</span><span class="comment">  matrix A.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  On input, X = [ c ].  On output, X = [ p ].
</span><span class="comment">*</span><span class="comment">                [ d ]                  [ q ]
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  This subroutine is designed for the condition number estimation
</span><span class="comment">*</span><span class="comment">  in routine <a name="STRSNA.45"></a><a href="strsna.f.html#STRSNA.1">STRSNA</a>.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  Arguments
</span><span class="comment">*</span><span class="comment">  =========
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  LTRAN   (input) LOGICAL
</span><span class="comment">*</span><span class="comment">          On entry, LTRAN specifies the option of conjugate transpose:
</span><span class="comment">*</span><span class="comment">             = .FALSE.,    op(T+i*B) = T+i*B,
</span><span class="comment">*</span><span class="comment">             = .TRUE.,     op(T+i*B) = (T+i*B)'.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  LREAL   (input) LOGICAL
</span><span class="comment">*</span><span class="comment">          On entry, LREAL specifies the input matrix structure:
</span><span class="comment">*</span><span class="comment">             = .FALSE.,    the input is complex
</span><span class="comment">*</span><span class="comment">             = .TRUE.,     the input is real
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  N       (input) INTEGER
</span><span class="comment">*</span><span class="comment">          On entry, N specifies the order of T+i*B. N &gt;= 0.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  T       (input) REAL array, dimension (LDT,N)
</span><span class="comment">*</span><span class="comment">          On entry, T contains a matrix in Schur canonical form.
</span><span class="comment">*</span><span class="comment">          If LREAL = .FALSE., then the first diagonal block of T must
</span><span class="comment">*</span><span class="comment">          be 1 by 1.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  LDT     (input) INTEGER
</span><span class="comment">*</span><span class="comment">          The leading dimension of the matrix T. LDT &gt;= max(1,N).
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  B       (input) REAL array, dimension (N)
</span><span class="comment">*</span><span class="comment">          On entry, B contains the elements to form the matrix
</span><span class="comment">*</span><span class="comment">          B as described above.
</span><span class="comment">*</span><span class="comment">          If LREAL = .TRUE., B is not referenced.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  W       (input) REAL
</span><span class="comment">*</span><span class="comment">          On entry, W is the diagonal element of the matrix B.
</span><span class="comment">*</span><span class="comment">          If LREAL = .TRUE., W is not referenced.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  SCALE   (output) REAL
</span><span class="comment">*</span><span class="comment">          On exit, SCALE is the scale factor.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  X       (input/output) REAL array, dimension (2*N)
</span><span class="comment">*</span><span class="comment">          On entry, X contains the right hand side of the system.
</span><span class="comment">*</span><span class="comment">          On exit, X is overwritten by the solution.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  WORK    (workspace) REAL array, dimension (N)
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  INFO    (output) INTEGER
</span><span class="comment">*</span><span class="comment">          On exit, INFO is set to
</span><span class="comment">*</span><span class="comment">             0: successful exit.
</span><span class="comment">*</span><span class="comment">               1: the some diagonal 1 by 1 block has been perturbed by
</span><span class="comment">*</span><span class="comment">                  a small number SMIN to keep nonsingularity.
</span><span class="comment">*</span><span class="comment">               2: the some diagonal 2 by 2 block has been perturbed by
</span><span class="comment">*</span><span class="comment">                  a small number in <a name="SLALN2.95"></a><a href="slaln2.f.html#SLALN2.1">SLALN2</a> to keep nonsingularity.
</span><span class="comment">*</span><span class="comment">          NOTE: In the interests of speed, this routine does not
</span><span class="comment">*</span><span class="comment">                check the inputs for errors.
</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, ONE
      PARAMETER          ( ZERO = 0.0E+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
      INTEGER            I, IERR, J, J1, J2, JNEXT, K, N1, N2
      REAL               BIGNUM, EPS, REC, SCALOC, SI, SMIN, SMINW,
     $                   SMLNUM, SR, TJJ, TMP, XJ, XMAX, XNORM, Z
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Local Arrays ..
</span>      REAL               D( 2, 2 ), V( 2, 2 )
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. External Functions ..
</span>      INTEGER            ISAMAX
      REAL               SASUM, SDOT, <a name="SLAMCH.116"></a><a href="slamch.f.html#SLAMCH.1">SLAMCH</a>, <a name="SLANGE.116"></a><a href="slange.f.html#SLANGE.1">SLANGE</a>
      EXTERNAL           ISAMAX, SASUM, SDOT, <a name="SLAMCH.117"></a><a href="slamch.f.html#SLAMCH.1">SLAMCH</a>, <a name="SLANGE.117"></a><a href="slange.f.html#SLANGE.1">SLANGE</a>
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. External Subroutines ..
</span>      EXTERNAL           SAXPY, <a name="SLADIV.120"></a><a href="sladiv.f.html#SLADIV.1">SLADIV</a>, <a name="SLALN2.120"></a><a href="slaln2.f.html#SLALN2.1">SLALN2</a>, SSCAL
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Intrinsic Functions ..
</span>      INTRINSIC          ABS, MAX
<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><span class="comment">*</span><span class="comment">     Do not test the input parameters for errors
</span><span class="comment">*</span><span class="comment">
</span>      NOTRAN = .NOT.LTRAN
      INFO = 0
<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">     Set constants to control overflow
</span><span class="comment">*</span><span class="comment">
</span>      EPS = <a name="SLAMCH.139"></a><a href="slamch.f.html#SLAMCH.1">SLAMCH</a>( <span class="string">'P'</span> )
      SMLNUM = <a name="SLAMCH.140"></a><a href="slamch.f.html#SLAMCH.1">SLAMCH</a>( <span class="string">'S'</span> ) / EPS
      BIGNUM = ONE / SMLNUM
<span class="comment">*</span><span class="comment">
</span>      XNORM = <a name="SLANGE.143"></a><a href="slange.f.html#SLANGE.1">SLANGE</a>( <span class="string">'M'</span>, N, N, T, LDT, D )
      IF( .NOT.LREAL )
     $   XNORM = MAX( XNORM, ABS( W ), <a name="SLANGE.145"></a><a href="slange.f.html#SLANGE.1">SLANGE</a>( <span class="string">'M'</span>, N, 1, B, N, D ) )
      SMIN = MAX( SMLNUM, EPS*XNORM )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Compute 1-norm of each column of strictly upper triangular
</span><span class="comment">*</span><span class="comment">     part of T to control overflow in triangular solver.
</span><span class="comment">*</span><span class="comment">
</span>      WORK( 1 ) = ZERO
      DO 10 J = 2, N
         WORK( J ) = SASUM( J-1, T( 1, J ), 1 )
   10 CONTINUE
<span class="comment">*</span><span class="comment">
</span>      IF( .NOT.LREAL ) THEN
         DO 20 I = 2, N
            WORK( I ) = WORK( I ) + ABS( B( I ) )
   20    CONTINUE
      END IF
<span class="comment">*</span><span class="comment">
</span>      N2 = 2*N
      N1 = N
      IF( .NOT.LREAL )
     $   N1 = N2
      K = ISAMAX( N1, X, 1 )
      XMAX = ABS( X( K ) )
      SCALE = ONE
<span class="comment">*</span><span class="comment">
</span>      IF( XMAX.GT.BIGNUM ) THEN
         SCALE = BIGNUM / XMAX
         CALL SSCAL( N1, SCALE, X, 1 )
         XMAX = BIGNUM
      END IF
<span class="comment">*</span><span class="comment">
</span>      IF( LREAL ) THEN
<span class="comment">*</span><span class="comment">
</span>         IF( NOTRAN ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Solve T*p = scale*c
</span><span class="comment">*</span><span class="comment">
</span>            JNEXT = N
            DO 30 J = N, 1, -1
               IF( J.GT.JNEXT )
     $            GO TO 30
               J1 = J
               J2 = J
               JNEXT = J - 1
               IF( J.GT.1 ) THEN
                  IF( T( J, J-1 ).NE.ZERO ) THEN
                     J1 = J - 1
                     JNEXT = J - 2
                  END IF
               END IF
<span class="comment">*</span><span class="comment">
</span>               IF( J1.EQ.J2 ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                 Meet 1 by 1 diagonal block
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                 Scale to avoid overflow when computing
</span><span class="comment">*</span><span class="comment">                     x(j) = b(j)/T(j,j)
</span><span class="comment">*</span><span class="comment">
</span>                  XJ = ABS( X( J1 ) )
                  TJJ = ABS( T( J1, J1 ) )
                  TMP = T( J1, J1 )
                  IF( TJJ.LT.SMIN ) THEN
                     TMP = SMIN
                     TJJ = SMIN
                     INFO = 1
                  END IF
<span class="comment">*</span><span class="comment">
123

⌨️ 快捷键说明

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