LAPACK  3.4.2
LAPACK: Linear Algebra PACKage
 All Files Functions Groups
chpev.f
Go to the documentation of this file.
1 *> \brief <b> CHPEV 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 CHPEV + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/chpev.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/chpev.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/chpev.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE CHPEV( JOBZ, UPLO, N, AP, W, Z, LDZ, WORK, RWORK,
22 * INFO )
23 *
24 * .. Scalar Arguments ..
25 * CHARACTER JOBZ, UPLO
26 * INTEGER INFO, LDZ, N
27 * ..
28 * .. Array Arguments ..
29 * REAL RWORK( * ), W( * )
30 * COMPLEX AP( * ), WORK( * ), Z( LDZ, * )
31 * ..
32 *
33 *
34 *> \par Purpose:
35 * =============
36 *>
37 *> \verbatim
38 *>
39 *> CHPEV computes all the eigenvalues and, optionally, eigenvectors of a
40 *> complex Hermitian matrix in packed storage.
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,out] AP
67 *> \verbatim
68 *> AP is COMPLEX array, dimension (N*(N+1)/2)
69 *> On entry, the upper or lower triangle of the Hermitian matrix
70 *> A, packed columnwise in a linear array. The j-th column of A
71 *> is stored in the array AP as follows:
72 *> if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j;
73 *> if UPLO = 'L', AP(i + (j-1)*(2*n-j)/2) = A(i,j) for j<=i<=n.
74 *>
75 *> On exit, AP is overwritten by values generated during the
76 *> reduction to tridiagonal form. If UPLO = 'U', the diagonal
77 *> and first superdiagonal of the tridiagonal matrix T overwrite
78 *> the corresponding elements of A, and if UPLO = 'L', the
79 *> diagonal and first subdiagonal of T overwrite the
80 *> corresponding elements of A.
81 *> \endverbatim
82 *>
83 *> \param[out] W
84 *> \verbatim
85 *> W is REAL array, dimension (N)
86 *> If INFO = 0, the eigenvalues in ascending order.
87 *> \endverbatim
88 *>
89 *> \param[out] Z
90 *> \verbatim
91 *> Z is COMPLEX array, dimension (LDZ, N)
92 *> If JOBZ = 'V', then if INFO = 0, Z contains the orthonormal
93 *> eigenvectors of the matrix A, with the i-th column of Z
94 *> holding the eigenvector associated with W(i).
95 *> If JOBZ = 'N', then Z is not referenced.
96 *> \endverbatim
97 *>
98 *> \param[in] LDZ
99 *> \verbatim
100 *> LDZ is INTEGER
101 *> The leading dimension of the array Z. LDZ >= 1, and if
102 *> JOBZ = 'V', LDZ >= max(1,N).
103 *> \endverbatim
104 *>
105 *> \param[out] WORK
106 *> \verbatim
107 *> WORK is COMPLEX array, dimension (max(1, 2*N-1))
108 *> \endverbatim
109 *>
110 *> \param[out] RWORK
111 *> \verbatim
112 *> RWORK is REAL array, dimension (max(1, 3*N-2))
113 *> \endverbatim
114 *>
115 *> \param[out] INFO
116 *> \verbatim
117 *> INFO is INTEGER
118 *> = 0: successful exit.
119 *> < 0: if INFO = -i, the i-th argument had an illegal value.
120 *> > 0: if INFO = i, the algorithm failed to converge; i
121 *> off-diagonal elements of an intermediate tridiagonal
122 *> form did not converge to zero.
123 *> \endverbatim
124 *
125 * Authors:
126 * ========
127 *
128 *> \author Univ. of Tennessee
129 *> \author Univ. of California Berkeley
130 *> \author Univ. of Colorado Denver
131 *> \author NAG Ltd.
132 *
133 *> \date November 2011
134 *
135 *> \ingroup complexOTHEReigen
136 *
137 * =====================================================================
138  SUBROUTINE chpev( JOBZ, UPLO, N, AP, W, Z, LDZ, WORK, RWORK,
139  $ info )
140 *
141 * -- LAPACK driver routine (version 3.4.0) --
142 * -- LAPACK is a software package provided by Univ. of Tennessee, --
143 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
144 * November 2011
145 *
146 * .. Scalar Arguments ..
147  CHARACTER jobz, uplo
148  INTEGER info, ldz, n
149 * ..
150 * .. Array Arguments ..
151  REAL rwork( * ), w( * )
152  COMPLEX ap( * ), work( * ), z( ldz, * )
153 * ..
154 *
155 * =====================================================================
156 *
157 * .. Parameters ..
158  REAL zero, one
159  parameter( zero = 0.0e0, one = 1.0e0 )
160 * ..
161 * .. Local Scalars ..
162  LOGICAL wantz
163  INTEGER iinfo, imax, inde, indrwk, indtau, indwrk,
164  $ iscale
165  REAL anrm, bignum, eps, rmax, rmin, safmin, sigma,
166  $ smlnum
167 * ..
168 * .. External Functions ..
169  LOGICAL lsame
170  REAL clanhp, slamch
171  EXTERNAL lsame, clanhp, slamch
172 * ..
173 * .. External Subroutines ..
174  EXTERNAL chptrd, csscal, csteqr, cupgtr, sscal, ssterf,
175  $ xerbla
176 * ..
177 * .. Intrinsic Functions ..
178  INTRINSIC sqrt
179 * ..
180 * .. Executable Statements ..
181 *
182 * Test the input parameters.
183 *
184  wantz = lsame( jobz, 'V' )
185 *
186  info = 0
187  IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
188  info = -1
189  ELSE IF( .NOT.( lsame( uplo, 'L' ) .OR. lsame( uplo, 'U' ) ) )
190  $ THEN
191  info = -2
192  ELSE IF( n.LT.0 ) THEN
193  info = -3
194  ELSE IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) ) THEN
195  info = -7
196  END IF
197 *
198  IF( info.NE.0 ) THEN
199  CALL xerbla( 'CHPEV ', -info )
200  return
201  END IF
202 *
203 * Quick return if possible
204 *
205  IF( n.EQ.0 )
206  $ return
207 *
208  IF( n.EQ.1 ) THEN
209  w( 1 ) = ap( 1 )
210  rwork( 1 ) = 1
211  IF( wantz )
212  $ z( 1, 1 ) = one
213  return
214  END IF
215 *
216 * Get machine constants.
217 *
218  safmin = slamch( 'Safe minimum' )
219  eps = slamch( 'Precision' )
220  smlnum = safmin / eps
221  bignum = one / smlnum
222  rmin = sqrt( smlnum )
223  rmax = sqrt( bignum )
224 *
225 * Scale matrix to allowable range, if necessary.
226 *
227  anrm = clanhp( 'M', uplo, n, ap, rwork )
228  iscale = 0
229  IF( anrm.GT.zero .AND. anrm.LT.rmin ) THEN
230  iscale = 1
231  sigma = rmin / anrm
232  ELSE IF( anrm.GT.rmax ) THEN
233  iscale = 1
234  sigma = rmax / anrm
235  END IF
236  IF( iscale.EQ.1 ) THEN
237  CALL csscal( ( n*( n+1 ) ) / 2, sigma, ap, 1 )
238  END IF
239 *
240 * Call CHPTRD to reduce Hermitian packed matrix to tridiagonal form.
241 *
242  inde = 1
243  indtau = 1
244  CALL chptrd( uplo, n, ap, w, rwork( inde ), work( indtau ),
245  $ iinfo )
246 *
247 * For eigenvalues only, call SSTERF. For eigenvectors, first call
248 * CUPGTR to generate the orthogonal matrix, then call CSTEQR.
249 *
250  IF( .NOT.wantz ) THEN
251  CALL ssterf( n, w, rwork( inde ), info )
252  ELSE
253  indwrk = indtau + n
254  CALL cupgtr( uplo, n, ap, work( indtau ), z, ldz,
255  $ work( indwrk ), iinfo )
256  indrwk = inde + n
257  CALL csteqr( jobz, n, w, rwork( inde ), z, ldz,
258  $ rwork( indrwk ), info )
259  END IF
260 *
261 * If matrix was scaled, then rescale eigenvalues appropriately.
262 *
263  IF( iscale.EQ.1 ) THEN
264  IF( info.EQ.0 ) THEN
265  imax = n
266  ELSE
267  imax = info - 1
268  END IF
269  CALL sscal( imax, one / sigma, w, 1 )
270  END IF
271 *
272  return
273 *
274 * End of CHPEV
275 *
276  END