1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152
|
SUBROUTINE CLAQHB( UPLO, N, KD, AB, LDAB, S, SCOND, AMAX, EQUED )
*
* -- LAPACK auxiliary routine (version 2.0) --
* Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd.,
* Courant Institute, Argonne National Lab, and Rice University
* September 30, 1994
*
* .. Scalar Arguments ..
CHARACTER EQUED, UPLO
INTEGER KD, LDAB, N
REAL AMAX, SCOND
* ..
* .. Array Arguments ..
REAL S( * )
COMPLEX AB( LDAB, * )
* ..
*
* Purpose
* =======
*
* CLAQHB equilibrates a symmetric band matrix A using the scaling
* factors in the vector S.
*
* Arguments
* =========
*
* UPLO (input) CHARACTER*1
* Specifies whether the upper or lower triangular part of the
* symmetric matrix A is stored.
* = 'U': Upper triangular
* = 'L': Lower triangular
*
* N (input) INTEGER
* The order of the matrix A. N >= 0.
*
* KD (input) INTEGER
* The number of super-diagonals of the matrix A if UPLO = 'U',
* or the number of sub-diagonals if UPLO = 'L'. KD >= 0.
*
* AB (input/output) COMPLEX array, dimension (LDAB,N)
* On entry, the upper or lower triangle of the symmetric band
* matrix A, stored in the first KD+1 rows of the array. The
* j-th column of A is stored in the j-th column of the array AB
* as follows:
* if UPLO = 'U', AB(kd+1+i-j,j) = A(i,j) for max(1,j-kd)<=i<=j;
* if UPLO = 'L', AB(1+i-j,j) = A(i,j) for j<=i<=min(n,j+kd).
*
* On exit, if INFO = 0, the triangular factor U or L from the
* Cholesky factorization A = U'*U or A = L*L' of the band
* matrix A, in the same storage format as A.
*
* LDAB (input) INTEGER
* The leading dimension of the array AB. LDAB >= KD+1.
*
* S (output) REAL array, dimension (N)
* The scale factors for A.
*
* SCOND (input) REAL
* Ratio of the smallest S(i) to the largest S(i).
*
* AMAX (input) REAL
* Absolute value of largest matrix entry.
*
* EQUED (output) CHARACTER*1
* Specifies whether or not equilibration was done.
* = 'N': No equilibration.
* = 'Y': Equilibration was done, i.e., A has been replaced by
* diag(S) * A * diag(S).
*
* Internal Parameters
* ===================
*
* THRESH is a threshold value used to decide if scaling should be done
* based on the ratio of the scaling factors. If SCOND < THRESH,
* scaling is done.
*
* LARGE and SMALL are threshold values used to decide if scaling should
* be done based on the absolute size of the largest matrix element.
* If AMAX > LARGE or AMAX < SMALL, scaling is done.
*
* =====================================================================
*
* .. Parameters ..
REAL ONE, THRESH
PARAMETER ( ONE = 1.0E+0, THRESH = 0.1E+0 )
* ..
* .. Local Scalars ..
INTEGER I, J
REAL CJ, LARGE, SMALL
* ..
* .. External Functions ..
LOGICAL LSAME
REAL SLAMCH
EXTERNAL LSAME, SLAMCH
* ..
* .. Intrinsic Functions ..
INTRINSIC MAX, MIN, REAL
* ..
* .. Executable Statements ..
*
* Quick return if possible
*
IF( N.LE.0 ) THEN
EQUED = 'N'
RETURN
END IF
*
* Initialize LARGE and SMALL.
*
SMALL = SLAMCH( 'Safe minimum' ) / SLAMCH( 'Precision' )
LARGE = ONE / SMALL
*
IF( SCOND.GE.THRESH .AND. AMAX.GE.SMALL .AND. AMAX.LE.LARGE ) THEN
*
* No equilibration
*
EQUED = 'N'
ELSE
*
* Replace A by diag(S) * A * diag(S).
*
IF( LSAME( UPLO, 'U' ) ) THEN
*
* Upper triangle of A is stored in band format.
*
DO 20 J = 1, N
CJ = S( J )
DO 10 I = MAX( 1, J-KD ), J - 1
AB( KD+1+I-J, J ) = CJ*S( I )*AB( KD+1+I-J, J )
10 CONTINUE
AB( KD+1, J ) = CJ*CJ*REAL( AB( KD+1, J ) )
20 CONTINUE
ELSE
*
* Lower triangle of A is stored.
*
DO 40 J = 1, N
CJ = S( J )
AB( 1, J ) = CJ*CJ*REAL( AB( 1, J ) )
DO 30 I = J + 1, MIN( N, J+KD )
AB( 1+I-J, J ) = CJ*S( I )*AB( 1+I-J, J )
30 CONTINUE
40 CONTINUE
END IF
EQUED = 'Y'
END IF
*
RETURN
*
* End of CLAQHB
*
END
|