cptcon.f.html

来自「famous linear algebra library (LAPACK) p」· HTML 代码 · 共 175 行

HTML
175 行
字号 
<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN">
<html>
 <head>
  <title>cptcon.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="CPTCON.1"></a><a href="cptcon.f.html#CPTCON.1">CPTCON</a>( N, D, E, ANORM, RCOND, RWORK, INFO )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  -- LAPACK 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>      INTEGER            INFO, N
      REAL               ANORM, RCOND
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Array Arguments ..
</span>      REAL               D( * ), RWORK( * )
      COMPLEX            E( * )
<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="CPTCON.19"></a><a href="cptcon.f.html#CPTCON.1">CPTCON</a> computes the reciprocal of the condition number (in the
</span><span class="comment">*</span><span class="comment">  1-norm) of a complex Hermitian positive definite tridiagonal matrix
</span><span class="comment">*</span><span class="comment">  using the factorization A = L*D*L**H or A = U**H*D*U computed by
</span><span class="comment">*</span><span class="comment">  <a name="CPTTRF.22"></a><a href="cpttrf.f.html#CPTTRF.1">CPTTRF</a>.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  Norm(inv(A)) is computed by a direct method, and the reciprocal of
</span><span class="comment">*</span><span class="comment">  the condition number is computed as
</span><span class="comment">*</span><span class="comment">                   RCOND = 1 / (ANORM * norm(inv(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">  N       (input) INTEGER
</span><span class="comment">*</span><span class="comment">          The order of the matrix A.  N &gt;= 0.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  D       (input) REAL array, dimension (N)
</span><span class="comment">*</span><span class="comment">          The n diagonal elements of the diagonal matrix D from the
</span><span class="comment">*</span><span class="comment">          factorization of A, as computed by <a name="CPTTRF.36"></a><a href="cpttrf.f.html#CPTTRF.1">CPTTRF</a>.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  E       (input) COMPLEX array, dimension (N-1)
</span><span class="comment">*</span><span class="comment">          The (n-1) off-diagonal elements of the unit bidiagonal factor
</span><span class="comment">*</span><span class="comment">          U or L from the factorization of A, as computed by <a name="CPTTRF.40"></a><a href="cpttrf.f.html#CPTTRF.1">CPTTRF</a>.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  ANORM   (input) REAL
</span><span class="comment">*</span><span class="comment">          The 1-norm of the original matrix A.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  RCOND   (output) REAL
</span><span class="comment">*</span><span class="comment">          The reciprocal of the condition number of the matrix A,
</span><span class="comment">*</span><span class="comment">          computed as RCOND = 1/(ANORM * AINVNM), where AINVNM is the
</span><span class="comment">*</span><span class="comment">          1-norm of inv(A) computed in this routine.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  RWORK   (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">          = 0:  successful exit
</span><span class="comment">*</span><span class="comment">          &lt; 0:  if INFO = -i, the i-th argument had an illegal value
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  Further Details
</span><span class="comment">*</span><span class="comment">  ===============
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  The method used is described in Nicholas J. Higham, &quot;Efficient
</span><span class="comment">*</span><span class="comment">  Algorithms for Computing the Condition Number of a Tridiagonal
</span><span class="comment">*</span><span class="comment">  Matrix&quot;, SIAM J. Sci. Stat. Comput., Vol. 7, No. 1, January 1986.
</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               ONE, ZERO
      PARAMETER          ( ONE = 1.0E+0, ZERO = 0.0E+0 )
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Local Scalars ..
</span>      INTEGER            I, IX
      REAL               AINVNM
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. External Functions ..
</span>      INTEGER            ISAMAX
      EXTERNAL           ISAMAX
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. External Subroutines ..
</span>      EXTERNAL           <a name="XERBLA.78"></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
<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">     Test the input arguments.
</span><span class="comment">*</span><span class="comment">
</span>      INFO = 0
      IF( N.LT.0 ) THEN
         INFO = -1
      ELSE IF( ANORM.LT.ZERO ) THEN
         INFO = -4
      END IF
      IF( INFO.NE.0 ) THEN
         CALL <a name="XERBLA.94"></a><a href="xerbla.f.html#XERBLA.1">XERBLA</a>( <span class="string">'<a name="CPTCON.94"></a><a href="cptcon.f.html#CPTCON.1">CPTCON</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>      RCOND = ZERO
      IF( N.EQ.0 ) THEN
         RCOND = ONE
         RETURN
      ELSE IF( ANORM.EQ.ZERO ) THEN
         RETURN
      END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Check that D(1:N) is positive.
</span><span class="comment">*</span><span class="comment">
</span>      DO 10 I = 1, N
         IF( D( I ).LE.ZERO )
     $      RETURN
   10 CONTINUE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Solve M(A) * x = e, where M(A) = (m(i,j)) is given by
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        m(i,j) =  abs(A(i,j)), i = j,
</span><span class="comment">*</span><span class="comment">        m(i,j) = -abs(A(i,j)), i .ne. j,
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     and e = [ 1, 1, ..., 1 ]'.  Note M(A) = M(L)*D*M(L)'.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Solve M(L) * x = e.
</span><span class="comment">*</span><span class="comment">
</span>      RWORK( 1 ) = ONE
      DO 20 I = 2, N
         RWORK( I ) = ONE + RWORK( I-1 )*ABS( E( I-1 ) )
   20 CONTINUE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Solve D * M(L)' * x = b.
</span><span class="comment">*</span><span class="comment">
</span>      RWORK( N ) = RWORK( N ) / D( N )
      DO 30 I = N - 1, 1, -1
         RWORK( I ) = RWORK( I ) / D( I ) + RWORK( I+1 )*ABS( E( I ) )
   30 CONTINUE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Compute AINVNM = max(x(i)), 1&lt;=i&lt;=n.
</span><span class="comment">*</span><span class="comment">
</span>      IX = ISAMAX( N, RWORK, 1 )
      AINVNM = ABS( RWORK( IX ) )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Compute the reciprocal condition number.
</span><span class="comment">*</span><span class="comment">
</span>      IF( AINVNM.NE.ZERO )
     $   RCOND = ( ONE / AINVNM ) / ANORM
<span class="comment">*</span><span class="comment">
</span>      RETURN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     End of <a name="CPTCON.148"></a><a href="cptcon.f.html#CPTCON.1">CPTCON</a>
</span><span class="comment">*</span><span class="comment">
</span>      END

</pre>

 </body>
</html>

cptcon.f.html - 源码说明

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

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

⌨️ 快捷键说明

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