LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine claein ( logical  RIGHTV,
logical  NOINIT,
integer  N,
complex, dimension( ldh, * )  H,
integer  LDH,
complex  W,
complex, dimension( * )  V,
complex, dimension( ldb, * )  B,
integer  LDB,
real, dimension( * )  RWORK,
real  EPS3,
real  SMLNUM,
integer  INFO 
)

CLAEIN computes a specified right or left eigenvector of an upper Hessenberg matrix by inverse iteration.

Download CLAEIN + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 CLAEIN uses inverse iteration to find a right or left eigenvector
 corresponding to the eigenvalue W of a complex upper Hessenberg
 matrix H.
Parameters
[in]RIGHTV
          RIGHTV is LOGICAL
          = .TRUE. : compute right eigenvector;
          = .FALSE.: compute left eigenvector.
[in]NOINIT
          NOINIT is LOGICAL
          = .TRUE. : no initial vector supplied in V
          = .FALSE.: initial vector supplied in V.
[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.
[in]LDH
          LDH is INTEGER
          The leading dimension of the array H.  LDH >= max(1,N).
[in]W
          W is COMPLEX
          The eigenvalue of H whose corresponding right or left
          eigenvector is to be computed.
[in,out]V
          V is COMPLEX array, dimension (N)
          On entry, if NOINIT = .FALSE., V must contain a starting
          vector for inverse iteration; otherwise V need not be set.
          On exit, V contains the computed eigenvector, normalized so
          that the component of largest magnitude has magnitude 1; here
          the magnitude of a complex number (x,y) is taken to be
          |x| + |y|.
[out]B
          B is COMPLEX array, dimension (LDB,N)
[in]LDB
          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).
[out]RWORK
          RWORK is REAL array, dimension (N)
[in]EPS3
          EPS3 is REAL
          A small machine-dependent value which is used to perturb
          close eigenvalues, and to replace zero pivots.
[in]SMLNUM
          SMLNUM is REAL
          A machine-dependent value close to the underflow threshold.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          = 1:  inverse iteration did not converge; V is set to the
                last iterate.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012

Definition at line 151 of file claein.f.

