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