00001 SUBROUTINE DPBEQU( UPLO, N, KD, AB, LDAB, S, SCOND, AMAX, INFO ) 00002 * 00003 * -- LAPACK routine (version 3.2.2) -- 00004 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00005 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00006 * June 2010 00007 * 00008 * .. Scalar Arguments .. 00009 CHARACTER UPLO 00010 INTEGER INFO, KD, LDAB, N 00011 DOUBLE PRECISION AMAX, SCOND 00012 * .. 00013 * .. Array Arguments .. 00014 DOUBLE PRECISION AB( LDAB, * ), S( * ) 00015 * .. 00016 * 00017 * Purpose 00018 * ======= 00019 * 00020 * DPBEQU computes row and column scalings intended to equilibrate a 00021 * symmetric positive definite band matrix A and reduce its condition 00022 * number (with respect to the two-norm). S contains the scale factors, 00023 * S(i) = 1/sqrt(A(i,i)), chosen so that the scaled matrix B with 00024 * elements B(i,j) = S(i)*A(i,j)*S(j) has ones on the diagonal. This 00025 * choice of S puts the condition number of B within a factor N of the 00026 * smallest possible condition number over all possible diagonal 00027 * scalings. 00028 * 00029 * Arguments 00030 * ========= 00031 * 00032 * UPLO (input) CHARACTER*1 00033 * = 'U': Upper triangular of A is stored; 00034 * = 'L': Lower triangular of A is stored. 00035 * 00036 * N (input) INTEGER 00037 * The order of the matrix A. N >= 0. 00038 * 00039 * KD (input) INTEGER 00040 * The number of superdiagonals of the matrix A if UPLO = 'U', 00041 * or the number of subdiagonals if UPLO = 'L'. KD >= 0. 00042 * 00043 * AB (input) DOUBLE PRECISION array, dimension (LDAB,N) 00044 * The upper or lower triangle of the symmetric band matrix A, 00045 * stored in the first KD+1 rows of the array. The j-th column 00046 * of A is stored in the j-th column of the array AB as follows: 00047 * if UPLO = 'U', AB(kd+1+i-j,j) = A(i,j) for max(1,j-kd)<=i<=j; 00048 * if UPLO = 'L', AB(1+i-j,j) = A(i,j) for j<=i<=min(n,j+kd). 00049 * 00050 * LDAB (input) INTEGER 00051 * The leading dimension of the array A. LDAB >= KD+1. 00052 * 00053 * S (output) DOUBLE PRECISION array, dimension (N) 00054 * If INFO = 0, S contains the scale factors for A. 00055 * 00056 * SCOND (output) DOUBLE PRECISION 00057 * If INFO = 0, S contains the ratio of the smallest S(i) to 00058 * the largest S(i). If SCOND >= 0.1 and AMAX is neither too 00059 * large nor too small, it is not worth scaling by S. 00060 * 00061 * AMAX (output) DOUBLE PRECISION 00062 * Absolute value of largest matrix element. If AMAX is very 00063 * close to overflow or very close to underflow, the matrix 00064 * should be scaled. 00065 * 00066 * INFO (output) INTEGER 00067 * = 0: successful exit 00068 * < 0: if INFO = -i, the i-th argument had an illegal value. 00069 * > 0: if INFO = i, the i-th diagonal element is nonpositive. 00070 * 00071 * ===================================================================== 00072 * 00073 * .. Parameters .. 00074 DOUBLE PRECISION ZERO, ONE 00075 PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 ) 00076 * .. 00077 * .. Local Scalars .. 00078 LOGICAL UPPER 00079 INTEGER I, J 00080 DOUBLE PRECISION SMIN 00081 * .. 00082 * .. External Functions .. 00083 LOGICAL LSAME 00084 EXTERNAL LSAME 00085 * .. 00086 * .. External Subroutines .. 00087 EXTERNAL XERBLA 00088 * .. 00089 * .. Intrinsic Functions .. 00090 INTRINSIC MAX, MIN, SQRT 00091 * .. 00092 * .. Executable Statements .. 00093 * 00094 * Test the input parameters. 00095 * 00096 INFO = 0 00097 UPPER = LSAME( UPLO, 'U' ) 00098 IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN 00099 INFO = -1 00100 ELSE IF( N.LT.0 ) THEN 00101 INFO = -2 00102 ELSE IF( KD.LT.0 ) THEN 00103 INFO = -3 00104 ELSE IF( LDAB.LT.KD+1 ) THEN 00105 INFO = -5 00106 END IF 00107 IF( INFO.NE.0 ) THEN 00108 CALL XERBLA( 'DPBEQU', -INFO ) 00109 RETURN 00110 END IF 00111 * 00112 * Quick return if possible 00113 * 00114 IF( N.EQ.0 ) THEN 00115 SCOND = ONE 00116 AMAX = ZERO 00117 RETURN 00118 END IF 00119 * 00120 IF( UPPER ) THEN 00121 J = KD + 1 00122 ELSE 00123 J = 1 00124 END IF 00125 * 00126 * Initialize SMIN and AMAX. 00127 * 00128 S( 1 ) = AB( J, 1 ) 00129 SMIN = S( 1 ) 00130 AMAX = S( 1 ) 00131 * 00132 * Find the minimum and maximum diagonal elements. 00133 * 00134 DO 10 I = 2, N 00135 S( I ) = AB( J, I ) 00136 SMIN = MIN( SMIN, S( I ) ) 00137 AMAX = MAX( AMAX, S( I ) ) 00138 10 CONTINUE 00139 * 00140 IF( SMIN.LE.ZERO ) THEN 00141 * 00142 * Find the first non-positive diagonal element and return. 00143 * 00144 DO 20 I = 1, N 00145 IF( S( I ).LE.ZERO ) THEN 00146 INFO = I 00147 RETURN 00148 END IF 00149 20 CONTINUE 00150 ELSE 00151 * 00152 * Set the scale factors to the reciprocals 00153 * of the diagonal elements. 00154 * 00155 DO 30 I = 1, N 00156 S( I ) = ONE / SQRT( S( I ) ) 00157 30 CONTINUE 00158 * 00159 * Compute SCOND = min(S(I)) / max(S(I)) 00160 * 00161 SCOND = SQRT( SMIN ) / SQRT( AMAX ) 00162 END IF 00163 RETURN 00164 * 00165 * End of DPBEQU 00166 * 00167 END
1.7.2