00001 SUBROUTINE ZLAQHB( UPLO, N, KD, AB, LDAB, S, SCOND, AMAX, EQUED ) 00002 * 00003 * -- LAPACK auxiliary routine (version 3.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 * November 2006 00007 * 00008 * .. Scalar Arguments .. 00009 CHARACTER EQUED, UPLO 00010 INTEGER KD, LDAB, N 00011 DOUBLE PRECISION AMAX, SCOND 00012 * .. 00013 * .. Array Arguments .. 00014 DOUBLE PRECISION S( * ) 00015 COMPLEX*16 AB( LDAB, * ) 00016 * .. 00017 * 00018 * Purpose 00019 * ======= 00020 * 00021 * ZLAQHB equilibrates a symmetric band matrix A using the scaling 00022 * factors in the vector S. 00023 * 00024 * Arguments 00025 * ========= 00026 * 00027 * UPLO (input) CHARACTER*1 00028 * Specifies whether the upper or lower triangular part of the 00029 * symmetric matrix A is stored. 00030 * = 'U': Upper triangular 00031 * = 'L': Lower triangular 00032 * 00033 * N (input) INTEGER 00034 * The order of the matrix A. N >= 0. 00035 * 00036 * KD (input) INTEGER 00037 * The number of super-diagonals of the matrix A if UPLO = 'U', 00038 * or the number of sub-diagonals if UPLO = 'L'. KD >= 0. 00039 * 00040 * AB (input/output) COMPLEX*16 array, dimension (LDAB,N) 00041 * On entry, the upper or lower triangle of the symmetric band 00042 * matrix A, stored in the first KD+1 rows of the array. The 00043 * j-th column of A is stored in the j-th column of the array AB 00044 * as follows: 00045 * if UPLO = 'U', AB(kd+1+i-j,j) = A(i,j) for max(1,j-kd)<=i<=j; 00046 * if UPLO = 'L', AB(1+i-j,j) = A(i,j) for j<=i<=min(n,j+kd). 00047 * 00048 * On exit, if INFO = 0, the triangular factor U or L from the 00049 * Cholesky factorization A = U'*U or A = L*L' of the band 00050 * matrix A, in the same storage format as A. 00051 * 00052 * LDAB (input) INTEGER 00053 * The leading dimension of the array AB. LDAB >= KD+1. 00054 * 00055 * S (output) DOUBLE PRECISION array, dimension (N) 00056 * The scale factors for A. 00057 * 00058 * SCOND (input) DOUBLE PRECISION 00059 * Ratio of the smallest S(i) to the largest S(i). 00060 * 00061 * AMAX (input) DOUBLE PRECISION 00062 * Absolute value of largest matrix entry. 00063 * 00064 * EQUED (output) CHARACTER*1 00065 * Specifies whether or not equilibration was done. 00066 * = 'N': No equilibration. 00067 * = 'Y': Equilibration was done, i.e., A has been replaced by 00068 * diag(S) * A * diag(S). 00069 * 00070 * Internal Parameters 00071 * =================== 00072 * 00073 * THRESH is a threshold value used to decide if scaling should be done 00074 * based on the ratio of the scaling factors. If SCOND < THRESH, 00075 * scaling is done. 00076 * 00077 * LARGE and SMALL are threshold values used to decide if scaling should 00078 * be done based on the absolute size of the largest matrix element. 00079 * If AMAX > LARGE or AMAX < SMALL, scaling is done. 00080 * 00081 * ===================================================================== 00082 * 00083 * .. Parameters .. 00084 DOUBLE PRECISION ONE, THRESH 00085 PARAMETER ( ONE = 1.0D+0, THRESH = 0.1D+0 ) 00086 * .. 00087 * .. Local Scalars .. 00088 INTEGER I, J 00089 DOUBLE PRECISION CJ, LARGE, SMALL 00090 * .. 00091 * .. External Functions .. 00092 LOGICAL LSAME 00093 DOUBLE PRECISION DLAMCH 00094 EXTERNAL LSAME, DLAMCH 00095 * .. 00096 * .. Intrinsic Functions .. 00097 INTRINSIC DBLE, MAX, MIN 00098 * .. 00099 * .. Executable Statements .. 00100 * 00101 * Quick return if possible 00102 * 00103 IF( N.LE.0 ) THEN 00104 EQUED = 'N' 00105 RETURN 00106 END IF 00107 * 00108 * Initialize LARGE and SMALL. 00109 * 00110 SMALL = DLAMCH( 'Safe minimum' ) / DLAMCH( 'Precision' ) 00111 LARGE = ONE / SMALL 00112 * 00113 IF( SCOND.GE.THRESH .AND. AMAX.GE.SMALL .AND. AMAX.LE.LARGE ) THEN 00114 * 00115 * No equilibration 00116 * 00117 EQUED = 'N' 00118 ELSE 00119 * 00120 * Replace A by diag(S) * A * diag(S). 00121 * 00122 IF( LSAME( UPLO, 'U' ) ) THEN 00123 * 00124 * Upper triangle of A is stored in band format. 00125 * 00126 DO 20 J = 1, N 00127 CJ = S( J ) 00128 DO 10 I = MAX( 1, J-KD ), J - 1 00129 AB( KD+1+I-J, J ) = CJ*S( I )*AB( KD+1+I-J, J ) 00130 10 CONTINUE 00131 AB( KD+1, J ) = CJ*CJ*DBLE( AB( KD+1, J ) ) 00132 20 CONTINUE 00133 ELSE 00134 * 00135 * Lower triangle of A is stored. 00136 * 00137 DO 40 J = 1, N 00138 CJ = S( J ) 00139 AB( 1, J ) = CJ*CJ*DBLE( AB( 1, J ) ) 00140 DO 30 I = J + 1, MIN( N, J+KD ) 00141 AB( 1+I-J, J ) = CJ*S( I )*AB( 1+I-J, J ) 00142 30 CONTINUE 00143 40 CONTINUE 00144 END IF 00145 EQUED = 'Y' 00146 END IF 00147 * 00148 RETURN 00149 * 00150 * End of ZLAQHB 00151 * 00152 END
1.7.2