151 *
152 * -- LAPACK auxiliary routine (version 3.4.2) --
153 * -- LAPACK is a software package provided by Univ. of Tennessee, --
154 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
155 * September 2012
156 *
157 * .. Scalar Arguments ..
158  LOGICAL noinit, rightv
159  INTEGER info, ldb, ldh, n
160  REAL eps3, smlnum
161  COMPLEX w
162 * ..
163 * .. Array Arguments ..
164  REAL rwork( * )
165  COMPLEX b( ldb, * ), h( ldh, * ), v( * )
166 * ..
167 *
168 * =====================================================================
169 *
170 * .. Parameters ..
171  REAL one, tenth
172  parameter ( one = 1.0e+0, tenth = 1.0e-1 )
173  COMPLEX zero
174  parameter ( zero = ( 0.0e+0, 0.0e+0 ) )
175 * ..
176 * .. Local Scalars ..
177  CHARACTER normin, trans
178  INTEGER i, ierr, its, j
179  REAL growto, nrmsml, rootn, rtemp, scale, vnorm
180  COMPLEX cdum, ei, ej, temp, x
181 * ..
182 * .. External Functions ..
183  INTEGER icamax
184  REAL scasum, scnrm2
185  COMPLEX cladiv
186  EXTERNAL icamax, scasum, scnrm2, cladiv
187 * ..
188 * .. External Subroutines ..
189  EXTERNAL clatrs, csscal
190 * ..
191 * .. Intrinsic Functions ..
192  INTRINSIC abs, aimag, max, REAL, sqrt
193 * ..
194 * .. Statement Functions ..
195  REAL cabs1
196 * ..
197 * .. Statement Function definitions ..
198  cabs1( cdum ) = abs( REAL( CDUM ) ) + abs( aimag( cdum ) )
199 * ..
200 * .. Executable Statements ..
201 *
202  info = 0
203 *
204 * GROWTO is the threshold used in the acceptance test for an
205 * eigenvector.
206 *
207  rootn = sqrt( REAL( N ) )
208  growto = tenth / rootn
209  nrmsml = max( one, eps3*rootn )*smlnum
210 *
211 * Form B = H - W*I (except that the subdiagonal elements are not
212 * stored).
213 *
214  DO 20 j = 1, n
215  DO 10 i = 1, j - 1
216  b( i, j ) = h( i, j )
217  10 CONTINUE
218  b( j, j ) = h( j, j ) - w
219  20 CONTINUE
220 *
221  IF( noinit ) THEN
222 *
223 * Initialize V.
224 *
225  DO 30 i = 1, n
226  v( i ) = eps3
227  30 CONTINUE
228  ELSE
229 *
230 * Scale supplied initial vector.
231 *
232  vnorm = scnrm2( n, v, 1 )
233  CALL csscal( n, ( eps3*rootn ) / max( vnorm, nrmsml ), v, 1 )
234  END IF
235 *
236  IF( rightv ) THEN
237 *
238 * LU decomposition with partial pivoting of B, replacing zero
239 * pivots by EPS3.
240 *
241  DO 60 i = 1, n - 1
242  ei = h( i+1, i )
243  IF( cabs1( b( i, i ) ).LT.cabs1( ei ) ) THEN
244 *
245 * Interchange rows and eliminate.
246 *
247  x = cladiv( b( i, i ), ei )
248  b( i, i ) = ei
249  DO 40 j = i + 1, n
250  temp = b( i+1, j )
251  b( i+1, j ) = b( i, j ) - x*temp
252  b( i, j ) = temp
253  40 CONTINUE
254  ELSE
255 *
256 * Eliminate without interchange.
257 *
258  IF( b( i, i ).EQ.zero )
259  $ b( i, i ) = eps3
260  x = cladiv( ei, b( i, i ) )
261  IF( x.NE.zero ) THEN
262  DO 50 j = i + 1, n
263  b( i+1, j ) = b( i+1, j ) - x*b( i, j )
264  50 CONTINUE
265  END IF
266  END IF
267  60 CONTINUE
268  IF( b( n, n ).EQ.zero )
269  $ b( n, n ) = eps3
270 *
271  trans = 'N'
272 *
273  ELSE
274 *
275 * UL decomposition with partial pivoting of B, replacing zero
276 * pivots by EPS3.
277 *
278  DO 90 j = n, 2, -1
279  ej = h( j, j-1 )
280  IF( cabs1( b( j, j ) ).LT.cabs1( ej ) ) THEN
281 *
282 * Interchange columns and eliminate.
283 *
284  x = cladiv( b( j, j ), ej )
285  b( j, j ) = ej
286  DO 70 i = 1, j - 1
287  temp = b( i, j-1 )
288  b( i, j-1 ) = b( i, j ) - x*temp
289  b( i, j ) = temp
290  70 CONTINUE
291  ELSE
292 *
293 * Eliminate without interchange.
294 *
295  IF( b( j, j ).EQ.zero )
296  $ b( j, j ) = eps3
297  x = cladiv( ej, b( j, j ) )
298  IF( x.NE.zero ) THEN
299  DO 80 i = 1, j - 1
300  b( i, j-1 ) = b( i, j-1 ) - x*b( i, j )
301  80 CONTINUE
302  END IF
303  END IF
304  90 CONTINUE
305  IF( b( 1, 1 ).EQ.zero )
306  $ b( 1, 1 ) = eps3
307 *
308  trans = 'C'
309 *
310  END IF
311 *
312  normin = 'N'
313  DO 110 its = 1, n
314 *
315 * Solve U*x = scale*v for a right eigenvector
316 * or U**H *x = scale*v for a left eigenvector,
317 * overwriting x on v.
318 *
319  CALL clatrs( 'Upper', trans, 'Nonunit', normin, n, b, ldb, v,
320  $ scale, rwork, ierr )
321  normin = 'Y'
322 *
323 * Test for sufficient growth in the norm of v.
324 *
325  vnorm = scasum( n, v, 1 )
326  IF( vnorm.GE.growto*scale )
327  $ GO TO 120
328 *
329 * Choose new orthogonal starting vector and try again.
330 *
331  rtemp = eps3 / ( rootn+one )
332  v( 1 ) = eps3
333  DO 100 i = 2, n
334  v( i ) = rtemp
335  100 CONTINUE
336  v( n-its+1 ) = v( n-its+1 ) - eps3*rootn
337  110 CONTINUE
338 *
339 * Failure to find eigenvector in N iterations.
340 *
341  info = 1
342 *
343  120 CONTINUE
344 *
345 * Normalize eigenvector.
346 *
347  i = icamax( n, v, 1 )
348  CALL csscal( n, one / cabs1( v( i ) ), v, 1 )
349 *
350  RETURN
351 *
352 * End of CLAEIN
353 *
real function scnrm2(N, X, INCX)
SCNRM2
Definition: scnrm2.f:56
real function scasum(N, CX, INCX)
SCASUM
Definition: scasum.f:54
subroutine clatrs(UPLO, TRANS, DIAG, NORMIN, N, A, LDA, X, SCALE, CNORM, INFO)
CLATRS solves a triangular system of equations with the scale factor set to prevent overflow...
Definition: clatrs.f:241
integer function icamax(N, CX, INCX)
ICAMAX
Definition: icamax.f:53
complex function cladiv(X, Y)
CLADIV performs complex division in real arithmetic, avoiding unnecessary overflow.
Definition: cladiv.f:66
subroutine csscal(N, SA, CX, INCX)
CSSCAL
Definition: csscal.f:54

Here is the call graph for this function:

Here is the caller graph for this function: