LAPACK 3.3.0

cpoequb.f

Go to the documentation of this file.
00001       SUBROUTINE CPOEQUB( N, A, LDA, S, SCOND, AMAX, INFO )
00002 *
00003 *     -- LAPACK routine (version 3.2)                                 --
00004 *     -- Contributed by James Demmel, Deaglan Halligan, Yozo Hida and --
00005 *     -- Jason Riedy of Univ. of California Berkeley.                 --
00006 *     -- November 2008                                                --
00007 *
00008 *     -- LAPACK is a software package provided by Univ. of Tennessee, --
00009 *     -- Univ. of California Berkeley and NAG Ltd.                    --
00010 *
00011       IMPLICIT NONE
00012 *     ..
00013 *     .. Scalar Arguments ..
00014       INTEGER            INFO, LDA, N
00015       REAL               AMAX, SCOND
00016 *     ..
00017 *     .. Array Arguments ..
00018       COMPLEX            A( LDA, * )
00019       REAL               S( * )
00020 *     ..
00021 *
00022 *  Purpose
00023 *  =======
00024 *
00025 *  CPOEQUB computes row and column scalings intended to equilibrate a
00026 *  symmetric positive definite matrix A and reduce its condition number
00027 *  (with respect to the two-norm).  S contains the scale factors,
00028 *  S(i) = 1/sqrt(A(i,i)), chosen so that the scaled matrix B with
00029 *  elements B(i,j) = S(i)*A(i,j)*S(j) has ones on the diagonal.  This
00030 *  choice of S puts the condition number of B within a factor N of the
00031 *  smallest possible condition number over all possible diagonal
00032 *  scalings.
00033 *
00034 *  Arguments
00035 *  =========
00036 *
00037 *  N       (input) INTEGER
00038 *          The order of the matrix A.  N >= 0.
00039 *
00040 *  A       (input) COMPLEX array, dimension (LDA,N)
00041 *          The N-by-N symmetric positive definite matrix whose scaling
00042 *          factors are to be computed.  Only the diagonal elements of A
00043 *          are referenced.
00044 *
00045 *  LDA     (input) INTEGER
00046 *          The leading dimension of the array A.  LDA >= max(1,N).
00047 *
00048 *  S       (output) REAL array, dimension (N)
00049 *          If INFO = 0, S contains the scale factors for A.
00050 *
00051 *  SCOND   (output) REAL
00052 *          If INFO = 0, S contains the ratio of the smallest S(i) to
00053 *          the largest S(i).  If SCOND >= 0.1 and AMAX is neither too
00054 *          large nor too small, it is not worth scaling by S.
00055 *
00056 *  AMAX    (output) REAL
00057 *          Absolute value of largest matrix element.  If AMAX is very
00058 *          close to overflow or very close to underflow, the matrix
00059 *          should be scaled.
00060 *
00061 *  INFO    (output) INTEGER
00062 *          = 0:  successful exit
00063 *          < 0:  if INFO = -i, the i-th argument had an illegal value
00064 *          > 0:  if INFO = i, the i-th diagonal element is nonpositive.
00065 *
00066 *  =====================================================================
00067 *
00068 *     .. Parameters ..
00069       REAL               ZERO, ONE
00070       PARAMETER          ( ZERO = 0.0E+0, ONE = 1.0E+0 )
00071 *     ..
00072 *     .. Local Scalars ..
00073       INTEGER            I
00074       REAL               SMIN, BASE, TMP
00075       COMPLEX            ZDUM
00076 *     ..
00077 *     .. External Functions ..
00078       REAL               SLAMCH
00079       EXTERNAL           SLAMCH
00080 *     ..
00081 *     .. External Subroutines ..
00082       EXTERNAL           XERBLA
00083 *     ..
00084 *     .. Intrinsic Functions ..
00085       INTRINSIC          MAX, MIN, SQRT, LOG, INT, REAL, AIMAG
00086 *     ..
00087 *     .. Statement Functions ..
00088       REAL               CABS1
00089 *     ..
00090 *     .. Statement Function Definitions ..
00091       CABS1( ZDUM ) = ABS( REAL( ZDUM ) ) + ABS( AIMAG( ZDUM ) )
00092 *     ..
00093 *     .. Executable Statements ..
00094 *
00095 *     Test the input parameters.
00096 *
00097 *     Positive definite only performs 1 pass of equilibration.
00098 *
00099       INFO = 0
00100       IF( N.LT.0 ) THEN
00101          INFO = -1
00102       ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
00103          INFO = -3
00104       END IF
00105       IF( INFO.NE.0 ) THEN
00106          CALL XERBLA( 'CPOEQUB', -INFO )
00107          RETURN
00108       END IF
00109 *
00110 *     Quick return if possible.
00111 *
00112       IF( N.EQ.0 ) THEN
00113          SCOND = ONE
00114          AMAX = ZERO
00115          RETURN
00116       END IF
00117 
00118       BASE = SLAMCH( 'B' )
00119       TMP = -0.5 / LOG ( BASE )
00120 *
00121 *     Find the minimum and maximum diagonal elements.
00122 *
00123       S( 1 ) = A( 1, 1 )
00124       SMIN = S( 1 )
00125       AMAX = S( 1 )
00126       DO 10 I = 2, N
00127          S( I ) = A( I, I )
00128          SMIN = MIN( SMIN, S( I ) )
00129          AMAX = MAX( AMAX, S( I ) )
00130    10 CONTINUE
00131 *
00132       IF( SMIN.LE.ZERO ) THEN
00133 *
00134 *        Find the first non-positive diagonal element and return.
00135 *
00136          DO 20 I = 1, N
00137             IF( S( I ).LE.ZERO ) THEN
00138                INFO = I
00139                RETURN
00140             END IF
00141    20    CONTINUE
00142       ELSE
00143 *
00144 *        Set the scale factors to the reciprocals
00145 *        of the diagonal elements.
00146 *
00147          DO 30 I = 1, N
00148             S( I ) = BASE ** INT( TMP * LOG( S( I ) ) )
00149    30    CONTINUE
00150 *
00151 *        Compute SCOND = min(S(I)) / max(S(I)).
00152 *
00153          SCOND = SQRT( SMIN ) / SQRT( AMAX )
00154       END IF
00155 *
00156       RETURN
00157 *
00158 *     End of CPOEQUB
00159 *
00160       END
 All Files Functions