 LAPACK  3.10.0 LAPACK: Linear Algebra PACKage

◆ chsein()

 subroutine chsein ( character SIDE, character EIGSRC, character INITV, logical, dimension( * ) SELECT, integer N, complex, dimension( ldh, * ) H, integer LDH, complex, dimension( * ) W, complex, dimension( ldvl, * ) VL, integer LDVL, complex, dimension( ldvr, * ) VR, integer LDVR, integer MM, integer M, complex, dimension( * ) WORK, real, dimension( * ) RWORK, integer, dimension( * ) IFAILL, integer, dimension( * ) IFAILR, integer INFO )

CHSEIN

Purpose:
CHSEIN uses inverse iteration to find specified right and/or left
eigenvectors of a complex upper Hessenberg matrix H.

The right eigenvector x and the left eigenvector y of the matrix H
corresponding to an eigenvalue w are defined by:

H * x = w * x,     y**h * H = w * y**h

where y**h denotes the conjugate transpose of the vector y.
Parameters
 [in] SIDE SIDE is CHARACTER*1 = 'R': compute right eigenvectors only; = 'L': compute left eigenvectors only; = 'B': compute both right and left eigenvectors. [in] EIGSRC EIGSRC is CHARACTER*1 Specifies the source of eigenvalues supplied in W: = 'Q': the eigenvalues were found using CHSEQR; thus, if H has zero subdiagonal elements, and so is block-triangular, then the j-th eigenvalue can be assumed to be an eigenvalue of the block containing the j-th row/column. This property allows CHSEIN to perform inverse iteration on just one diagonal block. = 'N': no assumptions are made on the correspondence between eigenvalues and diagonal blocks. In this case, CHSEIN must always perform inverse iteration using the whole matrix H. [in] INITV INITV is CHARACTER*1 = 'N': no initial vectors are supplied; = 'U': user-supplied initial vectors are stored in the arrays VL and/or VR. [in] SELECT SELECT is LOGICAL array, dimension (N) Specifies the eigenvectors to be computed. To select the eigenvector corresponding to the eigenvalue W(j), SELECT(j) must be set to .TRUE.. [in] N N is INTEGER The order of the matrix H. N >= 0. [in] H H is COMPLEX array, dimension (LDH,N) The upper Hessenberg matrix H. If a NaN is detected in H, the routine will return with INFO=-6. [in] LDH LDH is INTEGER The leading dimension of the array H. LDH >= max(1,N). [in,out] W W is COMPLEX array, dimension (N) On entry, the eigenvalues of H. On exit, the real parts of W may have been altered since close eigenvalues are perturbed slightly in searching for independent eigenvectors. [in,out] VL VL is COMPLEX array, dimension (LDVL,MM) On entry, if INITV = 'U' and SIDE = 'L' or 'B', VL must contain starting vectors for the inverse iteration for the left eigenvectors; the starting vector for each eigenvector must be in the same column in which the eigenvector will be stored. On exit, if SIDE = 'L' or 'B', the left eigenvectors specified by SELECT will be stored consecutively in the columns of VL, in the same order as their eigenvalues. If SIDE = 'R', VL is not referenced. [in] LDVL LDVL is INTEGER The leading dimension of the array VL. LDVL >= max(1,N) if SIDE = 'L' or 'B'; LDVL >= 1 otherwise. [in,out] VR VR is COMPLEX array, dimension (LDVR,MM) On entry, if INITV = 'U' and SIDE = 'R' or 'B', VR must contain starting vectors for the inverse iteration for the right eigenvectors; the starting vector for each eigenvector must be in the same column in which the eigenvector will be stored. On exit, if SIDE = 'R' or 'B', the right eigenvectors specified by SELECT will be stored consecutively in the columns of VR, in the same order as their eigenvalues. If SIDE = 'L', VR is not referenced. [in] LDVR LDVR is INTEGER The leading dimension of the array VR. LDVR >= max(1,N) if SIDE = 'R' or 'B'; LDVR >= 1 otherwise. [in] MM MM is INTEGER The number of columns in the arrays VL and/or VR. MM >= M. [out] M M is INTEGER The number of columns in the arrays VL and/or VR required to store the eigenvectors (= the number of .TRUE. elements in SELECT). [out] WORK WORK is COMPLEX array, dimension (N*N) [out] RWORK RWORK is REAL array, dimension (N) [out] IFAILL IFAILL is INTEGER array, dimension (MM) If SIDE = 'L' or 'B', IFAILL(i) = j > 0 if the left eigenvector in the i-th column of VL (corresponding to the eigenvalue w(j)) failed to converge; IFAILL(i) = 0 if the eigenvector converged satisfactorily. If SIDE = 'R', IFAILL is not referenced. [out] IFAILR IFAILR is INTEGER array, dimension (MM) If SIDE = 'R' or 'B', IFAILR(i) = j > 0 if the right eigenvector in the i-th column of VR (corresponding to the eigenvalue w(j)) failed to converge; IFAILR(i) = 0 if the eigenvector converged satisfactorily. If SIDE = 'L', IFAILR is not referenced. [out] INFO INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value > 0: if INFO = i, i is the number of eigenvectors which failed to converge; see IFAILL and IFAILR for further details.
