Annotation of rpl/lapack/lapack/zptcon.f, revision 1.9
1.9 ! bertrand 1: *> \brief \b ZPTCON
! 2: *
! 3: * =========== DOCUMENTATION ===========
! 4: *
! 5: * Online html documentation available at
! 6: * http://www.netlib.org/lapack/explore-html/
! 7: *
! 8: *> \htmlonly
! 9: *> Download ZPTCON + dependencies
! 10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zptcon.f">
! 11: *> [TGZ]</a>
! 12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zptcon.f">
! 13: *> [ZIP]</a>
! 14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zptcon.f">
! 15: *> [TXT]</a>
! 16: *> \endhtmlonly
! 17: *
! 18: * Definition:
! 19: * ===========
! 20: *
! 21: * SUBROUTINE ZPTCON( N, D, E, ANORM, RCOND, RWORK, INFO )
! 22: *
! 23: * .. Scalar Arguments ..
! 24: * INTEGER INFO, N
! 25: * DOUBLE PRECISION ANORM, RCOND
! 26: * ..
! 27: * .. Array Arguments ..
! 28: * DOUBLE PRECISION D( * ), RWORK( * )
! 29: * COMPLEX*16 E( * )
! 30: * ..
! 31: *
! 32: *
! 33: *> \par Purpose:
! 34: * =============
! 35: *>
! 36: *> \verbatim
! 37: *>
! 38: *> ZPTCON computes the reciprocal of the condition number (in the
! 39: *> 1-norm) of a complex Hermitian positive definite tridiagonal matrix
! 40: *> using the factorization A = L*D*L**H or A = U**H*D*U computed by
! 41: *> ZPTTRF.
! 42: *>
! 43: *> Norm(inv(A)) is computed by a direct method, and the reciprocal of
! 44: *> the condition number is computed as
! 45: *> RCOND = 1 / (ANORM * norm(inv(A))).
! 46: *> \endverbatim
! 47: *
! 48: * Arguments:
! 49: * ==========
! 50: *
! 51: *> \param[in] N
! 52: *> \verbatim
! 53: *> N is INTEGER
! 54: *> The order of the matrix A. N >= 0.
! 55: *> \endverbatim
! 56: *>
! 57: *> \param[in] D
! 58: *> \verbatim
! 59: *> D is DOUBLE PRECISION array, dimension (N)
! 60: *> The n diagonal elements of the diagonal matrix D from the
! 61: *> factorization of A, as computed by ZPTTRF.
! 62: *> \endverbatim
! 63: *>
! 64: *> \param[in] E
! 65: *> \verbatim
! 66: *> E is COMPLEX*16 array, dimension (N-1)
! 67: *> The (n-1) off-diagonal elements of the unit bidiagonal factor
! 68: *> U or L from the factorization of A, as computed by ZPTTRF.
! 69: *> \endverbatim
! 70: *>
! 71: *> \param[in] ANORM
! 72: *> \verbatim
! 73: *> ANORM is DOUBLE PRECISION
! 74: *> The 1-norm of the original matrix A.
! 75: *> \endverbatim
! 76: *>
! 77: *> \param[out] RCOND
! 78: *> \verbatim
! 79: *> RCOND is DOUBLE PRECISION
! 80: *> The reciprocal of the condition number of the matrix A,
! 81: *> computed as RCOND = 1/(ANORM * AINVNM), where AINVNM is the
! 82: *> 1-norm of inv(A) computed in this routine.
! 83: *> \endverbatim
! 84: *>
! 85: *> \param[out] RWORK
! 86: *> \verbatim
! 87: *> RWORK is DOUBLE PRECISION array, dimension (N)
! 88: *> \endverbatim
! 89: *>
! 90: *> \param[out] INFO
! 91: *> \verbatim
! 92: *> INFO is INTEGER
! 93: *> = 0: successful exit
! 94: *> < 0: if INFO = -i, the i-th argument had an illegal value
! 95: *> \endverbatim
! 96: *
! 97: * Authors:
! 98: * ========
! 99: *
! 100: *> \author Univ. of Tennessee
! 101: *> \author Univ. of California Berkeley
! 102: *> \author Univ. of Colorado Denver
! 103: *> \author NAG Ltd.
! 104: *
! 105: *> \date November 2011
! 106: *
! 107: *> \ingroup complex16OTHERcomputational
! 108: *
! 109: *> \par Further Details:
! 110: * =====================
! 111: *>
! 112: *> \verbatim
! 113: *>
! 114: *> The method used is described in Nicholas J. Higham, "Efficient
! 115: *> Algorithms for Computing the Condition Number of a Tridiagonal
! 116: *> Matrix", SIAM J. Sci. Stat. Comput., Vol. 7, No. 1, January 1986.
! 117: *> \endverbatim
! 118: *>
! 119: * =====================================================================
1.1 bertrand 120: SUBROUTINE ZPTCON( N, D, E, ANORM, RCOND, RWORK, INFO )
121: *
1.9 ! bertrand 122: * -- LAPACK computational routine (version 3.4.0) --
1.1 bertrand 123: * -- LAPACK is a software package provided by Univ. of Tennessee, --
124: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
1.9 ! bertrand 125: * November 2011
1.1 bertrand 126: *
127: * .. Scalar Arguments ..
128: INTEGER INFO, N
129: DOUBLE PRECISION ANORM, RCOND
130: * ..
131: * .. Array Arguments ..
132: DOUBLE PRECISION D( * ), RWORK( * )
133: COMPLEX*16 E( * )
134: * ..
135: *
136: * =====================================================================
137: *
138: * .. Parameters ..
139: DOUBLE PRECISION ONE, ZERO
140: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
141: * ..
142: * .. Local Scalars ..
143: INTEGER I, IX
144: DOUBLE PRECISION AINVNM
145: * ..
146: * .. External Functions ..
147: INTEGER IDAMAX
148: EXTERNAL IDAMAX
149: * ..
150: * .. External Subroutines ..
151: EXTERNAL XERBLA
152: * ..
153: * .. Intrinsic Functions ..
154: INTRINSIC ABS
155: * ..
156: * .. Executable Statements ..
157: *
158: * Test the input arguments.
159: *
160: INFO = 0
161: IF( N.LT.0 ) THEN
162: INFO = -1
163: ELSE IF( ANORM.LT.ZERO ) THEN
164: INFO = -4
165: END IF
166: IF( INFO.NE.0 ) THEN
167: CALL XERBLA( 'ZPTCON', -INFO )
168: RETURN
169: END IF
170: *
171: * Quick return if possible
172: *
173: RCOND = ZERO
174: IF( N.EQ.0 ) THEN
175: RCOND = ONE
176: RETURN
177: ELSE IF( ANORM.EQ.ZERO ) THEN
178: RETURN
179: END IF
180: *
181: * Check that D(1:N) is positive.
182: *
183: DO 10 I = 1, N
184: IF( D( I ).LE.ZERO )
185: $ RETURN
186: 10 CONTINUE
187: *
188: * Solve M(A) * x = e, where M(A) = (m(i,j)) is given by
189: *
190: * m(i,j) = abs(A(i,j)), i = j,
191: * m(i,j) = -abs(A(i,j)), i .ne. j,
192: *
1.8 bertrand 193: * and e = [ 1, 1, ..., 1 ]**T. Note M(A) = M(L)*D*M(L)**H.
1.1 bertrand 194: *
195: * Solve M(L) * x = e.
196: *
197: RWORK( 1 ) = ONE
198: DO 20 I = 2, N
199: RWORK( I ) = ONE + RWORK( I-1 )*ABS( E( I-1 ) )
200: 20 CONTINUE
201: *
1.8 bertrand 202: * Solve D * M(L)**H * x = b.
1.1 bertrand 203: *
204: RWORK( N ) = RWORK( N ) / D( N )
205: DO 30 I = N - 1, 1, -1
206: RWORK( I ) = RWORK( I ) / D( I ) + RWORK( I+1 )*ABS( E( I ) )
207: 30 CONTINUE
208: *
209: * Compute AINVNM = max(x(i)), 1<=i<=n.
210: *
211: IX = IDAMAX( N, RWORK, 1 )
212: AINVNM = ABS( RWORK( IX ) )
213: *
214: * Compute the reciprocal condition number.
215: *
216: IF( AINVNM.NE.ZERO )
217: $ RCOND = ( ONE / AINVNM ) / ANORM
218: *
219: RETURN
220: *
221: * End of ZPTCON
222: *
223: END
CVSweb interface <joel.bertrand@systella.fr>