Annotation of rpl/lapack/lapack/dpbcon.f, revision 1.9
1.9 ! bertrand 1: *> \brief \b DPBCON
! 2: *
! 3: * =========== DOCUMENTATION ===========
! 4: *
! 5: * Online html documentation available at
! 6: * http://www.netlib.org/lapack/explore-html/
! 7: *
! 8: *> \htmlonly
! 9: *> Download DPBCON + dependencies
! 10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dpbcon.f">
! 11: *> [TGZ]</a>
! 12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dpbcon.f">
! 13: *> [ZIP]</a>
! 14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dpbcon.f">
! 15: *> [TXT]</a>
! 16: *> \endhtmlonly
! 17: *
! 18: * Definition:
! 19: * ===========
! 20: *
! 21: * SUBROUTINE DPBCON( UPLO, N, KD, AB, LDAB, ANORM, RCOND, WORK,
! 22: * IWORK, INFO )
! 23: *
! 24: * .. Scalar Arguments ..
! 25: * CHARACTER UPLO
! 26: * INTEGER INFO, KD, LDAB, N
! 27: * DOUBLE PRECISION ANORM, RCOND
! 28: * ..
! 29: * .. Array Arguments ..
! 30: * INTEGER IWORK( * )
! 31: * DOUBLE PRECISION AB( LDAB, * ), WORK( * )
! 32: * ..
! 33: *
! 34: *
! 35: *> \par Purpose:
! 36: * =============
! 37: *>
! 38: *> \verbatim
! 39: *>
! 40: *> DPBCON estimates the reciprocal of the condition number (in the
! 41: *> 1-norm) of a real symmetric positive definite band matrix using the
! 42: *> Cholesky factorization A = U**T*U or A = L*L**T computed by DPBTRF.
! 43: *>
! 44: *> An estimate is obtained for norm(inv(A)), and the reciprocal of the
! 45: *> condition number is computed as RCOND = 1 / (ANORM * norm(inv(A))).
! 46: *> \endverbatim
! 47: *
! 48: * Arguments:
! 49: * ==========
! 50: *
! 51: *> \param[in] UPLO
! 52: *> \verbatim
! 53: *> UPLO is CHARACTER*1
! 54: *> = 'U': Upper triangular factor stored in AB;
! 55: *> = 'L': Lower triangular factor stored in AB.
! 56: *> \endverbatim
! 57: *>
! 58: *> \param[in] N
! 59: *> \verbatim
! 60: *> N is INTEGER
! 61: *> The order of the matrix A. N >= 0.
! 62: *> \endverbatim
! 63: *>
! 64: *> \param[in] KD
! 65: *> \verbatim
! 66: *> KD is INTEGER
! 67: *> The number of superdiagonals of the matrix A if UPLO = 'U',
! 68: *> or the number of subdiagonals if UPLO = 'L'. KD >= 0.
! 69: *> \endverbatim
! 70: *>
! 71: *> \param[in] AB
! 72: *> \verbatim
! 73: *> AB is DOUBLE PRECISION array, dimension (LDAB,N)
! 74: *> The triangular factor U or L from the Cholesky factorization
! 75: *> A = U**T*U or A = L*L**T of the band matrix A, stored in the
! 76: *> first KD+1 rows of the array. The j-th column of U or L is
! 77: *> stored in the j-th column of the array AB as follows:
! 78: *> if UPLO ='U', AB(kd+1+i-j,j) = U(i,j) for max(1,j-kd)<=i<=j;
! 79: *> if UPLO ='L', AB(1+i-j,j) = L(i,j) for j<=i<=min(n,j+kd).
! 80: *> \endverbatim
! 81: *>
! 82: *> \param[in] LDAB
! 83: *> \verbatim
! 84: *> LDAB is INTEGER
! 85: *> The leading dimension of the array AB. LDAB >= KD+1.
! 86: *> \endverbatim
! 87: *>
! 88: *> \param[in] ANORM
! 89: *> \verbatim
! 90: *> ANORM is DOUBLE PRECISION
! 91: *> The 1-norm (or infinity-norm) of the symmetric band matrix A.
! 92: *> \endverbatim
! 93: *>
! 94: *> \param[out] RCOND
! 95: *> \verbatim
! 96: *> RCOND is DOUBLE PRECISION
! 97: *> The reciprocal of the condition number of the matrix A,
! 98: *> computed as RCOND = 1/(ANORM * AINVNM), where AINVNM is an
! 99: *> estimate of the 1-norm of inv(A) computed in this routine.
! 100: *> \endverbatim
! 101: *>
! 102: *> \param[out] WORK
! 103: *> \verbatim
! 104: *> WORK is DOUBLE PRECISION array, dimension (3*N)
! 105: *> \endverbatim
! 106: *>
! 107: *> \param[out] IWORK
! 108: *> \verbatim
! 109: *> IWORK is INTEGER array, dimension (N)
! 110: *> \endverbatim
! 111: *>
! 112: *> \param[out] INFO
! 113: *> \verbatim
! 114: *> INFO is INTEGER
! 115: *> = 0: successful exit
! 116: *> < 0: if INFO = -i, the i-th argument had an illegal value
! 117: *> \endverbatim
! 118: *
! 119: * Authors:
! 120: * ========
! 121: *
! 122: *> \author Univ. of Tennessee
! 123: *> \author Univ. of California Berkeley
! 124: *> \author Univ. of Colorado Denver
! 125: *> \author NAG Ltd.
! 126: *
! 127: *> \date November 2011
! 128: *
! 129: *> \ingroup doubleOTHERcomputational
! 130: *
! 131: * =====================================================================
1.1 bertrand 132: SUBROUTINE DPBCON( UPLO, N, KD, AB, LDAB, ANORM, RCOND, WORK,
133: $ IWORK, INFO )
134: *
1.9 ! bertrand 135: * -- LAPACK computational routine (version 3.4.0) --
1.1 bertrand 136: * -- LAPACK is a software package provided by Univ. of Tennessee, --
137: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
1.9 ! bertrand 138: * November 2011
1.1 bertrand 139: *
140: * .. Scalar Arguments ..
141: CHARACTER UPLO
142: INTEGER INFO, KD, LDAB, N
143: DOUBLE PRECISION ANORM, RCOND
144: * ..
145: * .. Array Arguments ..
146: INTEGER IWORK( * )
147: DOUBLE PRECISION AB( LDAB, * ), WORK( * )
148: * ..
149: *
150: * =====================================================================
151: *
152: * .. Parameters ..
153: DOUBLE PRECISION ONE, ZERO
154: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
155: * ..
156: * .. Local Scalars ..
157: LOGICAL UPPER
158: CHARACTER NORMIN
159: INTEGER IX, KASE
160: DOUBLE PRECISION AINVNM, SCALE, SCALEL, SCALEU, SMLNUM
161: * ..
162: * .. Local Arrays ..
163: INTEGER ISAVE( 3 )
164: * ..
165: * .. External Functions ..
166: LOGICAL LSAME
167: INTEGER IDAMAX
168: DOUBLE PRECISION DLAMCH
169: EXTERNAL LSAME, IDAMAX, DLAMCH
170: * ..
171: * .. External Subroutines ..
172: EXTERNAL DLACN2, DLATBS, DRSCL, XERBLA
173: * ..
174: * .. Intrinsic Functions ..
175: INTRINSIC ABS
176: * ..
177: * .. Executable Statements ..
178: *
179: * Test the input parameters.
180: *
181: INFO = 0
182: UPPER = LSAME( UPLO, 'U' )
183: IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
184: INFO = -1
185: ELSE IF( N.LT.0 ) THEN
186: INFO = -2
187: ELSE IF( KD.LT.0 ) THEN
188: INFO = -3
189: ELSE IF( LDAB.LT.KD+1 ) THEN
190: INFO = -5
191: ELSE IF( ANORM.LT.ZERO ) THEN
192: INFO = -6
193: END IF
194: IF( INFO.NE.0 ) THEN
195: CALL XERBLA( 'DPBCON', -INFO )
196: RETURN
197: END IF
198: *
199: * Quick return if possible
200: *
201: RCOND = ZERO
202: IF( N.EQ.0 ) THEN
203: RCOND = ONE
204: RETURN
205: ELSE IF( ANORM.EQ.ZERO ) THEN
206: RETURN
207: END IF
208: *
209: SMLNUM = DLAMCH( 'Safe minimum' )
210: *
211: * Estimate the 1-norm of the inverse.
212: *
213: KASE = 0
214: NORMIN = 'N'
215: 10 CONTINUE
216: CALL DLACN2( N, WORK( N+1 ), WORK, IWORK, AINVNM, KASE, ISAVE )
217: IF( KASE.NE.0 ) THEN
218: IF( UPPER ) THEN
219: *
1.8 bertrand 220: * Multiply by inv(U**T).
1.1 bertrand 221: *
222: CALL DLATBS( 'Upper', 'Transpose', 'Non-unit', NORMIN, N,
223: $ KD, AB, LDAB, WORK, SCALEL, WORK( 2*N+1 ),
224: $ INFO )
225: NORMIN = 'Y'
226: *
227: * Multiply by inv(U).
228: *
229: CALL DLATBS( 'Upper', 'No transpose', 'Non-unit', NORMIN, N,
230: $ KD, AB, LDAB, WORK, SCALEU, WORK( 2*N+1 ),
231: $ INFO )
232: ELSE
233: *
234: * Multiply by inv(L).
235: *
236: CALL DLATBS( 'Lower', 'No transpose', 'Non-unit', NORMIN, N,
237: $ KD, AB, LDAB, WORK, SCALEL, WORK( 2*N+1 ),
238: $ INFO )
239: NORMIN = 'Y'
240: *
1.8 bertrand 241: * Multiply by inv(L**T).
1.1 bertrand 242: *
243: CALL DLATBS( 'Lower', 'Transpose', 'Non-unit', NORMIN, N,
244: $ KD, AB, LDAB, WORK, SCALEU, WORK( 2*N+1 ),
245: $ INFO )
246: END IF
247: *
248: * Multiply by 1/SCALE if doing so will not cause overflow.
249: *
250: SCALE = SCALEL*SCALEU
251: IF( SCALE.NE.ONE ) THEN
252: IX = IDAMAX( N, WORK, 1 )
253: IF( SCALE.LT.ABS( WORK( IX ) )*SMLNUM .OR. SCALE.EQ.ZERO )
254: $ GO TO 20
255: CALL DRSCL( N, SCALE, WORK, 1 )
256: END IF
257: GO TO 10
258: END IF
259: *
260: * Compute the estimate of the reciprocal condition number.
261: *
262: IF( AINVNM.NE.ZERO )
263: $ RCOND = ( ONE / AINVNM ) / ANORM
264: *
265: 20 CONTINUE
266: *
267: RETURN
268: *
269: * End of DPBCON
270: *
271: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to dpbcon.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack