1*> \brief <b> CHBEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for OTHER matrices</b>
2*
3*  =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6*            http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download CHBEV + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/chbev.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/chbev.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/chbev.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18*  Definition:
19*  ===========
20*
21*       SUBROUTINE CHBEV( JOBZ, UPLO, N, KD, AB, LDAB, W, Z, LDZ, WORK,
22*                         RWORK, INFO )
23*
24*       .. Scalar Arguments ..
25*       CHARACTER          JOBZ, UPLO
26*       INTEGER            INFO, KD, LDAB, LDZ, N
27*       ..
28*       .. Array Arguments ..
29*       REAL               RWORK( * ), W( * )
30*       COMPLEX            AB( LDAB, * ), WORK( * ), Z( LDZ, * )
31*       ..
32*
33*
34*> \par Purpose:
35*  =============
36*>
37*> \verbatim
38*>
39*> CHBEV computes all the eigenvalues and, optionally, eigenvectors of
40*> a complex Hermitian band matrix A.
41*> \endverbatim
42*
43*  Arguments:
44*  ==========
45*
46*> \param[in] JOBZ
47*> \verbatim
48*>          JOBZ is CHARACTER*1
49*>          = 'N':  Compute eigenvalues only;
50*>          = 'V':  Compute eigenvalues and eigenvectors.
51*> \endverbatim
52*>
53*> \param[in] UPLO
54*> \verbatim
55*>          UPLO is CHARACTER*1
56*>          = 'U':  Upper triangle of A is stored;
57*>          = 'L':  Lower triangle of A is stored.
58*> \endverbatim
59*>
60*> \param[in] N
61*> \verbatim
62*>          N is INTEGER
63*>          The order of the matrix A.  N >= 0.
64*> \endverbatim
65*>
66*> \param[in] KD
67*> \verbatim
68*>          KD is INTEGER
69*>          The number of superdiagonals of the matrix A if UPLO = 'U',
70*>          or the number of subdiagonals if UPLO = 'L'.  KD >= 0.
71*> \endverbatim
72*>
73*> \param[in,out] AB
74*> \verbatim
75*>          AB is COMPLEX array, dimension (LDAB, N)
76*>          On entry, the upper or lower triangle of the Hermitian band
77*>          matrix A, stored in the first KD+1 rows of the array.  The
78*>          j-th column of A is stored in the j-th column of the array AB
79*>          as follows:
80*>          if UPLO = 'U', AB(kd+1+i-j,j) = A(i,j) for max(1,j-kd)<=i<=j;
81*>          if UPLO = 'L', AB(1+i-j,j)    = A(i,j) for j<=i<=min(n,j+kd).
82*>
83*>          On exit, AB is overwritten by values generated during the
84*>          reduction to tridiagonal form.  If UPLO = 'U', the first
85*>          superdiagonal and the diagonal of the tridiagonal matrix T
86*>          are returned in rows KD and KD+1 of AB, and if UPLO = 'L',
87*>          the diagonal and first subdiagonal of T are returned in the
88*>          first two rows of AB.
89*> \endverbatim
90*>
91*> \param[in] LDAB
92*> \verbatim
93*>          LDAB is INTEGER
94*>          The leading dimension of the array AB.  LDAB >= KD + 1.
95*> \endverbatim
96*>
97*> \param[out] W
98*> \verbatim
99*>          W is REAL array, dimension (N)
100*>          If INFO = 0, the eigenvalues in ascending order.
101*> \endverbatim
102*>
103*> \param[out] Z
104*> \verbatim
105*>          Z is COMPLEX array, dimension (LDZ, N)
106*>          If JOBZ = 'V', then if INFO = 0, Z contains the orthonormal
107*>          eigenvectors of the matrix A, with the i-th column of Z
108*>          holding the eigenvector associated with W(i).
109*>          If JOBZ = 'N', then Z is not referenced.
110*> \endverbatim
111*>
112*> \param[in] LDZ
113*> \verbatim
114*>          LDZ is INTEGER
115*>          The leading dimension of the array Z.  LDZ >= 1, and if
116*>          JOBZ = 'V', LDZ >= max(1,N).
117*> \endverbatim
118*>
119*> \param[out] WORK
120*> \verbatim
121*>          WORK is COMPLEX array, dimension (N)
122*> \endverbatim
123*>
124*> \param[out] RWORK
125*> \verbatim
126*>          RWORK is REAL array, dimension (max(1,3*N-2))
127*> \endverbatim
128*>
129*> \param[out] INFO
130*> \verbatim
131*>          INFO is INTEGER
132*>          = 0:  successful exit.
133*>          < 0:  if INFO = -i, the i-th argument had an illegal value.
134*>          > 0:  if INFO = i, the algorithm failed to converge; i
135*>                off-diagonal elements of an intermediate tridiagonal
136*>                form did not converge to zero.
137*> \endverbatim
138*
139*  Authors:
140*  ========
141*
142*> \author Univ. of Tennessee
143*> \author Univ. of California Berkeley
144*> \author Univ. of Colorado Denver
145*> \author NAG Ltd.
146*
147*> \ingroup complexOTHEReigen
148*
149*  =====================================================================
150      SUBROUTINE CHBEV( JOBZ, UPLO, N, KD, AB, LDAB, W, Z, LDZ, WORK,
151     $                  RWORK, INFO )
152*
153*  -- LAPACK driver routine --
154*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
155*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
156*
157*     .. Scalar Arguments ..
158      CHARACTER          JOBZ, UPLO
159      INTEGER            INFO, KD, LDAB, LDZ, N
160*     ..
161*     .. Array Arguments ..
162      REAL               RWORK( * ), W( * )
163      COMPLEX            AB( LDAB, * ), WORK( * ), Z( LDZ, * )
164*     ..
165*
166*  =====================================================================
167*
168*     .. Parameters ..
169      REAL               ZERO, ONE
170      PARAMETER          ( ZERO = 0.0E0, ONE = 1.0E0 )
171*     ..
172*     .. Local Scalars ..
173      LOGICAL            LOWER, WANTZ
174      INTEGER            IINFO, IMAX, INDE, INDRWK, ISCALE
175      REAL               ANRM, BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA,
176     $                   SMLNUM
177*     ..
178*     .. External Functions ..
179      LOGICAL            LSAME
180      REAL               CLANHB, SLAMCH
181      EXTERNAL           LSAME, CLANHB, SLAMCH
182*     ..
183*     .. External Subroutines ..
184      EXTERNAL           CHBTRD, CLASCL, CSTEQR, SSCAL, SSTERF, XERBLA
185*     ..
186*     .. Intrinsic Functions ..
187      INTRINSIC          SQRT
188*     ..
189*     .. Executable Statements ..
190*
191*     Test the input parameters.
192*
193      WANTZ = LSAME( JOBZ, 'V' )
194      LOWER = LSAME( UPLO, 'L' )
195*
196      INFO = 0
197      IF( .NOT.( WANTZ .OR. LSAME( JOBZ, 'N' ) ) ) THEN
198         INFO = -1
199      ELSE IF( .NOT.( LOWER .OR. LSAME( UPLO, 'U' ) ) ) THEN
200         INFO = -2
201      ELSE IF( N.LT.0 ) THEN
202         INFO = -3
203      ELSE IF( KD.LT.0 ) THEN
204         INFO = -4
205      ELSE IF( LDAB.LT.KD+1 ) THEN
206         INFO = -6
207      ELSE IF( LDZ.LT.1 .OR. ( WANTZ .AND. LDZ.LT.N ) ) THEN
208         INFO = -9
209      END IF
210*
211      IF( INFO.NE.0 ) THEN
212         CALL XERBLA( 'CHBEV ', -INFO )
213         RETURN
214      END IF
215*
216*     Quick return if possible
217*
218      IF( N.EQ.0 )
219     $   RETURN
220*
221      IF( N.EQ.1 ) THEN
222         IF( LOWER ) THEN
223            W( 1 ) = REAL( AB( 1, 1 ) )
224         ELSE
225            W( 1 ) = REAL( AB( KD+1, 1 ) )
226         END IF
227         IF( WANTZ )
228     $      Z( 1, 1 ) = ONE
229         RETURN
230      END IF
231*
232*     Get machine constants.
233*
234      SAFMIN = SLAMCH( 'Safe minimum' )
235      EPS = SLAMCH( 'Precision' )
236      SMLNUM = SAFMIN / EPS
237      BIGNUM = ONE / SMLNUM
238      RMIN = SQRT( SMLNUM )
239      RMAX = SQRT( BIGNUM )
240*
241*     Scale matrix to allowable range, if necessary.
242*
243      ANRM = CLANHB( 'M', UPLO, N, KD, AB, LDAB, RWORK )
244      ISCALE = 0
245      IF( ANRM.GT.ZERO .AND. ANRM.LT.RMIN ) THEN
246         ISCALE = 1
247         SIGMA = RMIN / ANRM
248      ELSE IF( ANRM.GT.RMAX ) THEN
249         ISCALE = 1
250         SIGMA = RMAX / ANRM
251      END IF
252      IF( ISCALE.EQ.1 ) THEN
253         IF( LOWER ) THEN
254            CALL CLASCL( 'B', KD, KD, ONE, SIGMA, N, N, AB, LDAB, INFO )
255         ELSE
256            CALL CLASCL( 'Q', KD, KD, ONE, SIGMA, N, N, AB, LDAB, INFO )
257         END IF
258      END IF
259*
260*     Call CHBTRD to reduce Hermitian band matrix to tridiagonal form.
261*
262      INDE = 1
263      CALL CHBTRD( JOBZ, UPLO, N, KD, AB, LDAB, W, RWORK( INDE ), Z,
264     $             LDZ, WORK, IINFO )
265*
266*     For eigenvalues only, call SSTERF.  For eigenvectors, call CSTEQR.
267*
268      IF( .NOT.WANTZ ) THEN
269         CALL SSTERF( N, W, RWORK( INDE ), INFO )
270      ELSE
271         INDRWK = INDE + N
272         CALL CSTEQR( JOBZ, N, W, RWORK( INDE ), Z, LDZ,
273     $                RWORK( INDRWK ), INFO )
274      END IF
275*
276*     If matrix was scaled, then rescale eigenvalues appropriately.
277*
278      IF( ISCALE.EQ.1 ) THEN
279         IF( INFO.EQ.0 ) THEN
280            IMAX = N
281         ELSE
282            IMAX = INFO - 1
283         END IF
284         CALL SSCAL( IMAX, ONE / SIGMA, W, 1 )
285      END IF
286*
287      RETURN
288*
289*     End of CHBEV
290*
291      END
292