Further Details:
Each eigenvector is normalized so that the element of largest
magnitude has magnitude 1; here the magnitude of a complex number
(x,y) is taken to be |x|+|y|.

Definition at line 242 of file chsein.f.

245 *
246 * -- LAPACK computational routine --
247 * -- LAPACK is a software package provided by Univ. of Tennessee, --
248 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
249 *
250 * .. Scalar Arguments ..
251  CHARACTER EIGSRC, INITV, SIDE
252  INTEGER INFO, LDH, LDVL, LDVR, M, MM, N
253 * ..
254 * .. Array Arguments ..
255  LOGICAL SELECT( * )
256  INTEGER IFAILL( * ), IFAILR( * )
257  REAL RWORK( * )
258  COMPLEX H( LDH, * ), VL( LDVL, * ), VR( LDVR, * ),
259  \$ W( * ), WORK( * )
260 * ..
261 *
262 * =====================================================================
263 *
264 * .. Parameters ..
265  COMPLEX ZERO
266  parameter( zero = ( 0.0e+0, 0.0e+0 ) )
267  REAL RZERO
268  parameter( rzero = 0.0e+0 )
269 * ..
270 * .. Local Scalars ..
271  LOGICAL BOTHV, FROMQR, LEFTV, NOINIT, RIGHTV
272  INTEGER I, IINFO, K, KL, KLN, KR, KS, LDWORK
273  REAL EPS3, HNORM, SMLNUM, ULP, UNFL
274  COMPLEX CDUM, WK
275 * ..
276 * .. External Functions ..
277  LOGICAL LSAME, SISNAN
278  REAL CLANHS, SLAMCH
279  EXTERNAL lsame, clanhs, slamch, sisnan
280 * ..
281 * .. External Subroutines ..
282  EXTERNAL claein, xerbla
283 * ..
284 * .. Intrinsic Functions ..
285  INTRINSIC abs, aimag, max, real
286 * ..
287 * .. Statement Functions ..
288  REAL CABS1
289 * ..
290 * .. Statement Function definitions ..
291  cabs1( cdum ) = abs( real( cdum ) ) + abs( aimag( cdum ) )
292 * ..
293 * .. Executable Statements ..
294 *
295 * Decode and test the input parameters.
296 *
297  bothv = lsame( side, 'B' )
298  rightv = lsame( side, 'R' ) .OR. bothv
299  leftv = lsame( side, 'L' ) .OR. bothv
300 *
301  fromqr = lsame( eigsrc, 'Q' )
302 *
303  noinit = lsame( initv, 'N' )
304 *
305 * Set M to the number of columns required to store the selected
306 * eigenvectors.
307 *
308  m = 0
309  DO 10 k = 1, n
310  IF( SELECT( k ) )
311  \$ m = m + 1
312  10 CONTINUE
313 *
314  info = 0
315  IF( .NOT.rightv .AND. .NOT.leftv ) THEN
316  info = -1
317  ELSE IF( .NOT.fromqr .AND. .NOT.lsame( eigsrc, 'N' ) ) THEN
318  info = -2
319  ELSE IF( .NOT.noinit .AND. .NOT.lsame( initv, 'U' ) ) THEN
320  info = -3
321  ELSE IF( n.LT.0 ) THEN
322  info = -5
323  ELSE IF( ldh.LT.max( 1, n ) ) THEN
324  info = -7
325  ELSE IF( ldvl.LT.1 .OR. ( leftv .AND. ldvl.LT.n ) ) THEN
326  info = -10
327  ELSE IF( ldvr.LT.1 .OR. ( rightv .AND. ldvr.LT.n ) ) THEN
328  info = -12
329  ELSE IF( mm.LT.m ) THEN
330  info = -13
331  END IF
332  IF( info.NE.0 ) THEN
333  CALL xerbla( 'CHSEIN', -info )
334  RETURN
335  END IF
336 *
337 * Quick return if possible.
338 *
339  IF( n.EQ.0 )
340  \$ RETURN
341 *
342 * Set machine-dependent constants.
343 *
344  unfl = slamch( 'Safe minimum' )
345  ulp = slamch( 'Precision' )
346  smlnum = unfl*( n / ulp )
347 *
348  ldwork = n
349 *
350  kl = 1
351  kln = 0
352  IF( fromqr ) THEN
353  kr = 0
354  ELSE
355  kr = n
356  END IF
357  ks = 1
358 *
359  DO 100 k = 1, n
360  IF( SELECT( k ) ) THEN
361 *
362 * Compute eigenvector(s) corresponding to W(K).
363 *
364  IF( fromqr ) THEN
365 *
366 * If affiliation of eigenvalues is known, check whether
367 * the matrix splits.
368 *
369 * Determine KL and KR such that 1 <= KL <= K <= KR <= N
370 * and H(KL,KL-1) and H(KR+1,KR) are zero (or KL = 1 or
371 * KR = N).
372 *
373 * Then inverse iteration can be performed with the
374 * submatrix H(KL:N,KL:N) for a left eigenvector, and with
375 * the submatrix H(1:KR,1:KR) for a right eigenvector.
376 *
377  DO 20 i = k, kl + 1, -1
378  IF( h( i, i-1 ).EQ.zero )
379  \$ GO TO 30
380  20 CONTINUE
381  30 CONTINUE
382  kl = i
383  IF( k.GT.kr ) THEN
384  DO 40 i = k, n - 1
385  IF( h( i+1, i ).EQ.zero )
386  \$ GO TO 50
387  40 CONTINUE
388  50 CONTINUE
389  kr = i
390  END IF
391  END IF
392 *
393  IF( kl.NE.kln ) THEN
394  kln = kl
395 *
396 * Compute infinity-norm of submatrix H(KL:KR,KL:KR) if it
397 * has not ben computed before.
398 *
399  hnorm = clanhs( 'I', kr-kl+1, h( kl, kl ), ldh, rwork )
400  IF( sisnan( hnorm ) ) THEN
401  info = -6
402  RETURN
403  ELSE IF( (hnorm.GT.rzero) ) THEN
404  eps3 = hnorm*ulp
405  ELSE
406  eps3 = smlnum
407  END IF
408  END IF
409 *
410 * Perturb eigenvalue if it is close to any previous
411 * selected eigenvalues affiliated to the submatrix
412 * H(KL:KR,KL:KR). Close roots are modified by EPS3.
413 *
414  wk = w( k )
415  60 CONTINUE
416  DO 70 i = k - 1, kl, -1
417  IF( SELECT( i ) .AND. cabs1( w( i )-wk ).LT.eps3 ) THEN
418  wk = wk + eps3
419  GO TO 60
420  END IF
421  70 CONTINUE
422  w( k ) = wk
423 *
424  IF( leftv ) THEN
425 *
426 * Compute left eigenvector.
427 *
428  CALL claein( .false., noinit, n-kl+1, h( kl, kl ), ldh,
429  \$ wk, vl( kl, ks ), work, ldwork, rwork, eps3,
430  \$ smlnum, iinfo )
431  IF( iinfo.GT.0 ) THEN
432  info = info + 1
433  ifaill( ks ) = k
434  ELSE
435  ifaill( ks ) = 0
436  END IF
437  DO 80 i = 1, kl - 1
438  vl( i, ks ) = zero
439  80 CONTINUE
440  END IF
441  IF( rightv ) THEN
442 *
443 * Compute right eigenvector.
444 *
445  CALL claein( .true., noinit, kr, h, ldh, wk, vr( 1, ks ),
446  \$ work, ldwork, rwork, eps3, smlnum, iinfo )
447  IF( iinfo.GT.0 ) THEN
448  info = info + 1
449  ifailr( ks ) = k
450  ELSE
451  ifailr( ks ) = 0
452  END IF
453  DO 90 i = kr + 1, n
454  vr( i, ks ) = zero
455  90 CONTINUE
456  END IF
457  ks = ks + 1
458  END IF
459  100 CONTINUE
460 *
461  RETURN
462 *
463 * End of CHSEIN
464 *
logical function sisnan(SIN)
SISNAN tests input for NaN.
Definition: sisnan.f:59
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
real function clanhs(NORM, N, A, LDA, WORK)
CLANHS returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: clanhs.f:109
subroutine claein(RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, RWORK, EPS3, SMLNUM, INFO)
CLAEIN computes a specified right or left eigenvector of an upper Hessenberg matrix by inverse iterat...
Definition: claein.f:149
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:68
Here is the call graph for this function:
Here is the caller graph for this function: