diff options
author | julie <julielangou@users.noreply.github.com> | 2008-12-16 17:06:58 +0000 |
---|---|---|
committer | julie <julielangou@users.noreply.github.com> | 2008-12-16 17:06:58 +0000 |
commit | ff981f106bde4ce6a74aa4f4a572c943f5a395b2 (patch) | |
tree | a386cad907bcaefd6893535c31d67ec9468e693e /SRC/cpoequb.f | |
parent | e58b61578b55644f6391f3333262b72c1dc88437 (diff) | |
download | lapack-ff981f106bde4ce6a74aa4f4a572c943f5a395b2.tar.gz lapack-ff981f106bde4ce6a74aa4f4a572c943f5a395b2.tar.bz2 lapack-ff981f106bde4ce6a74aa4f4a572c943f5a395b2.zip |
Diffstat (limited to 'SRC/cpoequb.f')
-rw-r--r-- | SRC/cpoequb.f | 160 |
1 files changed, 160 insertions, 0 deletions
diff --git a/SRC/cpoequb.f b/SRC/cpoequb.f new file mode 100644 index 00000000..88a87b71 --- /dev/null +++ b/SRC/cpoequb.f @@ -0,0 +1,160 @@ + SUBROUTINE CPOEQUB( N, A, LDA, S, SCOND, AMAX, INFO ) +* +* -- LAPACK routine (version 3.2) -- +* -- Contributed by James Demmel, Deaglan Halligan, Yozo Hida and -- +* -- Jason Riedy of Univ. of California Berkeley. -- +* -- November 2008 -- +* +* -- LAPACK is a software package provided by Univ. of Tennessee, -- +* -- Univ. of California Berkeley and NAG Ltd. -- +* + IMPLICIT NONE +* .. +* .. Scalar Arguments .. + INTEGER INFO, LDA, N + REAL AMAX, SCOND +* .. +* .. Array Arguments .. + COMPLEX A( LDA, * ) + REAL S( * ) +* .. +* +* Purpose +* ======= +* +* CPOEQUB computes row and column scalings intended to equilibrate a +* symmetric positive definite matrix A and reduce its condition number +* (with respect to the two-norm). S contains the scale factors, +* S(i) = 1/sqrt(A(i,i)), chosen so that the scaled matrix B with +* elements B(i,j) = S(i)*A(i,j)*S(j) has ones on the diagonal. This +* choice of S puts the condition number of B within a factor N of the +* smallest possible condition number over all possible diagonal +* scalings. +* +* Arguments +* ========= +* +* N (input) INTEGER +* The order of the matrix A. N >= 0. +* +* A (input) COMPLEX array, dimension (LDA,N) +* The N-by-N symmetric positive definite matrix whose scaling +* factors are to be computed. Only the diagonal elements of A +* are referenced. +* +* LDA (input) INTEGER +* The leading dimension of the array A. LDA >= max(1,N). +* +* S (output) REAL array, dimension (N) +* If INFO = 0, S contains the scale factors for A. +* +* SCOND (output) REAL +* If INFO = 0, S contains the ratio of the smallest S(i) to +* the largest S(i). If SCOND >= 0.1 and AMAX is neither too +* large nor too small, it is not worth scaling by S. +* +* AMAX (output) REAL +* Absolute value of largest matrix element. If AMAX is very +* close to overflow or very close to underflow, the matrix +* should be scaled. +* +* INFO (output) INTEGER +* = 0: successful exit +* < 0: if INFO = -i, the i-th argument had an illegal value +* > 0: if INFO = i, the i-th diagonal element is nonpositive. +* +* ===================================================================== +* +* .. Parameters .. + REAL ZERO, ONE + PARAMETER ( ZERO = 0.0E+0, ONE = 1.0E+0 ) +* .. +* .. Local Scalars .. + INTEGER I + REAL SMIN, BASE, TMP + COMPLEX ZDUM +* .. +* .. External Functions .. + REAL SLAMCH + EXTERNAL SLAMCH +* .. +* .. External Subroutines .. + EXTERNAL XERBLA +* .. +* .. Intrinsic Functions .. + INTRINSIC MAX, MIN, SQRT, LOG, INT, REAL, AIMAG +* .. +* .. Statement Functions .. + REAL CABS1 +* .. +* .. Statement Function Definitions .. + CABS1( ZDUM ) = ABS( REAL( ZDUM ) ) + ABS( AIMAG( ZDUM ) ) +* .. +* .. Executable Statements .. +* +* Test the input parameters. +* +* Positive definite only performs 1 pass of equilibration. +* + INFO = 0 + IF( N.LT.0 ) THEN + INFO = -1 + ELSE IF( LDA.LT.MAX( 1, N ) ) THEN + INFO = -3 + END IF + IF( INFO.NE.0 ) THEN + CALL XERBLA( 'CPOEQUB', -INFO ) + RETURN + END IF +* +* Quick return if possible. +* + IF( N.EQ.0 ) THEN + SCOND = ONE + AMAX = ZERO + RETURN + END IF + + BASE = SLAMCH( 'B' ) + TMP = -0.5 / LOG ( BASE ) +* +* Find the minimum and maximum diagonal elements. +* + S( 1 ) = A( 1, 1 ) + SMIN = S( 1 ) + AMAX = S( 1 ) + DO 10 I = 2, N + S( I ) = A( I, I ) + SMIN = MIN( SMIN, S( I ) ) + AMAX = MAX( AMAX, S( I ) ) + 10 CONTINUE +* + IF( SMIN.LE.ZERO ) THEN +* +* Find the first non-positive diagonal element and return. +* + DO 20 I = 1, N + IF( S( I ).LE.ZERO ) THEN + INFO = I + RETURN + END IF + 20 CONTINUE + ELSE +* +* Set the scale factors to the reciprocals +* of the diagonal elements. +* + DO 30 I = 1, N + S( I ) = BASE ** INT( TMP * LOG( S( I ) ) ) + 30 CONTINUE +* +* Compute SCOND = min(S(I)) / max(S(I)). +* + SCOND = SQRT( SMIN ) / SQRT( AMAX ) + END IF +* + RETURN +* +* End of CPOEQUB +* + END |